Changeset 16295

Timestamp:

10/03/13 15:35:54 (11 years ago)

Author:

Mathieu Morlighem

Message:

CHG: preparing 2d horizontal mesh

Location:

issm/trunk-jpl/src/m/classes

Files:

• flowequation.m (modified) (1 diff)
• flowequation.py (modified) (1 diff)
• mesh2d.m (modified) (4 diffs)
• mesh2d.py (modified) (2 diffs)
• stressbalance.m (modified) (2 diffs)
• stressbalance.py (modified) (2 diffs)

issm/trunk-jpl/src/m/classes/flowequation.m

@@ -91 +91 @@
                                 md = checkfield(md,'flowequation.borderHO','size',[md.mesh.numberofvertices 1],'values',[0 1]);
                                 md = checkfield(md,'flowequation.borderFS','size',[md.mesh.numberofvertices 1],'values',[0 1]);
-                                if (md.mesh.dimension==2),
+                                if strcmp(meshtype(md.mesh),'2Dhorizontal')
                                         md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices 1],'values',[1:2]);
                                         md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements 1],'values',[1:2]);
-                                else
+                                elseif strcmp(meshtype(md.mesh),'3D'),
                                         md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices 1],'values',[0:8]);
                                         md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements 1],'values',[0:8]);
+                                else
+                                        error('Case not supported yet');
                                 end
                                 if ~(obj.isSIA || obj.isSSA || obj.isL1L2 || obj.isHO || obj.isFS),

issm/trunk-jpl/src/m/classes/flowequation.py

@@ -79 +79 @@
                         md = checkfield(md,'flowequation.borderHO','size',[md.mesh.numberofvertices],'values',[0,1])
                         md = checkfield(md,'flowequation.borderFS','size',[md.mesh.numberofvertices],'values',[0,1])
-                        if md.mesh.dimension==2:
+                        if strcmp(md.mesh.meshtype(),'2Dhorizontal'):
                                 md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',[1,2])
                                 md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements],'values',[1,2])
-                        else:
+                        elif strcmp(md.mesh.meshtype(),'3D'):
                                 md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',numpy.arange(0,8+1))
                                 md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements],'values',numpy.arange(0,8+1))
+                        else:
+                                raise RuntimeError('mesh type not supported yet')
                         if not (self.isSIA or self.isSSA or self.isL1L2 or self.isHO or self.isFS):
                                 md.checkmessage("no element types set for this model")

issm/trunk-jpl/src/m/classes/mesh2d.m

@@ -1 @@
-%MESH class definition
+%MESH2D class definition
 %
 %   Usage:
-%      mesh=mesh();
+%      mesh2d=mesh2d();
 
-classdef mesh
+classdef mesh2d
         properties (SetAccess=public) 
                 x                           = NaN;
@@ -30 +30 @@
         end
         methods
