source: issm/oecreview/Archive/25834-26739/ISSM-26554-26555.diff@ 26740

Last change on this file since 26740 was 26740, checked in by Mathieu Morlighem, 3 years ago

CHG: added 25834-26739
File size: 41.5 KB

  • ../trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

     
    2020//#define FSANALYTICAL 10
    2121//#define LATERALFRICTION 1
    2222//#define DISCSLOPE 1 //testing for SSA
    23 #define NOSPCSHEARVEL 1 //MLHO
    2423
    2524/*Model processing*/
    2625void StressbalanceAnalysis::CreateConstraints(Constraints* constraints,IoModel* iomodel){/*{{{*/
    … …  
    195194                                }
    196195                        }
    197196                        else{//MLHO
    198                                 #ifdef NOSPCSHEARVEL   
    199                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx",StressbalanceAnalysisEnum,finiteelement,0);
    200                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy",StressbalanceAnalysisEnum,finiteelement,2);
    201                                 #else
    202                                 /*Default: apply spcs to shear vx and shear vy*/
    203                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx",StressbalanceAnalysisEnum,finiteelement,0);
    204                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx",StressbalanceAnalysisEnum,finiteelement,1);
    205                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy",StressbalanceAnalysisEnum,finiteelement,2);
    206                                 IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy",StressbalanceAnalysisEnum,finiteelement,3);
    207                                 #endif
     197                                IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx_base",StressbalanceAnalysisEnum,finiteelement,0);
     198                                IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvx_shear",StressbalanceAnalysisEnum,finiteelement,1);
     199                                IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy_base",StressbalanceAnalysisEnum,finiteelement,2);
     200                                IoModelToConstraintsx(constraints,iomodel,"md.stressbalance.spcvy_shear",StressbalanceAnalysisEnum,finiteelement,3);
    208201                        }
    209202                }
    210203

  • ../trunk-jpl/src/m/boundaryconditions/SetMLHOBC.m

     
     1function md=SetMLHOBC(md)
     2%SETMLHOBC - Create the boundary conditions for stressbalance for MLHO: VxBase, VyBase, VxShear, VyShear
     3%
     4%   Usage:
     5%      md=SetMLHOBC(md)
     6%
     7
     8
     9%node on Dirichlet
     10if md.flowequation.isMLHO
     11        md.stressbalance.spcvx_base=md.stressbalance.spcvx;
     12        md.stressbalance.spcvy_base=md.stressbalance.spcvy;
     13
     14        md.stressbalance.spcvx_shear=NaN*ones(size(md.stressbalance.spcvx_base));
     15        md.stressbalance.spcvy_shear=NaN*ones(size(md.stressbalance.spcvy_base));
     16end

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

     
    6161                                if isfield(objstruct,'borderstokes'),  self.borderFS  = objstruct.borderstokes;   end;
    6262                        end
    6363
    64                         %Oct 8 2020
    65                         %Warning if MLHO
     64                        %Nov 6 2021
    6665                        if any(self.vertex_equation==4)
    67                                 warning(['Monolayer Higher-Order (MLHO) detected in md.flowequation, this is probably a mistake, you need to run setflowequation again']);
     66                                disp(['Monolayer Higher-Order (MLHO) detected in md.flowequation, this is still under development. Please double check your settings.']);
    6867                        end
    6968
    7069                end% }}}

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

     
    77        properties (SetAccess=public)
    88                spcvx                  = NaN;
    99                spcvy                  = NaN;
     10                spcvx_base             = NaN;
     11                spcvy_base             = NaN;
     12                spcvx_shear            = NaN;
     13                spcvy_shear            = NaN;
    1014                spcvz                  = NaN;
    1115                restol                 = 0;
    1216                reltol                 = 0;
    … …  
    3135                        self.referential=project3d(md,'vector',self.referential,'type','node');
    3236                        self.loadingforce=project3d(md,'vector',self.loadingforce,'type','node');
    3337
     38                        % for MLHO
     39                        if md.flowequation.isMLHO
     40                                self.spcvx_base=project3d(md,'vector',self.spcvx_base,'type','node');
     41                                self.spcvy_base=project3d(md,'vector',self.spcvy_base,'type','node');
     42                                self.spcvx_shear=project3d(md,'vector',self.spcvx_shear,'type','poly','degree',4);
     43                                self.spcvy_shear=project3d(md,'vector',self.spcvy_shear,'type','poly','degree',4);
     44                        end
    3445                end % }}}
    3546                function self = stressbalance(varargin) % {{{
    3647                        switch nargin
    … …  
    122133                                end
    123134                                md = checkfield(md,'fieldname','stressbalance.FSreconditioning','>',0);
    124135                        end
     136                        % CHECK THIS ONLY WORKS FOR MLHO
     137                        if md.flowequation.isMLHO
     138                                md = checkfield(md,'fieldname','stressbalance.spcvx_base','Inf',1,'timeseries',1);
     139                                md = checkfield(md,'fieldname','stressbalance.spcvy_base','Inf',1,'timeseries',1);
     140                                md = checkfield(md,'fieldname','stressbalance.spcvx_shear','Inf',1,'timeseries',1);
     141                                md = checkfield(md,'fieldname','stressbalance.spcvy_shear','Inf',1,'timeseries',1);
     142                        end
    125143                end % }}}
    126144                function list=defaultoutputs(self,md) % {{{
    127145
    … …  
    150168                        fielddisplay(self,'spcvy','y-axis velocity constraint (NaN means no constraint) [m/yr]');
    151169                        fielddisplay(self,'spcvz','z-axis velocity constraint (NaN means no constraint) [m/yr]');
    152170
     171                        disp(sprintf('\n      %s','MLHO boundary conditions:'));
     172                        fielddisplay(self,'spcvx_base','x-axis basal velocity constraint (NaN means no constraint) [m/yr]');
     173                        fielddisplay(self,'spcvy_base','y-axis basal velocity constraint (NaN means no constraint) [m/yr]');
     174                        fielddisplay(self,'spcvx_shear','x-axis shear velocity constraint (NaN means no constraint) [m/yr]');
     175                        fielddisplay(self,'spcvy_shear','y-axis shear velocity constraint (NaN means no constraint) [m/yr]');
     176
    153177                        disp(sprintf('\n      %s','Rift options:'));
    154178                        fielddisplay(self,'rift_penalty_threshold','threshold for instability of mechanical constraints');
    155179                        fielddisplay(self,'rift_penalty_lock','number of iterations before rift penalties are locked');
    … …  
    201225                                outputs      = [outputs defaultoutputs(self,md)]; %add defaults
    202226                        end
    203227                        WriteData(fid,prefix,'data',outputs,'name','md.stressbalance.requested_outputs','format','StringArray');
     228                        % for MLHO
     229                        if (md.flowequation.isMLHO)
     230                                WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvx_base','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
     231                                WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvy_base','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
     232                                WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvx_shear','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
     233                                WriteData(fid,prefix,'object',self,'class','stressbalance','fieldname','spcvy_shear','format','DoubleMat','mattype',1,'scale',1./yts,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
     234                        end
    204235                end % }}}
    205236                function savemodeljs(self,fid,modelname) % {{{
    206237
    … …  
    222253                        writejs2Darray(fid,[modelname '.stressbalance.loadingforce'],self.loadingforce);
    223254                        writejscellstring(fid,[modelname '.stressbalance.requested_outputs'],self.requested_outputs);
    224255
     256                        writejs1Darray(fid,[modelname '.stressbalance.spcvx_base'],self.spcvx_shear);
     257                        writejs1Darray(fid,[modelname '.stressbalance.spcvy_base'],self.spcvy_shear);
     258                        writejs1Darray(fid,[modelname '.stressbalance.spcvx_shear'],self.spcvx_shear);
     259                        writejs1Darray(fid,[modelname '.stressbalance.spcvy_shear'],self.spcvy_shear);
    225260                end % }}}
    226261        end
    227262end

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

     
    1919    def __init__(self):  # {{{
    2020        self.spcvx = np.nan
    2121        self.spcvy = np.nan
     22        self.spcvx_base = np.nan
     23        self.spcvy_base = np.nan
     24        self.spcvx_shear = np.nan
     25        self.spcvy_shear = np.nan
    2226        self.spcvz = np.nan
    2327        self.restol = 0
    2428        self.reltol = 0
    … …  
    5458        s += '{}\n'.format(fielddisplay(self, 'spcvy', 'y-axis velocity constraint (NaN means no constraint) [m / yr]'))
    5559        s += '{}\n'.format(fielddisplay(self, 'spcvz', 'z-axis velocity constraint (NaN means no constraint) [m / yr]'))
    5660        s += '{}\n'.format(fielddisplay(self, 'icefront', 'segments on ice front list (last column 0: Air, 1: Water, 2: Ice'))
     61        s += '      MLHO boundary conditions:\n'
     62        s += '{}\n'.format(fielddisplay(self, 'spcvx_base', 'x-axis basal velocity constraint (NaN means no constraint) [m / yr]'))
     63        s += '{}\n'.format(fielddisplay(self, 'spcvy_base', 'y-axis basal velocity constraint (NaN means no constraint) [m / yr]'))
     64        s += '{}\n'.format(fielddisplay(self, 'spcvx_shear', 'x-axis shear velocity constraint (NaN means no constraint) [m / yr]'))
     65        s += '{}\n'.format(fielddisplay(self, 'spcvy_shear', 'y-axis shear velocity constraint (NaN means no constraint) [m / yr]'))
    5766        s += '      Rift options:\n'
    5867        s += '{}\n'.format(fielddisplay(self, 'rift_penalty_threshold', 'threshold for instability of mechanical constraints'))
    5968        s += '{}\n'.format(fielddisplay(self, 'rift_penalty_lock', 'number of iterations before rift penalties are locked'))
    … …  
    7685        self.referential = project3d(md, 'vector', self.referential, 'type', 'node')
    7786        self.loadingforce = project3d(md, 'vector', self.loadingforce, 'type', 'node')
    7887
     88        if md.flowequation.isMLHO:
     89            self.spcvx_base = project3d(md, 'vector', self.spcvx_base, 'type', 'node')
     90            self.spcvy_base = project3d(md, 'vector', self.spcvy_base, 'type', 'node')
     91            self.spcvx_shear = project3d(md, 'vector', self.spcvx_shear, 'type', 'poly', 'degree', 4)
     92            self.spcvy_shear = project3d(md, 'vector', self.spcvy_shear, 'type', 'poly', 'degree', 4)
     93
    7994        return self
    8095    # }}}
    8196
    … …  
    160175            pos = np.nonzero(np.logical_and(md.mask.ocean_levelset, md.mesh.vertexonbase))
    161176            if np.any(np.logical_not(np.isnan(md.stressbalance.referential[pos, :]))):
    162177                md.checkmessage("no referential should be specified for basal vertices of grounded ice")
     178        if md.flowequation.isMLHO:
     179            md = checkfield(md, 'fieldname', 'stressbalance.spcvx_base', 'Inf', 1, 'timeseries', 1)
     180            md = checkfield(md, 'fieldname', 'stressbalance.spcvy_base', 'Inf', 1, 'timeseries', 1)
     181            md = checkfield(md, 'fieldname', 'stressbalance.spcvx_shear', 'Inf', 1, 'timeseries', 1)
     182            md = checkfield(md, 'fieldname', 'stressbalance.spcvy_shear', 'Inf', 1, 'timeseries', 1)
    163183        return md
    164184    # }}}
    165185
    … …  
    194214            outputscopy = outputs[0:max(0, indices[0] - 1)] + self.defaultoutputs(md) + outputs[indices[0] + 1:]
    195215            outputs = outputscopy
    196216        WriteData(fid, prefix, 'data', outputs, 'name', 'md.stressbalance.requested_outputs', 'format', 'StringArray')
     217        # MLHO
     218        if md.flowequation.isMLHO:
     219            WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvx_base', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
     220            WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvy_base', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
     221            WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvx_shear', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
     222            WriteData(fid, prefix, 'object', self, 'class', 'stressbalance', 'fieldname', 'spcvy_shear', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1. / yts, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', yts)
    197223    # }}}

  • ../trunk-jpl/src/m/extrusion/project3d.m

     
    66%   element vector of size (md.mesh.numberofelements2d,N/A).
    77%   arguments:
    88%      'vector': 2d vector
    9 %      'type': 'element' or 'node'.
     9%      'type': 'element' or 'node' or 'poly'
    1010%   options:
    1111%      'layer' a layer number where vector should keep its values. If not specified, all layers adopt the
    1212%             value of the 2d vector.
    1313%      'padding': default to 0 (value adopted by other 3d layers not being projected0
     14%                'degree': degree of polynomials when extrude from bottom to the top
    1415%
    1516%   Egs:
    1617%      extruded_vector=project3d(md,'vector',vector2d,'type','node','layer',1,'padding',NaN);
    … …  
    3233type         = getfieldvalue(options,'type');       %mandatory
    3334layer        = getfieldvalue(options,'layer',0);    %optional (do all layers otherwise)
    3435paddingvalue = getfieldvalue(options,'padding',0);  %0 by default
     36polyexponent = getfieldvalue(options,'degree',0);   %0 by default, 0-degree polynomial
    3537
    3638if length(vector2d)==1,
    3739        projected_vector=vector2d;
    … …  
    7981                projected_vector( ((layer-1)*md.mesh.numberofelements2d+1):(layer*md.mesh.numberofelements2d),:)=vector2d;
    8082        end
    8183
     84elseif strcmpi(type,'poly'), % interpolate values from 0 to 1 with a polynomial degree n
     85        %Initialize 3d vector
     86        if size(vector2d,1)==md.mesh.numberofvertices2d
     87                projected_vector=paddingvalue*ones(md.mesh.numberofvertices, size(vector2d,2));
     88        elseif size(vector2d,1)==md.mesh.numberofvertices2d+1
     89                projected_vector=paddingvalue*ones(md.mesh.numberofvertices+1,size(vector2d,2));
     90                projected_vector(end,:)=vector2d(end,:);
     91                vector2d=vector2d(1:end-1,:);
     92        else
     93                error('vector length not supported')
     94        end
     95
     96        polycoeff = [0:1./(md.mesh.numberoflayers-1):1];
     97
     98        %Fill in
     99        if layer==0,
     100                for i=1:md.mesh.numberoflayers,
     101                        projected_vector(((i-1)*md.mesh.numberofvertices2d+1):(i*md.mesh.numberofvertices2d),:)=vector2d*(1-(1-polycoeff(i)).^polyexponent);
     102                end
     103        else
     104                projected_vector(((layer-1)*md.mesh.numberofvertices2d+1):(layer*md.mesh.numberofvertices2d),:)=vector2d*(1-(1-polycoeff(layer)).^polyexponent);
     105        end
     106
    82107else
    83108        error('project3d error message: unknown projection type');
    84109end

  • ../trunk-jpl/test/NightlyRun/test127.m

     
    1717        massfluxatgate('name','MassFlux5','profilename',['../Exp/MassFlux5.exp'],'definitionstring','Outputdefinition5'),...
    1818        massfluxatgate('name','MassFlux6','profilename',['../Exp/MassFlux6.exp'],'definitionstring','Outputdefinition6')...
    1919        };
    20 
     20md=SetMLHOBC(md);
    2121md=solve(md,'Stressbalance');
    2222
    2323%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test128.m

     
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.transient.requested_outputs={'IceVolume','VxShear','VyShear','VxBase','VyBase','VxSurface','VySurface'};
    88
     9md=SetMLHOBC(md);
    910md=solve(md,'Transient');
    1011
    1112%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test128.py

     
    1515md.cluster = generic('name', gethostname(), 'np', 3)
    1616md.transient.requested_outputs = ['IceVolume','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase']
    1717
     18md = SetMLHOBC(md);
    1819md = solve(md, 'Transient')
    1920
    2021#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test129.m

     
    1313md.timestepping.final_time=19;
    1414md.settings.output_frequency=2;
    1515
     16md=SetMLHOBC(md);
    1617md=solve(md,'Transient');
    1718md2=solve(md,'Transient','restart',1);
    1819

  • ../trunk-jpl/test/NightlyRun/test129.py

     
    2222# time steps and resolution
    2323md.timestepping.final_time = 19
    2424md.settings.output_frequency = 2
     25md = SetMLHOBC(md);
    2526
    2627md = solve(md, 'Transient')
    2728md2 = copy.deepcopy(md)

  • ../trunk-jpl/test/NightlyRun/test248.m

     
    55md=setflowequation(md,'MLHO','all');
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
     8md=SetMLHOBC(md);
    89md=solve(md,'Stressbalance');
    910
    1011%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test248.py

     
    1515md = setflowequation(md, 'MLHO', 'all')
    1616md.cluster = generic('name', gethostname(), 'np', 3)
    1717md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
     18md = SetMLHOBC(md);
    1819md = solve(md, 'Stressbalance')
    1920
    2021#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test249.m

     
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.transient.requested_outputs={'default','FloatingArea','GroundedArea','TotalGroundedBmb','TotalFloatingBmb'};
    88md.basalforcings.floatingice_melting_rate(:)=1;
     9md=SetMLHOBC(md);
    910md=solve(md,'Transient');
    1011
    1112%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test249.py

     
    1616md = setflowequation(md, 'MLHO', 'all')
    1717md.cluster = generic('name', gethostname(), 'np', 3)
    1818md.transient.requested_outputs = ['default', 'FloatingArea', 'GroundedArea', 'TotalFloatingBmb', 'TotalGroundedBmb']
     19md = SetMLHOBC(md);
    1920md = solve(md, 'Transient')
    2021
    2122

  • ../trunk-jpl/test/NightlyRun/test254.m

     
    5757md.mask.ice_levelset=-1+nodeonicefront;
    5858
    5959md.stressbalance.requested_outputs={'default','VySurface','VyShear','VyBase'};
     60md=SetMLHOBC(md);
    6061md=solve(md,'Stressbalance');
    6162
    6263%create analytical solution: strain rate is constant = ((rho_ice*g*h)/4B)^3 (Paterson, 4th Edition, page 292.

  • ../trunk-jpl/test/NightlyRun/test254.py

     
    6969md.mask.ice_levelset = -1 + nodeonicefront
    7070
    7171md.stressbalance.requested_outputs = ['default', 'VySurface', 'VyShear', 'VyBase']
     72md = SetMLHOBC(md);
    7273md = solve(md, 'Stressbalance')
    7374
    7475# create analytical solution: strain rate is constant = ((rho_ice * g * h) / 4B)^3 (Paterson, 4th Edition, page 292.

  • ../trunk-jpl/test/NightlyRun/test255.m

     
    55md=setflowequation(md,'MLHO','all');
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.masstransport.hydrostatic_adjustment='Incremental';
     8md=SetMLHOBC(md);
    89md=solve(md,'Transient');
    910
    1011%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test255.py

     
    1414md = setflowequation(md, 'MLHO', 'all')
    1515md.cluster = generic('name', gethostname(), 'np', 3)
    1616md.masstransport.hydrostatic_adjustment = 'Incremental'
     17md = SetMLHOBC(md);
    1718md = solve(md, 'Transient')
    1819
    1920#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test256.m

     
    66md.geometry.base=md.geometry.base+50.; md.geometry.surface=md.geometry.surface+50.;
    77md.cluster=generic('name',oshostname(),'np',1);
    88md.masstransport.hydrostatic_adjustment='Incremental';
     9md=SetMLHOBC(md);
    910md=solve(md,'Transient');
    1011
    1112%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test256.py

     
    1616md.geometry.surface = md.geometry.surface + 50.
    1717md.cluster = generic('name', gethostname(), 'np', 1)
    1818md.masstransport.hydrostatic_adjustment = 'Incremental'
     19md = SetMLHOBC(md);
    1920md = solve(md, 'Transient')
    2021
    2122#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test330.m

     
    55md=setflowequation(md,'MLHO','all');
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
     8md=SetMLHOBC(md);
    89md=solve(md,'Stressbalance');
    910
    1011%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test330.py

     
    1414md = setflowequation(md, 'MLHO', 'all')
    1515md.cluster = generic('name', gethostname(), 'np', 3)
    1616md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
     17md = SetMLHOBC(md);
    1718md = solve(md, 'Stressbalance')
    1819
    1920#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test332.m

     
    44md=parameterize(md,'../Par/SquareSheetConstrained.par');
    55md=setflowequation(md,'MLHO','all');
    66md.cluster=generic('name',oshostname(),'np',3);
     7md=SetMLHOBC(md);
    78md=solve(md,'Transient');
    89
    910%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test332.py

     
    1313md = parameterize(md, '../Par/SquareSheetConstrained.py')
    1414md = setflowequation(md, 'MLHO', 'all')
    1515md.cluster = generic('name', gethostname(), 'np', 3)
     16md = SetMLHOBC(md);
    1617md = solve(md, 'Transient')
    1718
    1819#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test333.m

     
    1818md.inversion.vx_obs=md.initialization.vx; md.inversion.vy_obs=md.initialization.vy;
    1919
    2020md.cluster=generic('name',oshostname(),'np',3);
     21md=SetMLHOBC(md);
    2122md=solve(md,'Stressbalance');
    2223
    2324%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test333.py

     
    3131md.inversion.vy_obs = md.initialization.vy
    3232
    3333md.cluster = generic('name', gethostname(), 'np', 3)
     34md = SetMLHOBC(md);
    3435md = solve(md, 'Stressbalance')
    3536
    3637#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test334.m

     
    1919md.verbose.control=true;
    2020
    2121md.cluster=generic('name',oshostname(),'np',3);
     22md=SetMLHOBC(md);
    2223md=solve(md,'Stressbalance');
    2324
    2425%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test334.py

     
    3131
    3232
    3333md.cluster = generic('name', gethostname(), 'np', 3)
     34md = SetMLHOBC(md);
    3435md = solve(md, 'Stressbalance')
    3536
    3637

  • ../trunk-jpl/test/NightlyRun/test335.m

     
    66md = extrude(md, 5, 1);
    77md.cluster=generic('name',oshostname(),'np',3);
    88md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
     9md=SetMLHOBC(md);
    910md=solve(md,'Stressbalance');
    1011
    1112%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test335.py

     
    1212md = setmask(md, '', '')
    1313md = parameterize(md, '../Par/SquareSheetConstrained.py')
    1414md = setflowequation(md, 'MLHO', 'all')
     15md = SetMLHOBC(md);
    1516md.extrude(5, 1.)
    1617md.cluster = generic('name', gethostname(), 'np', 3)
    1718md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']

  • ../trunk-jpl/test/NightlyRun/test446.m

     
    55md=setflowequation(md,'MLHO','all');
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
     8md=SetMLHOBC(md);
    89md=solve(md,'Stressbalance');
    910
    1011%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test446.py

     
    1414md = setflowequation(md, 'MLHO', 'all')
    1515md.cluster = generic('name', gethostname(), 'np', 3)
    1616md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
     17md = SetMLHOBC(md);
    1718md = solve(md, 'Stressbalance')
    1819
    1920#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test447.m

     
    1515md=setflowequation(md,'MLHO','all');
    1616md.cluster=generic('name',oshostname(),'np',3);
    1717md.transient.requested_outputs={'default','GroundedArea','FloatingArea','TotalFloatingBmb','TotalGroundedBmb','TotalSmb'};
     18md=SetMLHOBC(md);
    1819md=solve(md,'Transient');
    1920
    2021%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test447.py

     
    2626md = setflowequation(md, 'MLHO', 'all')
    2727md.cluster = generic('name', gethostname(), 'np', 3)
    2828md.transient.requested_outputs = ['default', 'GroundedArea', 'FloatingArea', 'TotalFloatingBmb', 'TotalGroundedBmb', 'TotalSmb']
     29md = SetMLHOBC(md);
    2930md = solve(md, 'Transient')
    3031
    3132#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test448.m

     
    2424md.transient.isstressbalance=1;
    2525md.transient.isgroundingline=1;
    2626
     27md=SetMLHOBC(md);
    2728%test different grounding line dynamics.
    2829md.groundingline.migration='AggressiveMigration';
    2930md=solve(md,'Transient');

  • ../trunk-jpl/test/NightlyRun/test448.py

     
    3737
    3838#test different grounding line dynamics.
    3939md.groundingline.migration = 'AggressiveMigration'
     40md = SetMLHOBC(md);
    4041md = solve(md, 'Transient')
    4142element_on_iceshelf_agressive = md.results.TransientSolution[0].MaskOceanLevelset
    4243vel_agressive = md.results.TransientSolution[0].Vel

  • ../trunk-jpl/test/NightlyRun/test449.m

     
    4444md.timestepping.time_step=10;
    4545
    4646md.cluster=generic('name',oshostname(),'np',3);
     47md=SetMLHOBC(md);
    4748md=solve(md,'Transient');
    4849
    4950%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test449.py

     
    5555md.timestepping.time_step = 10
    5656
    5757md.cluster = generic('name', gethostname(), 'np', 3)
     58md = SetMLHOBC(md);
    5859md = solve(md, 'Transient')
    5960#print md.results.TransientSolution[0].BasalforcingsFloatingiceMeltingRate
    6061#print md.results.TransientSolution[1].BasalforcingsFloatingiceMeltingRate

  • ../trunk-jpl/test/NightlyRun/test518.m

     
    55md=setflowequation(md,'MLHO','all');
    66md.cluster=generic('name',oshostname(),'np',3);
    77md.stressbalance.requested_outputs={'default','VxSurface','VySurface','VxShear','VyShear','VxBase','VyBase'};
     8md=SetMLHOBC(md);
    89md=solve(md,'Stressbalance');
    910
    1011%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test518.py

     
    1414md = setflowequation(md, 'MLHO', 'all')
    1515md.cluster = generic('name', gethostname(), 'np', 3)
    1616md.stressbalance.requested_outputs = ['default', 'VxSurface', 'VySurface', 'VxShear', 'VyShear', 'VxBase', 'VyBase']
     17md = SetMLHOBC(md);
    1718md = solve(md, 'Stressbalance')
    1819
    1920#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test519.m

     
    88[md.mesh.lat,md.mesh.long] = xy2ll(md.mesh.x,md.mesh.y,-1);
    99md.mesh.scale_factor=0.9*ones(md.mesh.numberofvertices,1);
    1010md.transient.requested_outputs={'default','IceVolume','IceVolumeScaled','GroundedArea','GroundedAreaScaled','FloatingArea','FloatingAreaScaled','TotalSmb','TotalSmbScaled','TotalFloatingBmb','TotalFloatingBmbScaled'};
     11md=SetMLHOBC(md);
    1112md=solve(md,'Transient');
    1213
    1314%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test519.py

     
    1515md.mesh.scale_factor = 0.9 * np.ones((md.mesh.numberofvertices))
    1616md.transient.requested_outputs = ['default', 'IceVolume', 'IceVolumeScaled', 'GroundedArea', 'GroundedAreaScaled', 'FloatingArea', 'FloatingAreaScaled', 'TotalSmb', 'TotalSmbScaled', 'TotalFloatingBmb', 'TotalFloatingBmbScaled']
    1717md.cluster = generic('name', gethostname(), 'np', 3)
     18md = SetMLHOBC(md);
    1819md = solve(md, 'Transient')
    1920
    2021# Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test810.m

     
    1414md.transient.isthermal=0;
    1515md.transient.isgroundingline=1;
    1616
     17md=SetMLHOBC(md);
    1718md=solve(md,'Transient');
    1819
    1920%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test810.py

     
    2121md.transient.issmb = True
    2222md.transient.isthermal = False
    2323md.transient.isgroundingline = True
     24md = SetMLHOBC(md);
    2425
    2526md = solve(md, 'Transient')
    2627

  • ../trunk-jpl/test/NightlyRun/test811.m

     
    1717md.frontalforcings.meltingrate=zeros(md.mesh.numberofvertices,1);
    1818md.levelset.migration_max = 1e10;
    1919
     20md=SetMLHOBC(md);
    2021md=solve(md,'Transient');
    2122
    2223%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test811.py

     
    2727md.frontalforcings.meltingrate = np.zeros((md.mesh.numberofvertices))
    2828md.levelset.migration_max = 1e10
    2929
     30md = SetMLHOBC(md);
    3031md = solve(md, 'Transient')
    3132
    3233#Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test812.m

     
    3434
    3535md.transient.requested_outputs={'default','StrainRateparallel','StrainRateperpendicular','Calvingratex','Calvingratey','CalvingCalvingrate'};
    3636
     37md=SetMLHOBC(md);
    3738md=solve(md,'Transient');
    3839
    3940%Fields and tolerances to track changes

  • ../trunk-jpl/test/NightlyRun/test812.py

     
    4343
    4444md.transient.requested_outputs = ['default', 'StrainRateparallel', 'StrainRateperpendicular', 'Calvingratex', 'Calvingratey', 'CalvingCalvingrate']
    4545
     46md = SetMLHOBC(md);
    4647md = solve(md, 'Transient')
    4748
    4849#Fields and tolerances to track changes

  • ../trunk-jpl/src/m/boundaryconditions/SetMLHOBC.py

     
     1import os
     2import MatlabFuncs as m
     3
     4
     5def SetMLHOBC(md):
     6    """
     7    SETMLHOBC - Create the boundary conditions for stressbalance for MLHO: VxBase, VyBase, VxShear, VyShear
     8
     9       Usage:
     10          md = SetIceShelfBC(md, varargin)
     11
     12       Example:
     13          md = SetIceShelfBC(md)
     14
     15    """
     16
     17    #node on Dirichlet (boundary and ~icefront)
     18    md.stressbalance.spcvx_base = md.stressbalance.spcvx
     19    md.stressbalance.spcvy_base = md.stressbalance.spcvy
     20    md.stressbalance.spcvx_shear = float('nan') * md.stressbalance.spcvx
     21    md.stressbalance.spcvy_shear = float('nan') * md.stressbalance.spcvy

  • ../trunk-jpl/src/m/extrusion/project3d.py

     
    1414
    1515        arguments:
    1616            'vector': 2d vector
    17             'type': 'element' or 'node'
     17            'type': 'element' or 'node' or 'poly'
    1818
    1919        options:
    2020            'layer'     a layer number where vector should keep its values. If
    … …  
    2222                        vector.
    2323            'padding':  default to 0 (value adopted by other 3d layers not
    2424                        being projected.
     25            'degree':   degree of polynomials when extrude from bottom to the top
    2526
    2627        Examples:
    2728            extruded_vector = project3d(md, 'vector', vector2d, 'type', 'node', 'layer', 1, 'padding', NaN)
    … …  
    4142    vectype = options.getfieldvalue('type')  #mandatory
    4243    layer = options.getfieldvalue('layer', 0)  #optional (do all layers otherwise)
    4344    paddingvalue = options.getfieldvalue('padding', 0)  #0 by default
     45    polyexponent = options.getfieldvalue('degree', 0)  #0 by default, 0-degree polynomial
    4446
    4547    #Handle special case where vector2d is single element (differs from representation in MATLAB)
    4648    if isinstance(vector2d, (bool, int, float)):
    … …  
    114116                    projected_vector[(i * md.mesh.numberofelements2d):((i + 1) * md.mesh.numberofelements2d), :] = vector2d
    115117            else:
    116118                projected_vector[((layer - 1) * md.mesh.numberofelements2d):(layer * md.mesh.numberofelements2d), :] = vector2d
     119    elif vectype.lower() == 'poly':
     120        #Initialize 3d vector
     121        if np.ndim(vector2d) == 1:
     122            if vector2d.shape[0] == md.mesh.numberofelements2d:
     123                projected_vector = (paddingvalue * np.ones((md.mesh.numberofelements))).astype(vector2d.dtype)
     124            elif vector2d.shape[0] == md.mesh.numberofelements2d + 1:
     125                projected_vector = (paddingvalue * np.ones((md.mesh.numberofelements + 1))).astype(vector2d.dtype)
     126                projected_vector[-1] = vector2d[-1]
     127                vector2d = vector2d[:-1]
     128            else:
     129                raise TypeError("vector length not supported")
     130            #Fill in
     131            if layer == 0:
     132                for i in range(md.mesh.numberoflayers - 1):
     133                    projected_vector[(i * md.mesh.numberofelements2d):((i + 1) * md.mesh.numberofelements2d)] = vector2d*(1.0-(1.0-i/(md.mesh.numberoflayers - 1.0))**polyexponent)
     134            else:
     135                projected_vector[((layer - 1) * md.mesh.numberofelements2d):(layer * md.mesh.numberofelements2d)] = vector2d*(1.0-(1.0-layer/(md.mesh.numberoflayers - 1.0))**polyexponent)
     136        else:
     137            if vector2d.shape[0] == md.mesh.numberofelements2d:
     138                projected_vector = (paddingvalue * np.ones((md.mesh.numberofelements, np.size(vector2d, axis=1)))).astype(vector2d.dtype)
     139            elif vector2d.shape[0] == md.mesh.numberofelements2d + 1:
     140                projected_vector = (paddingvalue * np.ones((md.mesh.numberofelements + 1, np.size(vector2d, axis=1)))).astype(vector2d.dtype)
     141                projected_vector[-1, :] = vector2d[-1, :]
     142                vector2d = vector2d[:-1, :]
     143            else:
     144                raise TypeError("vector length not supported")
     145            #Fill in
     146            if layer == 0:
     147                for i in range(md.mesh.numberoflayers - 1):
     148                    projected_vector[(i * md.mesh.numberofelements2d):((i + 1) * md.mesh.numberofelements2d), :] = vector2d*(1.0-(1.0-i/(md.mesh.numberoflayers - 1.0))**polyexponent)
     149            else:
    117150
    118151    else:
    119152        raise TypeError("project3d error message: unknown projection type")

  • ../trunk-jpl/test/NightlyRun/test127.py

     
    2424                                   massfluxatgate('name', 'MassFlux4', 'profilename', '../Exp/MassFlux4.exp', 'definitionstring', 'Outputdefinition4'),
    2525                                   massfluxatgate('name', 'MassFlux5', 'profilename', '../Exp/MassFlux5.exp', 'definitionstring', 'Outputdefinition5'),
    2626                                   massfluxatgate('name', 'MassFlux6', 'profilename', '../Exp/MassFlux6.exp', 'definitionstring', 'Outputdefinition6')]
    27 
     27md = SetMLHOBC(md);
    2828md = solve(md, 'Stressbalance')
    2929
    3030#Fields and tolerances to track changes
Note: See TracBrowser for help on using the repository browser.