Changeset 5442 for issm/trunk/src/m/classes/public/modelextract.m

Timestamp:

08/20/10 09:05:55 (15 years ago)

Author:

Mathieu Morlighem

Message:

Better checks in ismodelselfconsistent (problem with matrices) and fixed edges extraction

File:

• issm/trunk/src/m/classes/public/modelextract.m (modified) (1 diff)

issm/trunk/src/m/classes/public/modelextract.m

@@ -148 +148 @@
         if size(md2.edges,2)>1, %do not use ~isnan because there are some NaNs...
                 %renumber first two columns
-                md2.edges(:,1)=Pgrid(md2.edges(:,1)); 
-                md2.edges(:,2)=Pgrid(md2.edges(:,2)); 
-                md2.edges(:,3)=Pelem(md2.edges(:,3));
-                %remove column where the first element is not in the domain
-                md2.edges=md2.edges(find(md2.edges(:,1) & md2.edges(:,2) & md2.edges(:,3)),:);
-                %now we must renumber the second column and put NaN if the element is not in the domain
-                A=md2.edges(:,4);
-                A(find(~isnan(A)))=Pelem(A(find(~isnan(A))));
-                A(find(A==0))=NaN;
-                md2.edges(:,4)=A; clear A;
-        end
-
+                pos=find(~isnan(md2.edges(:,4)));
+                md2.edges(:  ,1)=Pgrid(md2.edges(:,1)); 
+                md2.edges(:  ,2)=Pgrid(md2.edges(:,2)); 
+                md2.edges(:  ,3)=Pelem(md2.edges(:,3));
+                md2.edges(pos,4)=Pelem(md2.edges(pos,4));
+                %remove edges when the 2 vertices are not in the domain.
+                md2.edges=md2.edges(find(md2.edges(:,1) & md2.edges(:,2)),:);
+                %Replace all zeros by NaN in the last two columns;
+                pos=find(md2.edges(:,3)==0);
+                md2.edges(pos,3)=NaN;
+                pos=find(md2.edges(:,4)==0);
+                md2.edges(pos,4)=NaN;
+                %Invert NaN of the third column with last column
+                pos=find(isnan(md2.edges(:,3)));
+                md2.edges(pos,3)=md2.edges(pos,4);
+                md2.edges(pos,4)=NaN;
+                %Finally remove edges that do not belong to any element
+                pos=find(isnan(md2.edges(:,3)) & isnan(md2.edges(:,4)));
+                md2.edges(pos,:)=[];
+        end
 
         %Penalties