Data Structure: Boundary Conditions

We use bdFlag(1:NT,1:3) to record the type of three edges of each triangle. Similarly in 3-D, we use bdFlag(1:NT,1:4) to record the type of four faces of each tetrahedron. The value is the type of boundary condition.

0: non-boundary, i.e., an interior edge or face.
1: first type, i.e., a Dirichlet boundary edge or face.
2: second type, i.e., a Neumann boundary edge or face.
3: third type, i.e., a Robin boundary edge or face.

Local labeling of edges and faces
Set up boundary conditions
Example: Crack Domain
Example: Prism Domain
Remark

Local labeling of edges and faces

We label three edges of a triangle such that bdFlag(t,i) is the edge opposite to the i-th vertex. Similarly bdFlag(t,i) is the face opposite to the i-th vertex.

node = [1,0; 1,1; 0,0];
elem = [1 2 3];
edge = [2 3; 1 3; 1 2];
showmesh(node,elem);
findnode(node);
findedge(node,edge,'all','vec');

Set up boundary conditions

The function setboundary is to set up the bdFlag matrix for a 2-D triangulation and setboundary3 for a 3-D triangulation.

help setboundary

  SETBOUNDARY set type of boundary edges.
 
   bdFlag = SETBOUNDARY(node,elem,'Dirichlet') set all boundary edges to
   Dirichlet type. 
 
   bdFlag = SETBOUNDARY(node,elem,'Neumann') set all boundary edges to
   Neumann type. 
 
   bdFlag = SETBOUNDARY(node,elem,'Robin') set all boundary edges to
   Robin type. 
 
   bdFlag = SETBOUNDARY(node,elem,'Dirichlet','(x==1) | (x==-1)') set
   Dirichlet boundary condition on x=1 and x=-1. Other edges are
   homongenous Neumann boundary condition.
 
   bdFlag = SETBOUNDARY(node,elem,'Dirichlet','(x==1) | ...
   (x==-1)','Neumann','(y==1) | (y==-1)') set
   Dirichlet boundary condition on x=1 or x=-1 and Neumann boundary
   condition on y=1 or y=-1.
 
   bdFlag = SETBOUNDARY(node,elem,'Dirichlet','(x==1) | ...
   (x==-1)','Neumann','y==1', 'Robin',' y==-1') set
   Dirichlet boundary condition on x=1 or x=-1 and Neumann boundary
   condition on y=1, and Robin boundary condition on y=-1.
 
   bdFlag = SETBOUNDARY(node,elem,'Dirichlet','all','Neumann','y==1') set
   Neumann boundary condition on y=1 and others are Dirichlet boundary condition.
 
  Example
    
       node = [0,0; 1,0; 1,1; 0,1];
       elem = [2,3,1; 4,1,3];
       bdFlag = setboundary(node,elem,'Dirichlet','all','Neumann','y==1');
       [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
       [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
       showmesh(node,elem);
       allEdge = [elem(:,[2,3]); elem(:,[3,1]); elem(:,[1,2])];
       Dirichlet = allEdge((bdFlag(:) == 1),:);
       Neumann = allEdge((bdFlag(:) == 2) | (bdFlag(:) == 3),:);
       findedge(node,Dirichlet,[],'noindex','LineWidth',4,'Color','r');
       findedge(node,Neumann,[],'noindex','LineWidth',4,'Color','b');
 
  See also setboundary3
 
  Copyright (C) Long Chen. See COPYRIGHT.txt for details.

Note that if the i-th edge of t is on the boundary but bdFlag(t,i)=0, it is equivalent to use homogenous Neumann boundary condition (zero flux).

Example: Crack Domain

node = [1,0; 0,1; -1,0; 0,-1; 0,0; 1,0];        % nodes
elem = [5,1,2; 5,2,3; 5,3,4; 5,4,6];            % elements
elem = label(node,elem);                        % label the mesh
figure;
showmesh(node,elem);                            % plot mesh
findelem(node,elem);                            % plot element indices
findnode(node,2:6);                             % plot node indices
text(node(6,1),node(6,2)+0.075,int2str(1),'FontSize',16,'FontWeight','bold');
hold on;
plot([node(1,1),node(5,1)], [node(1,2),node(5,2)],'r-', 'LineWidth',3);
bdFlag = setboundary(node,elem,'Dirichlet');                % Dirichlet boundary condition
display(elem)
display(bdFlag)

elem =

     5     1     2
     5     2     3
     5     3     4
     5     4     6


bdFlag =

    1    0    1
    1    0    0
    1    0    0
    1    1    0

node 1 and node 6 are the same point (1,0)

bdFlag = setboundary(node,elem,'Dirichlet','abs(x) + abs(y) == 1','Neumann','y==0');
display(bdFlag)

bdFlag =

    1    0    2
    1    0    0
    1    0    0
    1    2    0

Example: Prism Domain

node = [-1,-1,-1; 1,-1,-1; 1,1,-1; -1,1,-1; -1,-1,1; 1,-1,1; 1,1,1; -1,1,1];
elem = [1,2,3,7; 1,6,2,7; 1,5,6,7];
elem = label3(node,elem);
figure;
showmesh3(node,elem);
view([-53,8]);
findnode3(node,[1 2 3 5 6 7]);
findelem3(node,elem);
bdFlag = setboundary3(node,elem,'Dirichlet','(z==1) | (z==-1)');
display(elem)
display(bdFlag)

elem =

     1     7     2     3
     1     7     6     2
     1     7     5     6


bdFlag =

    0    1    0    0
    0    0    0    0
    1    0    0    0

The top and bottom of the prism is set as Dirichlet boundary condition and other faces are zero flux boundary condition.

Remark