-                function obj = mesh(varargin) % {{{
+                function obj = mesh2d(varargin) % {{{
                         switch nargin
                                 case 0
@@ -58 +58 @@
                         md = checkfield(md,'mesh.numberofelements','>',0);
                         md = checkfield(md,'mesh.numberofvertices','>',0);
-                        md = checkfield(md,'mesh.vertexonbed','size',[md.mesh.numberofvertices 1],'values',[0 1]);
                         md = checkfield(md,'mesh.average_vertex_connectivity','>=',9,'message','''mesh.average_vertex_connectivity'' should be at least 9 in 2d');
 
@@ -96 +95 @@
                 end % }}}
                 function marshall(obj,md,fid) % {{{
-                        WriteData(fid,'object',obj,'fieldname','x','format','DoubleMat','mattype',1);
-                        WriteData(fid,'object',obj,'fieldname','y','format','DoubleMat','mattype',1);
-                        WriteData(fid,'object',obj,'fieldname','elements','format','DoubleMat','mattype',2);
-                        WriteData(fid,'object',obj,'fieldname','numberofelements','format','Integer');
-                        WriteData(fid,'object',obj,'fieldname','numberofvertices','format','Integer');
-                        WriteData(fid,'object',obj,'fieldname','average_vertex_connectivity','format','Integer');
+                        WriteData(fid,'enum',MeshTypeEnum(),'data',StringToEnum(['Mesh' meshtype(obj)]),'format','Integer');
+                        WriteData(fid,'object',obj,'class','mesh','fieldname','x','format','DoubleMat','mattype',1);
+                        WriteData(fid,'object',obj,'class','mesh','fieldname','y','format','DoubleMat','mattype',1);
+                        WriteData(fid,'enum',MeshZEnum(),'data',zeros(obj.numberofvertices,1),'format','DoubleMat','mattype',1);
+                        WriteData(fid,'object',obj,'class','mesh','fieldname','elements','format','DoubleMat','mattype',2);
+                        WriteData(fid,'object',obj,'class','mesh','fieldname','numberofelements','format','Integer');
+                        WriteData(fid,'object',obj,'class','mesh','fieldname','numberofvertices','format','Integer');
+                        WriteData(fid,'object',obj,'class','mesh','fieldname','average_vertex_connectivity','format','Integer');
+                end % }}}
+                function t = meshtype(obj) % {{{
+                        t = '2Dhorizontal';
                 end % }}}
         end

issm/trunk-jpl/src/m/classes/mesh2d.py

@@ -5 +5 @@
 from MatlabFuncs import *
 
-class mesh(object):
+class mesh2d(object):
         """
-        MESH class definition
+        MESH2D class definition
 
            Usage:
-              mesh=mesh();
+              mesh2d=mesh2d();
         """
 
@@ -98 +98 @@
                 return md
         # }}}
+        def meshtype(self): # {{{
+                return "2Dhorizontal"
+        #}}}
         def marshall(self,md,fid):    # {{{
-                WriteData(fid,'object',self,'fieldname','x','format','DoubleMat','mattype',1)
-                WriteData(fid,'object',self,'fieldname','y','format','DoubleMat','mattype',1)
-                WriteData(fid,'object',self,'fieldname','elements','format','DoubleMat','mattype',2)
-                WriteData(fid,'object',self,'fieldname','numberofelements','format','Integer')
-                WriteData(fid,'object',self,'fieldname','numberofvertices','format','Integer')
-                WriteData(fid,'object',self,'fieldname','average_vertex_connectivity','format','Integer')
+                WriteData(fid,'enum',MeshTypeEnum(),'data',StringToEnum("Mesh"+self.meshtype())[0],'format','Integer');
+                WriteData(fid,'object',self,'class','mesh','fieldname','x','format','DoubleMat','mattype',1)
+                WriteData(fid,'object',self,'class','mesh','fieldname','y','format','DoubleMat','mattype',1)
+                WriteData(fid,'enum',MeshZEnum(),'data',numpy.zeros(self.numberofvertices),'format','DoubleMat','mattype',1);
+                WriteData(fid,'object',self,'class','mesh','fieldname','elements','format','DoubleMat','mattype',2)
+                WriteData(fid,'object',self,'class','mesh','fieldname','numberofelements','format','Integer')
+                WriteData(fid,'object',self,'class','mesh','fieldname','numberofvertices','format','Integer')
+                WriteData(fid,'object',self,'class','mesh','fieldname','average_vertex_connectivity','format','Integer')
         # }}}

issm/trunk-jpl/src/m/classes/stressbalance.m

@@ -81 +81 @@
                         md = checkfield(md,'stressbalance.spcvx','forcing',1);
                         md = checkfield(md,'stressbalance.spcvy','forcing',1);
-                        if md.mesh.dimension==3, md = checkfield(md,'stressbalance.spcvz','forcing',1); end
+                        if strcmp(meshtype(md.mesh),'3D'), md = checkfield(md,'stressbalance.spcvz','forcing',1); end
                         md = checkfield(md,'stressbalance.restol','size',[1 1],'>',0,'NaN',1);
                         md = checkfield(md,'stressbalance.reltol','size',[1 1]);
@@ -111 +111 @@
                         end
                         %CHECK THAT NO rotation specified for FS Grounded ice at base
-                        if md.mesh.dimension==3 & md.flowequation.isFS,
+                        if strcmp(meshtype(md.mesh),'3D') & md.flowequation.isFS,
                                 pos=find(md.mask.groundedice_levelset>0. & md.mesh.vertexonbed);
                                 if any(~isnan(md.stressbalance.referential(pos,:))),

issm/trunk-jpl/src/m/classes/stressbalance.py

@@ -112 +112 @@
                 md = checkfield(md,'stressbalance.spcvx','forcing',1)
                 md = checkfield(md,'stressbalance.spcvy','forcing',1)
-                if md.mesh.dimension==3:
+                if strcmp(md.mesh.meshtype(),'3D'):
                         md = checkfield(md,'stressbalance.spcvz','forcing',1)
                 md = checkfield(md,'stressbalance.restol','size',[1],'>',0)
@@ -145 +145 @@
                                         md.checkmessage("Vectors in stressbalance.referential (columns 1 to 3 and 4 to 6) must be orthogonal")
                 #CHECK THAT NO rotation specified for FS Grounded ice at base
-#               if md.mesh.dimension==3 & md.flowequation.isFS,
-                if md.mesh.dimension==3 and md.flowequation.isFS:
+                if strcmp(md.mesh.meshtype(),'3D') and md.flowequation.isFS:
                         pos=numpy.nonzero(numpy.logical_and(md.mask.groundedice_levelset,md.mesh.vertexonbed))
                         if numpy.any(numpy.logical_not(numpy.isnan(md.stressbalance.referential[pos,:]))):