Changeset 21626


Ignore:
Timestamp:
03/23/17 15:47:16 (8 years ago)
Author:
youngmc3
Message:

CHG: added calving stress thresholds and calving max parameters

Location:
issm/trunk-jpl/src
Files:
12 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp

    r20983 r21626  
    6161                        break;
    6262                case CalvingDevEnum:
    63                         iomodel->FetchDataToInput(elements,"md.calving.coeff",CalvingdevCoeffEnum);
    6463                        iomodel->FetchDataToInput(elements,"md.calving.meltingrate",CalvingMeltingrateEnum);
    6564                        break;
    … …  
    7574        parameters->AddObject(iomodel->CopyConstantObject("md.levelset.stabilization",LevelsetStabilizationEnum));
    7675        parameters->AddObject(iomodel->CopyConstantObject("md.levelset.reinit_frequency",LevelsetReinitFrequencyEnum));
     76        parameters->AddObject(iomodel->CopyConstantObject("md.levelset.calving_max",CalvingMaxEnum));
    7777        int  calvinglaw;
    7878        iomodel->FindConstant(&calvinglaw,"md.calving.law");
    … …  
    8080                case DefaultCalvingEnum:
    8181                case CalvingLevermannEnum:
     82                        break;
    8283                case CalvingDevEnum:
     84                        parameters->AddObject(iomodel->CopyConstantObject("md.calving.stress_threshold_groundedice",CalvingStressThresholdGroundediceEnum));
     85                        parameters->AddObject(iomodel->CopyConstantObject("md.calving.stress_threshold_floatingice",CalvingStressThresholdFloatingiceEnum));
    8386                        break;
    8487                case CalvingMinthicknessEnum:
    … …  
    138141        IssmDouble h,hx,hy,hz;
    139142        IssmDouble vel;
    140         IssmDouble norm_dlsf, norm_calving, calvingrate, meltingrate;
     143        IssmDouble norm_dlsf, norm_calving, calvingrate, meltingrate, groundedice;
     144        IssmDouble calvingmax;
    141145        IssmDouble* xyz_list = NULL;
    142146
    … …  
    152156        }
    153157
     158        /*Calving threshold*/
     159       
    154160        /*Fetch number of nodes and dof for this finite element*/
    155161        int numnodes    = basalelement->GetNumberOfNodes();
    … …  
    170176        basalelement->GetVerticesCoordinates(&xyz_list);
    171177        basalelement->FindParam(&dt,TimesteppingTimeStepEnum);
     178        basalelement->FindParam(&calvingmax,CalvingMaxEnum);
    172179        Input* vx_input           = NULL;
    173180        Input* vy_input           = NULL;
    … …  
    178185        Input* calvingrate_input  = NULL;
    179186        Input* meltingrate_input  = NULL;
     187        Input* gr_input           = NULL;
    180188
    181189        /*Load velocities*/
    … …  
    187195                        vx_input=basalelement->GetInput(VxEnum); _assert_(vx_input);
    188196                        vy_input=basalelement->GetInput(VyEnum); _assert_(vy_input);
     197                        gr_input=basalelement->GetInput(MaskGroundediceLevelsetEnum); _assert_(gr_input);
    189198                        break;
    190199                case Domain3DEnum:
    191200                        vx_input=basalelement->GetInput(VxAverageEnum); _assert_(vx_input);
    192201                        vy_input=basalelement->GetInput(VyAverageEnum); _assert_(vy_input);
     202                        gr_input=basalelement->GetInput(MaskGroundediceLevelsetEnum); _assert_(gr_input);
    193203                        break;
    194204                default: _error_("mesh "<<EnumToStringx(domaintype)<<" not supported yet");
    … …  
    253263                vx_input->GetInputValue(&v[0],gauss);
    254264                vy_input->GetInputValue(&v[1],gauss);
     265                gr_input->GetInputValue(&groundedice,gauss);
    255266
    256267                /*Get calving speed*/
    … …  
    262273                                calvingrate_input->GetInputValue(&calvingrate,gauss);
    263274                                meltingrate_input->GetInputValue(&meltingrate,gauss);
     275
     276                                /*Limit calving rate to c <= v + 3 km/yr */
     277                                vel=sqrt(v[0]*v[0] + v[1]*v[1]);
     278                                if(calvingrate>calvingmax+vel) calvingrate = vel+calvingmax;
     279                                if(groundedice<0) meltingrate = 0;
    264280
    265281                                norm_dlsf=0.;

  • issm/trunk-jpl/src/c/classes/Elements/Penta.cpp

    r21583 r21626  
    186186        IssmDouble  calvingrate[NUMVERTICES];
    187187        IssmDouble  lambda1,lambda2,ex,ey,vx,vy,vel;
    188         IssmDouble  sigma_vm,sigma_max,epse_2,groundedice;
     188        IssmDouble  sigma_vm,sigma_max,sigma_max_floating,sigma_max_grounded;
     189        IssmDouble  epse_2,groundedice;
    189190
    190191        /* Get node coordinates and dof list: */
    … …  
    202203        IssmDouble  B   = this->GetMaterialParameter(MaterialsRheologyBbarEnum);
    203204        IssmDouble  n   = this->GetMaterialParameter(MaterialsRheologyNEnum);
     205        this->parameters->FindParam(&sigma_max_floating,CalvingStressThresholdFloatingiceEnum);
     206        this->parameters->FindParam(&sigma_max_grounded,CalvingStressThresholdGroundediceEnum);
    204207
    205208        /* Start looping on the number of vertices: */
    … …  
    229232                epse_2    = 1./2. *(lambda1*lambda1 + lambda2*lambda2);
    230233                sigma_vm  = sqrt(3.) * B * pow(epse_2,1./(2.*n));
    231                 sigma_max = 1000.e+3;
    232 
    233                 if(groundedice<0) sigma_max=200.e+3;
     234
     235                /*Tensile stress threshold*/
     236                if(groundedice<0)
     237                 sigma_max = sigma_max_floating;
     238                else
     239                 sigma_max = sigma_max_grounded;
    234240
    235241                /*Assign values*/

  • issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r21620 r21626  
    216216        IssmDouble  calvingrate[NUMVERTICES];
    217217        IssmDouble  lambda1,lambda2,ex,ey,vx,vy,vel;
    218         IssmDouble  sigma_vm,sigma_max,epse_2,groundedice;
     218        IssmDouble  sigma_vm,sigma_max,sigma_max_floating,sigma_max_grounded;
     219        IssmDouble  epse_2,groundedice;
    219220
    220221        /* Get node coordinates and dof list: */
    … …  
    227228        IssmDouble  B   = this->GetMaterialParameter(MaterialsRheologyBbarEnum);
    228229        IssmDouble  n   = this->GetMaterialParameter(MaterialsRheologyNEnum);
     230        this->parameters->FindParam(&sigma_max_floating,CalvingStressThresholdFloatingiceEnum);
     231        this->parameters->FindParam(&sigma_max_grounded,CalvingStressThresholdGroundediceEnum);
    229232
    230233        /* Start looping on the number of vertices: */
    … …  
    254257                epse_2    = 1./2. *(lambda1*lambda1 + lambda2*lambda2);
    255258                sigma_vm  = sqrt(3.) * B * pow(epse_2,1./(2.*n));
    256                 //sigma_max = 125.e+3;
    257                 sigma_max = 350.e+3;
    258                 sigma_max = 450.e+3;
    259                 sigma_max = 800.e+3; //too much
    260                 //sigma_max = 700.e+3;
    261                 //sigma_max = 670.e+3;
    262                 //sigma_max = 550.e+3;
    263                 sigma_max = 750.e+3; //too high
    264                 sigma_max = 850.e+3; //too low
    265                 sigma_max = 800.e+3; //IUGG previous test
    266                 sigma_max = 1000.e+3; //850 seems small
    267 
    268                 if(groundedice<0) sigma_max=200.e+3;
     259
     260                /*OLD (keep for a little bit)*/
     261                //sigma_max = 800.e+3; //IUGG previous test
     262                //sigma_max = 1000.e+3; //GRL
     263                //if(groundedice<0) sigma_max=150.e+3;
     264
     265                /*Tensile stress threshold*/
     266                if(groundedice<0)
     267                 sigma_max = sigma_max_floating;
     268                else
     269                 sigma_max = sigma_max_grounded;
    269270
    270271                /*Assign values*/
    … …  
    294295        IssmDouble  calvingratey[NUMVERTICES];
    295296        IssmDouble  calvingrate[NUMVERTICES];
    296 
    297297
    298298        /* Get node coordinates and dof list: */

  • issm/trunk-jpl/src/c/cores/levelsetfunctionslope_core.cpp

    r17700 r21626  
    3333        }
    3434        if(domaintype==Domain2DverticalEnum){
    35               femmodel->parameters->SetParam(LevelsetfunctionSlopeXEnum,InputToExtrudeEnum);
     35                femmodel->parameters->SetParam(LevelsetfunctionSlopeXEnum,InputToExtrudeEnum);
    3636                extrudefrombase_core(femmodel);
    3737        }

  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r21620 r21626  
    259259        CalvingratexAverageEnum,
    260260        CalvingrateyAverageEnum,
     261        CalvingStressThresholdGroundediceEnum,
     262        CalvingStressThresholdFloatingiceEnum,
     263        CalvingMaxEnum,
    261264        StrainRateparallelEnum,
    262265        StrainRateperpendicularEnum,

  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r21620 r21626  
    265265                case CalvingratexAverageEnum : return "CalvingratexAverage";
    266266                case CalvingrateyAverageEnum : return "CalvingrateyAverage";
     267                case CalvingStressThresholdGroundediceEnum : return "CalvingStressThresholdGroundedice";
     268                case CalvingStressThresholdFloatingiceEnum : return "CalvingStressThresholdFloatingice";
     269                case CalvingMaxEnum : return "CalvingMax";
    267270                case StrainRateparallelEnum : return "StrainRateparallel";
    268271                case StrainRateperpendicularEnum : return "StrainRateperpendicular";

  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r21620 r21626  
    271271              else if (strcmp(name,"CalvingratexAverage")==0) return CalvingratexAverageEnum;
    272272              else if (strcmp(name,"CalvingrateyAverage")==0) return CalvingrateyAverageEnum;
     273              else if (strcmp(name,"CalvingStressThresholdGroundedice")==0) return CalvingStressThresholdGroundediceEnum;
     274              else if (strcmp(name,"CalvingStressThresholdFloatingice")==0) return CalvingStressThresholdFloatingiceEnum;
     275              else if (strcmp(name,"CalvingMax")==0) return CalvingMaxEnum;
    273276              else if (strcmp(name,"StrainRateparallel")==0) return StrainRateparallelEnum;
    274277              else if (strcmp(name,"StrainRateperpendicular")==0) return StrainRateperpendicularEnum;
    … …  
    380383              else if (strcmp(name,"SmbRequestedOutputs")==0) return SmbRequestedOutputsEnum;
    381384              else if (strcmp(name,"SmbIsInitialized")==0) return SmbIsInitializedEnum;
    382               else if (strcmp(name,"SmbDzini")==0) return SmbDziniEnum;
    383               else if (strcmp(name,"SmbDini")==0) return SmbDiniEnum;
    384               else if (strcmp(name,"SmbReini")==0) return SmbReiniEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"SmbGdnini")==0) return SmbGdniniEnum;
     388              if (strcmp(name,"SmbDzini")==0) return SmbDziniEnum;
     389              else if (strcmp(name,"SmbDini")==0) return SmbDiniEnum;
     390              else if (strcmp(name,"SmbReini")==0) return SmbReiniEnum;
     391              else if (strcmp(name,"SmbGdnini")==0) return SmbGdniniEnum;
    389392              else if (strcmp(name,"SmbGspini")==0) return SmbGspiniEnum;
    390393              else if (strcmp(name,"SmbECini")==0) return SmbECiniEnum;
    … …  
    503506              else if (strcmp(name,"SurfaceAbsVelMisfit")==0) return SurfaceAbsVelMisfitEnum;
    504507              else if (strcmp(name,"SurfaceArea")==0) return SurfaceAreaEnum;
    505               else if (strcmp(name,"SurfaceAverageVelMisfit")==0) return SurfaceAverageVelMisfitEnum;
    506               else if (strcmp(name,"SurfaceLogVelMisfit")==0) return SurfaceLogVelMisfitEnum;
    507               else if (strcmp(name,"SurfaceLogVxVyMisfit")==0) return SurfaceLogVxVyMisfitEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"SurfaceRelVelMisfit")==0) return SurfaceRelVelMisfitEnum;
     511              if (strcmp(name,"SurfaceAverageVelMisfit")==0) return SurfaceAverageVelMisfitEnum;
     512              else if (strcmp(name,"SurfaceLogVelMisfit")==0) return SurfaceLogVelMisfitEnum;
     513              else if (strcmp(name,"SurfaceLogVxVyMisfit")==0) return SurfaceLogVxVyMisfitEnum;
     514              else if (strcmp(name,"SurfaceRelVelMisfit")==0) return SurfaceRelVelMisfitEnum;
    512515              else if (strcmp(name,"SurfaceSlopeX")==0) return SurfaceSlopeXEnum;
    513516              else if (strcmp(name,"SurfaceSlopeY")==0) return SurfaceSlopeYEnum;
    … …  
    626629              else if (strcmp(name,"Outputdefinition31")==0) return Outputdefinition31Enum;
    627630              else if (strcmp(name,"Outputdefinition32")==0) return Outputdefinition32Enum;
    628               else if (strcmp(name,"Outputdefinition33")==0) return Outputdefinition33Enum;
    629               else if (strcmp(name,"Outputdefinition34")==0) return Outputdefinition34Enum;
    630               else if (strcmp(name,"Outputdefinition35")==0) return Outputdefinition35Enum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"Outputdefinition36")==0) return Outputdefinition36Enum;
     634              if (strcmp(name,"Outputdefinition33")==0) return Outputdefinition33Enum;
     635              else if (strcmp(name,"Outputdefinition34")==0) return Outputdefinition34Enum;
     636              else if (strcmp(name,"Outputdefinition35")==0) return Outputdefinition35Enum;
     637              else if (strcmp(name,"Outputdefinition36")==0) return Outputdefinition36Enum;
    635638              else if (strcmp(name,"Outputdefinition37")==0) return Outputdefinition37Enum;
    636639              else if (strcmp(name,"Outputdefinition38")==0) return Outputdefinition38Enum;
    … …  
    749752              else if (strcmp(name,"Scaled")==0) return ScaledEnum;
    750753              else if (strcmp(name,"Separate")==0) return SeparateEnum;
    751               else if (strcmp(name,"Sset")==0) return SsetEnum;
    752               else if (strcmp(name,"Dense")==0) return DenseEnum;
    753               else if (strcmp(name,"MpiDense")==0) return MpiDenseEnum;
    754754         else stage=7;
    755755   }
    756756   if(stage==7){
    757               if (strcmp(name,"MpiSparse")==0) return MpiSparseEnum;
     757              if (strcmp(name,"Sset")==0) return SsetEnum;
     758              else if (strcmp(name,"Dense")==0) return DenseEnum;
     759              else if (strcmp(name,"MpiDense")==0) return MpiDenseEnum;
     760              else if (strcmp(name,"MpiSparse")==0) return MpiSparseEnum;
    758761              else if (strcmp(name,"Seq")==0) return SeqEnum;
    759762              else if (strcmp(name,"Mpi")==0) return MpiEnum;
    … …  
    872875              else if (strcmp(name,"Penta")==0) return PentaEnum;
    873876              else if (strcmp(name,"PentaInput")==0) return PentaInputEnum;
    874               else if (strcmp(name,"Vertex")==0) return VertexEnum;
    875               else if (strcmp(name,"VertexPId")==0) return VertexPIdEnum;
    876               else if (strcmp(name,"VertexSId")==0) return VertexSIdEnum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"Option")==0) return OptionEnum;
     880              if (strcmp(name,"Vertex")==0) return VertexEnum;
     881              else if (strcmp(name,"VertexPId")==0) return VertexPIdEnum;
     882              else if (strcmp(name,"VertexSId")==0) return VertexSIdEnum;
     883              else if (strcmp(name,"Option")==0) return OptionEnum;
    881884              else if (strcmp(name,"GenericOption")==0) return GenericOptionEnum;
    882885              else if (strcmp(name,"OptionCell")==0) return OptionCellEnum;
    … …  
    995998              else if (strcmp(name,"LACrouzeixRaviart")==0) return LACrouzeixRaviartEnum;
    996999              else if (strcmp(name,"Closed")==0) return ClosedEnum;
    997               else if (strcmp(name,"Free")==0) return FreeEnum;
    998               else if (strcmp(name,"Open")==0) return OpenEnum;
    999               else if (strcmp(name,"Air")==0) return AirEnum;
    10001000         else stage=9;
    10011001   }
    10021002   if(stage==9){
    1003               if (strcmp(name,"Ice")==0) return IceEnum;
     1003              if (strcmp(name,"Free")==0) return FreeEnum;
     1004              else if (strcmp(name,"Open")==0) return OpenEnum;
     1005              else if (strcmp(name,"Air")==0) return AirEnum;
     1006              else if (strcmp(name,"Ice")==0) return IceEnum;
    10041007              else if (strcmp(name,"Melange")==0) return MelangeEnum;
    10051008              else if (strcmp(name,"Water")==0) return WaterEnum;

  • issm/trunk-jpl/src/m/classes/calvingdev.m

    r21049 r21626  
    66classdef calvingdev
    77        properties (SetAccess=public)
    8                 coeff         = NaN;
     8                stress_threshold_groundedice = 0.;
     9                stress_threshold_floatingice = 0.;
    910                meltingrate   = NaN;
    1011        end
    … …  
    2930                end % }}}
    3031                function self = extrude(self,md) % {{{
    31                         self.coeff=project3d(md,'vector',self.coeff,'type','node');
    3232                        self.meltingrate=project3d(md,'vector',self.meltingrate,'type','node');
    3333                end % }}}
    … …  
    3535
    3636                        %Proportionality coefficient in Pi model
    37                         self.coeff=2e13;
     37                        %self.coeff=2e13;
     38
     39                        %Default sigma max
     40                        self.stress_threshold_groundedice = 1e6;
     41                        self.stress_threshold_floatingice = 150e3;
    3842                end % }}}
    3943                function md = checkconsistency(self,md,solution,analyses) % {{{
    … …  
    4145                        if (~strcmp(solution,'TransientSolution') | md.transient.ismovingfront==0), return; end
    4246
    43                         md = checkfield(md,'fieldname','calving.coeff','>',0,'size',[md.mesh.numberofvertices 1]);
     47                        md = checkfield(md,'fieldname','calving.stress_threshold_groundedice','>',0,'numel',1,'nan',1,'Inf',1);
     48                        md = checkfield(md,'fieldname','calving.stress_threshold_floatingice','>',0,'numel',1,'nan',1,'Inf',1);
    4449                        md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'Inf',1,'timeseries',1,'>=',0);
    4550                end % }}}
    4651                function disp(self) % {{{
    4752                        disp(sprintf('   Calving Pi parameters:'));
    48                         fielddisplay(self,'coeff','proportionality coefficient in Pi model');
     53                        fielddisplay(self,'stress_threshold_groundedice','sigma_max applied to grounded ice only [Pa]');
     54                        fielddisplay(self,'stress_threshold_floatingice','sigma_max applied to floating ice only [Pa]');
    4955                        fielddisplay(self,'meltingrate','melting rate at given location [m/a]');
    5056
    … …  
    5359                        yts=md.constants.yts;
    5460                        WriteData(fid,prefix,'name','md.calving.law','data',2,'format','Integer');
    55                         WriteData(fid,prefix,'object',self,'fieldname','coeff','format','DoubleMat','mattype',1);
    56                         WriteData(fid,prefix,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
     61                        WriteData(fid,prefix,'object',self,'fieldname','stress_threshold_groundedice','format','DoubleMat','mattype',1);
     62                        WriteData(fid,prefix,'object',self,'fieldname','stress_threshold_floatingice','format','DoubleMat','mattype',1);
     63                        WriteData(fid,prefix,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'scale',1./yts);
    5764                end % }}}
    5865        end

  • issm/trunk-jpl/src/m/classes/calvinglevermann.m

    r21049 r21626  
    4242
    4343                        md = checkfield(md,'fieldname','calving.coeff','>',0,'size',[md.mesh.numberofvertices 1]);
    44                         md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1],'>=',0);
     44                        md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'Inf',1,'timeseries',1,'>=',0);
    4545                end % }}}
    4646                function disp(self) % {{{
    … …  
    5454                        WriteData(fid,prefix,'name','md.calving.law','data',3,'format','Integer');
    5555                        WriteData(fid,prefix,'object',self,'fieldname','coeff','format','DoubleMat','mattype',1);
    56                         WriteData(fid,prefix,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts,'scale',1./yts);
     56                        WriteData(fid,prefix,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'scale',1./yts);
    5757                end % }}}
    5858        end

  • issm/trunk-jpl/src/m/classes/calvingminthickness.m

    r21049 r21626  
    5353                        WriteData(fid,prefix,'name','md.calving.law','data',4,'format','Integer');
    5454                        WriteData(fid,prefix,'object',self,'fieldname','min_thickness','format','Double');
    55                         WriteData(fid,prefix,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
     55                        WriteData(fid,prefix,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'scale',1./yts);
    5656                end % }}}
    5757        end

  • issm/trunk-jpl/src/m/classes/levelset.m

    r21049 r21626  
    99                spclevelset                     = NaN;
    1010                reinit_frequency        = 5;
     11                calving_max       = 0.;
    1112        end
    1213        methods
    … …  
    3637
    3738                        %stabilization = 2 by default
    38                         self.stabilization = 2;
     39                        self.stabilization    = 2;
    3940                        self.reinit_frequency = 5;
     41                        self.calving_max      = 3000.;
    4042
    4143                end % }}}
    … …  
    4648                        md = checkfield(md,'fieldname','levelset.spclevelset','Inf',1,'timeseries',1);
    4749                        md = checkfield(md,'fieldname','levelset.stabilization','values',[0 1 2]);
     50                        md = checkfield(md,'fieldname','levelset.calving_max','numel',1,'NaN',1,'Inf',1,'>',0);
    4851                end % }}}
    4952                function disp(self) % {{{
    … …  
    5255                        fielddisplay(self,'spclevelset','Levelset constraints (NaN means no constraint)');
    5356                        fielddisplay(self,'reinit_frequency','Amount of time steps after which the levelset function in re-initialized');
     57                        fielddisplay(self,'calving_max','maximum allowed calving rate (m/a)');
    5458                end % }}}
    5559                function marshall(self,prefix,md,fid) % {{{
     60
     61                        yts=md.constants.yts;
     62
    5663                        WriteData(fid,prefix,'object',self,'fieldname','stabilization','format','Integer');
    5764                        WriteData(fid,prefix,'object',self,'fieldname','spclevelset','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
    5865                        WriteData(fid,prefix,'object',self,'fieldname','reinit_frequency','format','Integer');
     66                        WriteData(fid,prefix,'object',self,'fieldname','calving_max','format','Double','scale',1./yts);
    5967                end % }}}
    6068                function savemodeljs(self,fid,modelname) % {{{
    … …  
    6371                        writejs1Darray(fid,[modelname '.levelset.spclevelset'],self.spclevelset);
    6472                        writejs1Darray(fid,[modelname '.levelset.reinit_frequency'],self.reinit_frequency);
     73                        writejsdouble(fid,[modelname '.levelset.calving_max'],self.calving_max);
    6574
    6675                end % }}}

  • issm/trunk-jpl/src/m/classes/levelset.py

    r21049 r21626  
    1414        def __init__(self): # {{{
    1515
    16                 self.stabilization = 0
    17                 self.spclevelset   = float('NaN')
     16                self.stabilization    = 0
     17                self.spclevelset      = float('NaN')
    1818                self.reinit_frequency = 0
     19                self.calving_max      = 0.
    1920
    2021                #set defaults
    … …  
    2728                string="%s\n%s"%(string,fielddisplay(self,'spclevelset','levelset constraints (NaN means no constraint)'))
    2829                string="%s\n%s"%(string,fielddisplay(self,'reinit_frequency','Amount of time steps after which the levelset function in re-initialized'))
     30                string="%s\n%s"%(string,fielddisplay(self,'calving_max','maximum allowed calving rate (m/a)'))
    2931
    3032                return string
    … …  
    3941                self.stabilization = 2
    4042                self.reinit_frequency = 5
     43                self.calving_max      = 3000
    4144
    4245                return self
    … …  
    5053                md = checkfield(md,'fieldname','levelset.spclevelset','Inf',1,'timeseries',1)
    5154                md = checkfield(md,'fieldname','levelset.stabilization','values',[0,1,2]);
     55                md = checkfield(md,'fieldname','levelset.calving_max','numel',1,'NaN',1,'Inf',1,'>',0);
    5256
    5357                return md
    … …  
    5559        def marshall(self,prefix,md,fid):    # {{{
    5660
     61                yts=md.constants.yts;
     62
    5763                WriteData(fid,prefix,'object',self,'fieldname','stabilization','format','Integer');
    5864                WriteData(fid,prefix,'object',self,'fieldname','spclevelset','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
    5965                WriteData(fid,prefix,'object',self,'fieldname','reinit_frequency','format','Integer');
     66                WriteData(fid,prefix,'object',self,'fieldname','calving_max','format','Double','scale',1./yts);
    6067        # }}}
Note: See TracChangeset for help on using the changeset viewer.