Adding example edge and midpoint matrices

get_edge_connectivity.m

@@ -0,0 +1,62 @@
+function [t_ed] = get_edge_connectivity(gd,sf,ns,mesh)
+% GET_EDGE_CONNECTIVITY  build edge connectivitity matrix for given mesh level
+%
+% Syntax:
+%     [t_ed] = get_edges(gd,sf,ns,mesh)
+%
+% Inputs:
+%     gd,sf,ns - outputs of pdepoly specifying domain
+%     mesh - max mesh level
+%
+% Outputs:
+%    t_ed - a 3xNumTriangles matrix representing the which edges
+%         correspond to which triangles. t_ed(i,T) represents the ith edge
+%         in triangle T.
+%
+% Usage Example:
+%    mesh = 8;
+%    pdepoly([0,1,1,0], [0,0,1,1]);
+%       (OR) [gd,sf,ns] = get_gd_sf_ns([0,1,1,0],[0,0,1,1]);
+%    [t_ed] = get_edge_connectivity(gd,sf,ns,mesh)
+%
+% Dependencies:
+%    find_edges.m
+%
+% Author: Nicole Stock
+% Date: Fall 2020
+
+addpath('data')
+addpath('edge_resources')
+
+model=createpde(1);
+g=decsg(gd,sf,ns);
+geometryFromEdges(model,g);
+[p,e,t]=initmesh(g,'hmax',inf);
+tr = triangulation(t(1:3,:)',p');
+
+if mesh > 1
+    % To ensure we refine every triangle the same
+    [~,num_node]=size(p);
+    it=zeros(1,num_node);
+    for i=1:num_node
+        it(i)=i;
+    end   
+
+    for i = 2:mesh
+        % Refine mesh to next level
+        [p,e,t]=refinemesh(g,p,e,t,it,'regular');
+        tr = triangulation(t(1:3,:)',p');
+    end
+
+    % get edge information for highest mesh level
+    ed = edges(tr);
+    t_ed = find_edge_connectivity(t,ed);
+end
+% get edge information for highest mesh level
+ed = edges(tr);
+t_ed = find_edge_connectivity(t,ed);
+% mesh level must be greater than 1
+
+% end main
+end
+

get_midpoints.m

@@ -0,0 +1,55 @@
+function [p2,t2] = get_midpoints(gd,sf,ns,mesh)
+%GET_MIDPOINTS
+%
+% Syntax:
+%     [p2,t2] = get_midpoints(gd,sf,ns,mesh)
+%
+% Inputs:
+%     gd,sf,ns - outputs of pdepoly specifying domain
+%     mesh - max mesh level
+%
+% Outputs:
+%     p2 - a 2xNumNodes matrix representing midpoint nodal coordinates.
+%     t2 - a 3xNumTriangles matrix representing the element connectivity in 
+%         terms of node IDs. The three node IDs in a column are the three
+%         midpoints of the node IDS in corresponding column in t.
+%
+% Usage Exampled:
+%    mesh = 8;
+%    pdepoly([0,1,1,0], [0,0,1,1]);
+%       (OR) [gd,sf,ns] = get_gd_sf_ns([0,1,1,0],[0,0,1,1]);
+%    [p2,t2] = get_midpoints(gd,sf,ns,mesh)
+%
+% Dependencies:
+%    find_midpoints.m
+%
+% Author: Nicole Stock
+% Date: Fall 2020
+
+model=createpde(1);
+g=decsg(gd,sf,ns);
+geometryFromEdges(model,g);
+[p,e,t]=initmesh(g,'hmax',inf);
+%pdemesh(p,e,t, 'NodeLabels','on', 'ElementLabels','on');
+
+if mesh > 1
+    % To ensure we refine every triangle the same
+    [~,num_node]=size(p);
+    it=zeros(1,num_node);
+    for i=1:num_node
+        it(i)=i;
+    end
+
+    for i = 2:mesh
+        % Refine mesh to next level
+        [p,e,t]=refinemesh(g,p,e,t,it,'regular');
+        %pdemesh(p,e,t, 'NodeLabels','on', 'ElementLabels','on');
+    end
+
+    % Find the midpoints for P2 nodal points for the highest mesh level
+    [p2,t2] = find_midpoints(p,t);
+end
+% mesh level must be greater than 1
+
+% end main
+end