Changeset 15796

Timestamp:

08/12/13 14:18:56 (12 years ago)

Author:

seroussi

Message:

CHG: trying to add boundary condition change in python

Location:

issm/trunk-jpl/src/m/boundaryconditions

Files:

• SetIceShelfBC.m (modified) (1 diff)
• SetMarineIceSheetBC.m (modified) (2 diffs)
• SetMarineIceSheetBC.py (modified) (2 diffs)

issm/trunk-jpl/src/m/boundaryconditions/SetIceShelfBC.m

@@ -46 +46 @@
 segments=find(sum(segmentsfront,2)~=numbernodesfront);
 %Find all nodes for these segments and spc them
-pos=[md.mesh.segments(segments,1:numbernodesfront)];
+pos=md.mesh.segments(segments,1:numbernodesfront);
 md.stressbalance.spcvx(pos(:))=0;
 md.stressbalance.spcvy(pos(:))=0;

issm/trunk-jpl/src/m/boundaryconditions/SetMarineIceSheetBC.m

@@ -41 +41 @@
 md.stressbalance.loadingforce=0*ones(md.mesh.numberofvertices,3);
 
-
 %Position of ice front
 md.mask.icelevelset(find(vertexonicefront))=0;
@@ -54 +53 @@
 segments=find(sum(segmentsfront,2)~=numbernodesfront);
 %Find all nodes for these segments and spc them
-pos=[md.mesh.segments(segments,1:numbernodesfront)];
+pos=md.mesh.segments(segments,1:numbernodesfront);
 md.stressbalance.spcvx(pos(:))=0;
 md.stressbalance.spcvy(pos(:))=0;

issm/trunk-jpl/src/m/boundaryconditions/SetMarineIceSheetBC.py

@@ -43 +43 @@
         md.stressbalance.spcvy=float('nan')*numpy.ones((md.mesh.numberofvertices,1))
         md.stressbalance.spcvz=float('nan')*numpy.ones((md.mesh.numberofvertices,1))
-        md.stressbalance.spcvx[pos]=0
-        md.stressbalance.spcvy[pos]=0
-        md.stressbalance.spcvz[pos]=0
         md.stressbalance.referential=float('nan')*numpy.ones((md.mesh.numberofvertices,6))
         md.stressbalance.loadingforce=0*numpy.ones((md.mesh.numberofvertices,3))
+
+        #Position of ice front
+        pos=numpy.nonzero(vertexonicefront)[0]
+        md.mask.icelevelset[pos]=0
+
+        #First find segments that are not completely on the front
+        if md.mesh.dimension==2:
+                numbernodesfront=2
+        else:
+                numbernodesfront=4
+        segmentsfront=numpy.nonzeros(md.mask.icelevelset(md.mesh.segments(:,0:-1)))
+        segments=numpy.sum(segmentsfront,axis=1)!=numbernodesfront).astype(int) 
+        #Find all nodes for these segments and spc them
+        pos=md.mesh.segments[segments,0:-1]
+        md.stressbalance.spcvx[pos[:]]=0
+        md.stressbalance.spcvy[pos[:]]=0
+        md.stressbalance.spcvz[pos[:]]=0
 
         #Dirichlet Values
@@ -60 +74 @@
         pos=numpy.nonzero(md.mesh.vertexonboundary)[0]
         md.hydrology.spcwatercolumn[pos,0]=1
-
-        #Position of ice front
-        pos=numpy.nonzero(vertexonicefront)[0]
-        md.mask.icelevelset[pos]=0
 
         #Create zeros basalforcings and surfaceforcings