MATLAB Examples

Nested Function for Geometry with Additional Parameters

This example shows how to include additional parameters into a function for creating a 2-D geometry.

When a 2-D geometry function requires additional parameters, you cannot use a standard anonymous function approach because geometry functions return a varying number of arguments. Instead, you can use global variables or nested functions. In most cases, the recommended approach is to use nested functions.

The example solves a Poisson's equation with zero Dirichlet boundary conditions on all boundaries. The geometry is a cardioid with an elliptic hole that has a center at (1,-1) and variable semiaxes. To set up and solve the PDE model with this geometry, use a nested function. Here, the parent function accepts the lengths of the semiaxes, rr and ss, as input parameters. The reason to nest cardioidWithEllipseGeom within cardioidWithEllipseModel is that nested functions share the workspace of their parent functions. Therefore, the cardioidWithEllipseGeom function can access the values of rr and ss that you pass to cardioidWithEllipseModel.

`function cardioidWithEllipseModel(rr,ss)`
```if (rr > 0) & (ss > 0)
model = createpde();
geometryFromEdges(model,@cardioidWithEllipseGeom);
pdegplot(model,'EdgeLabels','on','FaceLabels','on')
axis equal```
```  applyBoundaryCondition(model,'dirichlet','Edge',1:8,'u',0);
specifyCoefficients(model,'m',0,'d',0,'c',1,'a',0,'f',1);```
```  generateMesh(model);
u = solvepde(model);
figure
pdeplot(model,'XYData',u.NodalSolution)
axis equal```
```else
display('Semiaxes values must be positive numbers.')
end```
`function [x,y] = cardioidWithEllipseGeom(bs,s)`
```if nargin == 0
x = 8; % eight segments in boundary
return
end```
```if nargin == 1
% Cardioid
dlc = [ 0   4   8  12
4   8  12  16
1   1   1   1
0   0   0   0];
% Ellipse
dle = [0      pi/2   pi       3*pi/2
pi/2   pi     3*pi/2   2*pi
0      0      0        0
1      1      1        1];
% Combine the edge matrices
dl = [dlc,dle];
x = dl(:,bs);
return
end```
```x = zeros(size(s));
y = zeros(size(s));
if numel(bs) == 1 % Does bs need scalar expansion?
bs = bs*ones(size(s)); % Expand bs
end```
```cbs = find(bs < 3); % Upper half of cardioid
x(cbs) = s(cbs).^4/512 - 3*s(cbs).^2/16 + 4;
y(cbs) = s(cbs).*(64 - s(cbs).^2).^(3/2)/512;
cbs = find(bs >= 3 & bs <= 4); % Lower half of cardioid
s(cbs) = 16 - s(cbs);
x(cbs) = s(cbs).^4/512 - 3*s(cbs).^2/16 + 4;
y(cbs) = -s(cbs).*(64 - s(cbs).^2).^(3/2)/512;
cbs = find(bs > 4); % Inner ellipse center (1,-1) axes rr and ss
x(cbs) = 1 + rr*cos(s(cbs));
y(cbs) = -1 + ss*sin(s(cbs));```
`end`
`end`

When calling cardioidWithEllipseModel, ensure that the semiaxes values are small enough, so that the elliptic hole appears entirely within the outer cardioid. Otherwise, the geometry becomes invalid.

For example, call the function for the ellipse with the major semiaxis rr = 0.5 and the minor semiaxis ss = 0.25. This function call returns the following geometry and the solution.

```cardioidWithEllipseModel(0.5,0.25) ```