Changeset 6705

Timestamp:

12/09/10 11:13:45 (14 years ago)

Author:

seroussi

Message:

include macayeal/stokes in setelementstype

File:

• issm/trunk/src/m/model/setelementstype.m (modified) (4 diffs)

issm/trunk/src/m/model/setelementstype.m

@@ -67 +67 @@
 
 %Stokes can only be used alone for now:
-if any(stokesflag) &any(hutterflag+macayealflag),
+if any(stokesflag) &any(hutterflag),
         error('setelementstype error message: stokes cannot be used with any other model for now, put stokes everywhere')
 end
@@ -123 +123 @@
 gridonmacayealpattyn=zeros(md.numberofgrids,1);
 gridonpattynstokes=zeros(md.numberofgrids,1);
+gridonmacayealstokes=zeros(md.numberofgrids,1);
 if strcmpi(coupling_method,'penalties'),
         %Create the border grids between Pattyn and MacAyeal and extrude them
@@ -177 +178 @@
                 gridonpattynstokes=zeros(md.numberofgrids,1);
                 gridonpattynstokes(md.elements(find(pattynstokesflag),:))=1;
-        elseif any(stokesflag) & any(macayealflag+hutterflag),
+        elseif any(stokesflag) & any(macayealflag),
+                %Find grid at the border
+                gridonmacayealstokes(find(gridonmacayeal & gridonstokes))=1;
+                %Stokes elements in contact with this layer become MacAyealStokes elements
+                matrixelements=ismember(md.elements,find(gridonmacayealstokes));
+                commonelements=sum(matrixelements,2)~=0;
+                commonelements(find(macayealflag))=0; %only one layer: the elements previously in macayeal
+                stokesflag(find(commonelements))=0; %these elements are now macayealmacayealelements
+                macayealstokesflag=zeros(md.numberofelements,1);
+                macayealstokesflag(find(commonelements))=1;
+                gridonstokes=zeros(md.numberofgrids,1);
+                gridonstokes(md.elements(find(stokesflag),:))=1;
+                md.gridonstokes=gridonstokes;
+
+                %Create MacaAyealMacAyealApproximation where needed
+                md.elements_type(find(macayealstokesflag))=MacAyealStokesApproximationEnum();
+
+                %Now recreate gridonmacayealstokes
+                gridonmacayealstokes=zeros(md.numberofgrids,1);
+                gridonmacayealstokes(md.elements(find(macayealstokesflag),:))=1;
+        elseif any(stokesflag) & any(hutterflag),
                 error('type of coupling not supported yet');
         end
@@ -206 +227 @@
 pos=find(gridonpattynstokes);
 md.vertices_type(pos)=PattynStokesApproximationEnum();
+pos=find(gridonmacayealstokes);
+md.vertices_type(pos)=MacAyealStokesApproximationEnum();
 
 %figure out solution types