% ------ These codes show you how to construct and mesh complex geometry
% ------ using command lines, very useful for solving
% ------ PDEs in continous domains
% ------ by Hongxue Cai @ Northwestern University -------------->
% These code need MATLAB PDE Toolbox
%x, y-coordinates of the geometrical elements which dfining the geometrical
%domain. These values come from an application which models the tectorial
%membrane in the inner ear. These should be different for different
%applications, therefore readers should not pay much attention to the
%individual numbers.
xy = [
-0.01430000000000 0.00255766422613
-0.00954483106568 0.00776466551413
-0.00411900762330 0.00903776158615
0.00152629824713 0.00978098020364
0.00236820179834 0.00774844523211
-0.01760000000000 0.00187878393617
-0.01760000000000 -0.00013871129459
-0.01474101729286 0.00047264657699
-0.01122101729286 0.00031358294595
];
% x, y- coordinates and radius of the center point for circular arcs
xycR = [
-0.01760000000000 0.01023878393617 0.00836000000000
-0.00605407636840 0.00508611522649 0.00440000000000
-0.00268321950888 0.01994365506126 0.01100000000000
0.00194725002273 0.00876471271787 0.00110000000000
-0.01354587812765 -0.01210697797173 0.01263626967942
];
x = xy(:, 1);
y = xy(:, 2);
xc = xycR(:, 1);
yc = xycR(:, 2);
Rc = xycR(:, 3);
%geometry matrix, in which, each colum defines a geometrical element that
%can be either a line or a circular arc. (Note that any complex geometry can
%be composed by line and circular arts.) In each colum, the number in row#1
%indicates the type of the geometrica element, 2 for line element, 2 for
%circular arts. Then, row#2 is the x-coordinate of the starting point;
%row#3, that of the ending point; Row #4 and #5 are, respectively the
%y-coordinates of the the starting and ending point. The number in the 6th
%anf 7th
%rows indicate the location of geometric domain, (0,1) on the left side of this
%geometric element, (1,0) on th right. The rows 8, 9, 10 are zero for line,
%and repesent the x, y-coordinates and the radius of the center point for a
%circular arc.
gTM = [
2 1 1 1 2 2 1 1 1
x(1) x(3) x(3) x(5) x(9) x(7) x(6) x(9) x(8)
x(2) x(2) x(4) x(4) x(5) x(6) x(1) x(8) x(7)
y(1) y(3) y(3) y(5) y(9) y(7) y(6) y(9) y(8)
y(2) y(2) y(4) y(4) y(5) y(6) y(1) y(8) y(7)
0 1 0 1 1 0 0 0 0
1 0 1 0 0 1 1 1 1
0 xc(2) xc(3) xc(4) 0 0 xc(1) xc(5) xc(5)
0 yc(2) yc(3) yc(4) 0 0 yc(1) yc(5) yc(5)
0 Rc(2) Rc(3) Rc(4) 0 0 Rc(1) Rc(5) Rc(5)
];
%plotting geometry
figure
pdegplot(gTM)
%meshing
[pTM,eTM,tTM] = initmesh(gTM);
[pTM,eTM,tTM] = refinemesh(gTM,pTM,eTM,tTM);
%plotting mesh
figure
pdemesh(pTM,eTM,tTM);