Add a getBoundaryOperator method in multiblock.DiffOp that uses the the diffop get_boundary_operator method rather than its properties. Werpers might not like this.

function z = gaussian(x,x0,d)
    z = exp(-sum((x-x0).^2,2)/d^2);
end