Changeset 23189 for issm/trunk/src/m/mesh/augment2dmesh.m

Timestamp:

08/28/18 09:45:51 (7 years ago)

Author:

Mathieu Morlighem

Message:

merged trunk-jpl and trunk for revision 23187

Location:

issm/trunk

Files:

• . (modified) (1 prop)
• src (modified) (1 prop)
• src/m/mesh/augment2dmesh.m (modified) (6 diffs)

issm/trunk/src/m/mesh/augment2dmesh.m

@@ -1 @@
-function mh=augment2dmesh(mh,mhband,varargin)
+function [mh,internalsegments,internalelements]=augment2dmesh(mh,mhband,varargin)
 %AUGMENT2DMESH - augment mh mesh with a band around it (provided by mhband)
 %
 %   Usage:
 %      mh=augment2dmesh(mh,mhband);
+%      [mh,internalsegments]=augment2dmesh(mh,mhband);
 %
 %   Example: 
@@ -17 +18 @@
 mhband.segments(:,3)=mhband.segments(:,3)+mh.numberofelements;
 
+bandsegments=mhband.segments; %keep  a copy around for when user requests this be returned. 
+
 %The innner segments of mhband and the outer segments of mh are identical. Go into  the elements of 
 %mhband and set them to their md1 equivalent: 
-tol=1; %1 meter 
+tol=getfieldvalue(options,'tolerance',1); %1 meter 
 for i=1:length(mhband.segments),
         node2=mhband.segments(i,1);
@@ -29 +32 @@
                         %go into the mesh of mhband, and replace by node1.
                         pos=find(mhband.elements==node2); mhband.elements(pos)=node1;
-                        segs=mhband.segments(:,1:2); pos=find(segs==node2); segs(pos)=node1; mhband.segments(:,1:2)=segs;
+                        segs=mhband.segments(:,1:2); pos=find(segs==node2); segs(pos)=node1; mhband.segments(:,1:2)=segs; bandsegments(:,1:2)=segs;
+                        %flag the segment for future erasal: 
+                        pos=find(mhband.segments(:,1)==node1 | mhband.segments(:,2)==node1); mhband.segments(pos,3)=NaN; 
                         break;
                 end
@@ -36 +41 @@
 
 %Do the merge: 
+internalelements=mh.elements;
 mh.elements=[mh.elements;mhband.elements];
 mh.x=[mh.x;mhband.x];
 mh.y=[mh.y;mhband.y];
-mh.lat=[mh.lat;mhband.lat];
-mh.long=[mh.long;mhband.long];
-mh.segments=[mh.segments;mhband.segments];
+if ~isnan(mh.lat), mh.lat=[mh.lat;mhband.lat]; end
+if ~isnan(mh.long), mh.long=[mh.long;mhband.long];  end
+mh.segments=[mhband.segments]; 
+
+%segments that are internal mght have been requested: 
+if strcmpi(getfieldvalue(options,'internalsegments','off'),'on'),
+        pos=find(isnan(mh.segments(:,3)));
+        internalsegments=bandsegments(pos,:);
+end
+
+%get rid of segments that are internal: 
+pos=find(~isnan(mh.segments(:,3)));
+mh.segments=mh.segments(pos,:); 
 
 %Remove orphans:
@@ -53 +69 @@
         x=[x(1:orphan(i)-(i-1)-1); x(orphan(i)-(i-1)+1:end)];
         y=[y(1:orphan(i)-(i-1)-1); y(orphan(i)-(i-1)+1:end)];
-        lat=[lat(1:orphan(i)-(i-1)-1); lat(orphan(i)-(i-1)+1:end)];
-        long=[long(1:orphan(i)-(i-1)-1); long(orphan(i)-(i-1)+1:end)];
+        if ~isnan(lat), lat=[lat(1:orphan(i)-(i-1)-1); lat(orphan(i)-(i-1)+1:end)]; end
+                if ~isnan(long), long=[long(1:orphan(i)-(i-1)-1); long(orphan(i)-(i-1)+1:end)]; end
         %update elements
         pos=find(elements>orphan(i)-(i-1));
@@ -73 +89 @@
 mh.numberofelements=length(mh.elements);
 mh.numberofvertices=length(mh.x);
+
+%vertex on boundaries: 
+mh.vertexonboundary=zeros(mh.numberofvertices,1);
+mh.vertexonboundary(mh.segments(:,1:2))=1;
+
+%connectivities: 
+mh.vertexconnectivity=NodeConnectivity(mh.elements,mh.numberofvertices);
+mh.elementconnectivity=ElementConnectivity(mh.elements,mh.vertexconnectivity);
+
+%return: 
+if nargout==1,
+        varargout{1}=mh;
+elseif nargout==2,
+        varargout{1}=mh;
+        varargout{2}=internalsegments;
+elseif nargout==3,
+        varargout{1}=mh;
+        varargout{2}=internalsegments;
+        varargout{3}=internalelements;
+else 
+        error('unknow number of output arguments')
+end