Support for Neumann and Robin boundary conditions

[danfortunato]

This feature request would allow the user to specify boundary conditions as either Dirichlet, Neumann, or Robin. Syntax would need to be introduced in order to specify which conditions and boundary data correspond to which boundaries of the given Domain. The boundary conditions would also need to conform to the element boundaries to avoid hanging nodes.

This could be implemented at the top level of the ultraSEM hierarchy by mapping the given boundary data (whether Dirichlet, Neumann, of Robin) to pure Dirichlet data. This involves solving an O(p) x O(p) linear system using the top-level Dirichlet-to-Neumann map.

[danfortunato]

Here is an example of imposing a Neumann boundary condition using the top-level Dirichlet-to-Neumann map. The result is a vector dir of coefficients of pure Dirichlet data:

n = 20;
rhs = randnfun2(1, [-2 2 -1 1]);
pdo = {1,0,0};
dom = ultraSEM.rectangle([-2 2 -1 1], 1, 2);
S = ultraSEM(dom, pdo, rhs, n);

% Boundary data (in coefficient form):
dir_cfs = [1; zeros(n-1,1)]; % Dirichlet data: u = 1
neu_cfs = zeros(n,1);        % Neumann data: du/dn = 0

% Build the ultraSEM so that we can get the top-level D2N map:
build(S)
D2N = S.patches{1}.D2N;

% Let's impose Neumann conditions on the right side
side = 3*n+1:4*n;
bc = repmat(dir_cfs, 6, 1);
bc(side) = neu_cfs - D2N(side,end);

% Map the mixed boundary data to pure Dirichlet data:
A = speye(6*n);
A(side,:) = D2N(side,1:end-1);
dir = A \ bc;

[danfortunato]

The PR in #37 partly implements this. It does not support pure Neumann BCs with a mean zero condition yet.