Changeset 26126


Ignore:
Timestamp:
03/22/21 09:39:18 (4 years ago)
Author:
Eric.Larour
Message:

CHG: fixing bugs in optimized sea level core. Running side by side with optim flag.

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

Legend:

Unmodified
Added
Removed

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

    r26121 r26126  
    6969        InputUpdateFromConstantx(inputs,elements,0.,BedNorthEnum);
    7070    iomodel->FetchDataToInput(inputs,elements,"md.initialization.sealevel",SealevelEnum);
     71
     72        /*Initialize loads:*/
     73        iomodel->ConstantToInput(inputs,elements,0.,DeltaTwsEnum,P1Enum);
     74        iomodel->ConstantToInput(inputs,elements,0.,DeltaIceThicknessEnum,P1Enum);
     75        iomodel->ConstantToInput(inputs,elements,0.,DeltaBottomPressureEnum,P1Enum);
    7176
    7277
    … …  
    140145        parameters->AddObject(iomodel->CopyConstantObject("md.solidearth.settings.glfraction",SolidearthSettingsGlfractionEnum));
    141146        parameters->AddObject(iomodel->CopyConstantObject("md.solidearth.settings.cross_section_shape",SolidearthSettingsCrossSectionShapeEnum));
     147        parameters->AddObject(iomodel->CopyConstantObject("md.solidearth.settings.optim",SolidearthSettingsOptimEnum));
    142148        parameters->AddObject(new DoubleParam(CumBslcEnum,0.0));
    143149        parameters->AddObject(new DoubleParam(CumBslcIceEnum,0.0));

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

    r26121 r26126  
    65386538
    65396539        /*recover some parameters:*/
    6540         rho_ice=FindParam(MaterialsRhoIceEnum);
    6541         rho_water=FindParam(MaterialsRhoSeawaterEnum);
     6540        this->parameters->FindParam(&rho_ice,MaterialsRhoIceEnum);
     6541        this->parameters->FindParam(&rho_water,MaterialsRhoSeawaterEnum);
     6542        this->parameters->FindParam(&rho_freshwater,MaterialsRhoFreshwaterEnum);
    65426543        this->parameters->FindParam(&computerigid,SolidearthSettingsRigidEnum);
    65436544        this->parameters->FindParam(&glfraction,SolidearthSettingsGlfractionEnum);
    … …  
    66886689
    66896690        /*compute ocean average over element:*/
    6690         OceanAverageOptim(&oceanaverage,&oceanarea,SealevelGRD,masks);
     6691        OceanAverageOptim(&oceanaverage,&oceanarea,&SealevelGRD[0],masks);
    66916692       
    66926693        /*add ocean average in the global sealevelloads vector:*/
    … …  
    68706871                total_weight+=gauss->weight;
    68716872        }
    6872         Sg_avg=Sg_avg/total_weight;
     6873        if(total_weight) Sg_avg=Sg_avg/total_weight;
    68736874        delete gauss;
    68746875
    … …  
    68816882void       Tria::SealevelchangeMomentOfInertiaElement(IssmDouble* dI_list, IssmDouble* loads, IssmDouble* sealevelloads){/*{{{*/
    68826883               
    6883         IssmDouble S;
     6884        IssmDouble S=0;
    68846885
    68856886        /*Compute area of element:*/
    … …  
    69346935
    69356936        /*recover total load: */
    6936         S=sealevelloads[this->Sid()]+loads[this->Sid()];
     6937        if(loads) S+=loads[this->Sid()];
     6938        if(sealevelloads) S+=sealevelloads[this->Sid()];
    69376939
    69386940        /* Perturbation terms for moment of inertia (moi_list):

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

    r26121 r26126  
    2929        /*Parameters, variables:*/
    3030        bool save_results;
     31        int  optim=0;
    3132
    3233        /*Retrieve parameters:*/
    3334        femmodel->parameters->FindParam(&save_results,SaveResultsEnum);
     35        femmodel->parameters->FindParam(&optim,SolidearthSettingsOptimEnum);
    3436       
    3537        /*Verbose: */
    … …  
    4951
    5052        /*Run geodetic:*/
    51         grd_core(femmodel);
     53        if(!optim) grd_core(femmodel);
     54        else grd_core_optim(femmodel);
    5255
    5356        /*Run steric core for sure:*/
    … …  
    413416        IssmDouble           eps_rel;
    414417        IssmDouble           eps_abs;
     418        int                  max_nonlinear_iterations;
     419        int                  iterations=0;
    415420        int                  step;
    416421        IssmDouble           time;
    … …  
    421426        int  count,frequency,iscoupling;
    422427        int  grd=0;
    423         int computesealevel=0;
     428        int  computesealevel=0;
     429
    424430        /*}}}*/
    425431
    … …  
    432438        femmodel->parameters->FindParam(&count,SealevelchangeRunCountEnum);
    433439        femmodel->parameters->FindParam(&computesealevel,SolidearthSettingsComputesealevelchangeEnum);
     440        femmodel->parameters->FindParam(&max_nonlinear_iterations,SolidearthSettingsMaxiterEnum);
    434441        /*}}}*/
    435442
    … …  
    462469        /*call masks core: */
    463470        masks=sealevel_masks(femmodel);
     471        if(VerboseSolution()) _printf0_("         starting  GRD convolutions\n");
    464472       
    465473        /*buildup loads: */
    … …  
    501509        //compute rotation axis motion:
    502510        RotationAxisMotion(&rotationaxismotionvector[0],femmodel,allloads,allsealevelloads);
    503 
     511       
     512        if(VerboseSolution()) _printf0_("         converging ocean GRD convolutions\n");
    504513        for(;;){
    505514                       
    … …  
    525534                slcconvergence(&converged,sealevelloads,oldsealevelloads,eps_rel,eps_abs);
    526535                if (converged)break;
     536
     537                //early return?
     538                if(iterations>max_nonlinear_iterations)break;
     539                else iterations++;
    527540        }
    528541
    529542        deformation:
     543       
     544        if(VerboseSolution()) _printf0_("         deformation GRD convolutions\n");
    530545
    531546        /*convolve loads and sea level loads to get the deformation:*/
    … …  
    534549                element->SealevelchangeDeformationConvolution(allsealevelloads, allloads, rotationaxismotionvector,masks);
    535550        }
     551       
     552        if(VerboseSolution()) _printf0_("         updating GRD fields\n");
    536553
    537554        /*Update bedrock motion and geoid:*/
    … …  
    543560                barycontrib->Save(femmodel->results,femmodel->parameters,oceanarea);
    544561        }
    545                
    546         femmodel->inputs->AXPY(SealevelEnum,SealevelGRDEnum,1);
    547         femmodel->inputs->AXPY(BedEnum,BedGRDEnum,1);
     562
     563        femmodel->inputs->AXPY(1,SealevelEnum,SealevelGRDEnum);
     564        femmodel->inputs->AXPY(1,BedEnum,BedGRDEnum);
    548565        if(horiz){
    549                 femmodel->inputs->AXPY(BedEastEnum,BedEastGRDEnum,1);
    550                 femmodel->inputs->AXPY(BedNorthEnum,BedNorthGRDEnum,1);
     566                femmodel->inputs->AXPY(1,BedEastEnum,BedEastGRDEnum);
     567                femmodel->inputs->AXPY(1,BedNorthEnum,BedNorthGRDEnum);
    551568        }
    552569
    … …  
    594611        int  horiz;
    595612        bool geometrydone = false;
    596         bool optim=false;
     613        int  optim;
    597614        int  grdmodel=0;
    598615
    … …  
    602619        femmodel->parameters->FindParam(&geometrydone,SealevelchangeGeometryDoneEnum);
    603620        femmodel->parameters->FindParam(&grdmodel,GrdModelEnum);
     621        femmodel->parameters->FindParam(&optim,SolidearthSettingsOptimEnum);
    604622       
    605623        /*early return?:*/
    … …  
    623641        for(Object* & object : femmodel->elements->objects){
    624642                Element*   element=xDynamicCast<Element*>(object);
    625                 element->SealevelchangeGeometry(latitude,longitude,radius,xx,yy,zz,xxe,yye,zze);
    626                 //if(optim) element->SealevelchangeGeometryOptim(latitude,longitude,radius,xx,yy,zz,xxe,yye,zze);
     643                if(optim) element->SealevelchangeGeometryOptim(latitude,longitude,radius,xx,yy,zz,xxe,yye,zze);
     644                else      element->SealevelchangeGeometry(latitude,longitude,radius,xx,yy,zz,xxe,yye,zze);
    627645        }
    628646
    … …  
    12401258        return sealevelloadsaverage/oceanarea;
    12411259} /*}}}*/
    1242 void RotationAxisMotion(IssmDouble* m, FemModel* femmodel,IssmDouble* loads, IssmDouble* sealeveloads){ /*{{{*/
     1260void RotationAxisMotion(IssmDouble* m, FemModel* femmodel,IssmDouble* loads, IssmDouble* sealevelloads){ /*{{{*/
    12431261
    12441262        IssmDouble  moi_list[3]={0,0,0};
    … …  
    12611279        for(Object* & object : femmodel->elements->objects){
    12621280                Element* element = xDynamicCast<Element*>(object);
    1263                 element->SealevelchangeMomentOfInertiaElement(&moi_list[0],loads,sealeveloads);
     1281                element->SealevelchangeMomentOfInertiaElement(&moi_list[0],loads,sealevelloads);
    12641282                moi_list_cpu[0] += moi_list[0];
    12651283                moi_list_cpu[1] += moi_list[1];

  • issm/trunk-jpl/src/c/shared/Enum/Enum.vim

    r26121 r26126  
    350350syn keyword cConstant SolidearthSettingsAbstolEnum
    351351syn keyword cConstant SolidearthSettingsCrossSectionShapeEnum
     352syn keyword cConstant SolidearthSettingsOptimEnum
    352353syn keyword cConstant RotationalAngularVelocityEnum
    353354syn keyword cConstant SolidearthSettingsElasticEnum
    354355syn keyword cConstant SealevelchangeGeometryDoneEnum
    355 syn keyword cConstant SealevelGrotm1Enum
    356 syn keyword cConstant SealevelGrotm2Enum
    357 syn keyword cConstant SealevelGrotm3Enum
    358356syn keyword cConstant RotationalEquatorialMoiEnum
    359357syn keyword cConstant TidalLoveHEnum
    … …  
    758756syn keyword cConstant BslcRateEnum
    759757syn keyword cConstant GmtslcEnum
     758syn keyword cConstant SealevelGrotm1Enum
     759syn keyword cConstant SealevelGrotm2Enum
     760syn keyword cConstant SealevelGrotm3Enum
    760761syn keyword cConstant SealevelRSLBarystaticEnum
    761762syn keyword cConstant SealevelRSLRateEnum

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

    r26121 r26126  
    344344        SolidearthSettingsAbstolEnum,
    345345        SolidearthSettingsCrossSectionShapeEnum,
     346        SolidearthSettingsOptimEnum,
    346347        RotationalAngularVelocityEnum,
    347348        SolidearthSettingsElasticEnum,
    348349        SealevelchangeGeometryDoneEnum,
    349         SealevelGrotm1Enum,
    350         SealevelGrotm2Enum,
    351         SealevelGrotm3Enum,
    352350        RotationalEquatorialMoiEnum,
    353351        TidalLoveHEnum,
    … …  
    754752        BslcRateEnum,
    755753        GmtslcEnum,
     754        SealevelGrotm1Enum,
     755        SealevelGrotm2Enum,
     756        SealevelGrotm3Enum,
    756757        SealevelRSLBarystaticEnum,
    757758        SealevelRSLRateEnum,

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

    r26121 r26126  
    352352                case SolidearthSettingsAbstolEnum : return "SolidearthSettingsAbstol";
    353353                case SolidearthSettingsCrossSectionShapeEnum : return "SolidearthSettingsCrossSectionShape";
     354                case SolidearthSettingsOptimEnum : return "SolidearthSettingsOptim";
    354355                case RotationalAngularVelocityEnum : return "RotationalAngularVelocity";
    355356                case SolidearthSettingsElasticEnum : return "SolidearthSettingsElastic";
    356357                case SealevelchangeGeometryDoneEnum : return "SealevelchangeGeometryDone";
    357                 case SealevelGrotm1Enum : return "SealevelGrotm1";
    358                 case SealevelGrotm2Enum : return "SealevelGrotm2";
    359                 case SealevelGrotm3Enum : return "SealevelGrotm3";
    360358                case RotationalEquatorialMoiEnum : return "RotationalEquatorialMoi";
    361359                case TidalLoveHEnum : return "TidalLoveH";
    … …  
    760758                case BslcRateEnum : return "BslcRate";
    761759                case GmtslcEnum : return "Gmtslc";
     760                case SealevelGrotm1Enum : return "SealevelGrotm1";
     761                case SealevelGrotm2Enum : return "SealevelGrotm2";
     762                case SealevelGrotm3Enum : return "SealevelGrotm3";
    762763                case SealevelRSLBarystaticEnum : return "SealevelRSLBarystatic";
    763764                case SealevelRSLRateEnum : return "SealevelRSLRate";

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

    r26121 r26126  
    358358              else if (strcmp(name,"SolidearthSettingsAbstol")==0) return SolidearthSettingsAbstolEnum;
    359359              else if (strcmp(name,"SolidearthSettingsCrossSectionShape")==0) return SolidearthSettingsCrossSectionShapeEnum;
     360              else if (strcmp(name,"SolidearthSettingsOptim")==0) return SolidearthSettingsOptimEnum;
    360361              else if (strcmp(name,"RotationalAngularVelocity")==0) return RotationalAngularVelocityEnum;
    361362              else if (strcmp(name,"SolidearthSettingsElastic")==0) return SolidearthSettingsElasticEnum;
    362363              else if (strcmp(name,"SealevelchangeGeometryDone")==0) return SealevelchangeGeometryDoneEnum;
    363               else if (strcmp(name,"SealevelGrotm1")==0) return SealevelGrotm1Enum;
    364               else if (strcmp(name,"SealevelGrotm2")==0) return SealevelGrotm2Enum;
    365               else if (strcmp(name,"SealevelGrotm3")==0) return SealevelGrotm3Enum;
    366364              else if (strcmp(name,"RotationalEquatorialMoi")==0) return RotationalEquatorialMoiEnum;
    367365              else if (strcmp(name,"TidalLoveH")==0) return TidalLoveHEnum;
    … …  
    383381              else if (strcmp(name,"SolidearthSettingsOceanAreaScaling")==0) return SolidearthSettingsOceanAreaScalingEnum;
    384382              else if (strcmp(name,"RotationalPolarMoi")==0) return RotationalPolarMoiEnum;
     383              else if (strcmp(name,"SolidearthSettingsReltol")==0) return SolidearthSettingsReltolEnum;
     384              else if (strcmp(name,"SealevelchangeRequestedOutputs")==0) return SealevelchangeRequestedOutputsEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"SolidearthSettingsReltol")==0) return SolidearthSettingsReltolEnum;
    389               else if (strcmp(name,"SealevelchangeRequestedOutputs")==0) return SealevelchangeRequestedOutputsEnum;
    390               else if (strcmp(name,"SolidearthSettingsRigid")==0) return SolidearthSettingsRigidEnum;
     388              if (strcmp(name,"SolidearthSettingsRigid")==0) return SolidearthSettingsRigidEnum;
    391389              else if (strcmp(name,"SolidearthSettingsRotation")==0) return SolidearthSettingsRotationEnum;
    392390              else if (strcmp(name,"SealevelchangeRunCount")==0) return SealevelchangeRunCountEnum;
    … …  
    506504              else if (strcmp(name,"TransientIsesa")==0) return TransientIsesaEnum;
    507505              else if (strcmp(name,"TransientIsgia")==0) return TransientIsgiaEnum;
     506              else if (strcmp(name,"TransientIsgroundingline")==0) return TransientIsgroundinglineEnum;
     507              else if (strcmp(name,"TransientIshydrology")==0) return TransientIshydrologyEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"TransientIsgroundingline")==0) return TransientIsgroundinglineEnum;
    512               else if (strcmp(name,"TransientIshydrology")==0) return TransientIshydrologyEnum;
    513               else if (strcmp(name,"TransientIsmasstransport")==0) return TransientIsmasstransportEnum;
     511              if (strcmp(name,"TransientIsmasstransport")==0) return TransientIsmasstransportEnum;
    514512              else if (strcmp(name,"TransientIsoceantransport")==0) return TransientIsoceantransportEnum;
    515513              else if (strcmp(name,"TransientIsmovingfront")==0) return TransientIsmovingfrontEnum;
    … …  
    629627              else if (strcmp(name,"Enthalpy")==0) return EnthalpyEnum;
    630628              else if (strcmp(name,"EnthalpyPicard")==0) return EnthalpyPicardEnum;
     629              else if (strcmp(name,"EplHead")==0) return EplHeadEnum;
     630              else if (strcmp(name,"EplHeadOld")==0) return EplHeadOldEnum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"EplHead")==0) return EplHeadEnum;
    635               else if (strcmp(name,"EplHeadOld")==0) return EplHeadOldEnum;
    636               else if (strcmp(name,"EplHeadSlopeX")==0) return EplHeadSlopeXEnum;
     634              if (strcmp(name,"EplHeadSlopeX")==0) return EplHeadSlopeXEnum;
    637635              else if (strcmp(name,"EplHeadSlopeY")==0) return EplHeadSlopeYEnum;
    638636              else if (strcmp(name,"EplHeadSubstep")==0) return EplHeadSubstepEnum;
    … …  
    752750              else if (strcmp(name,"RadarAttenuationMacGregor")==0) return RadarAttenuationMacGregorEnum;
    753751              else if (strcmp(name,"RadarAttenuationWolff")==0) return RadarAttenuationWolffEnum;
     752              else if (strcmp(name,"RadarIcePeriod")==0) return RadarIcePeriodEnum;
     753              else if (strcmp(name,"RadarPowerMacGregor")==0) return RadarPowerMacGregorEnum;
    754754         else stage=7;
    755755   }
    756756   if(stage==7){
    757               if (strcmp(name,"RadarIcePeriod")==0) return RadarIcePeriodEnum;
    758               else if (strcmp(name,"RadarPowerMacGregor")==0) return RadarPowerMacGregorEnum;
    759               else if (strcmp(name,"RadarPowerWolff")==0) return RadarPowerWolffEnum;
     757              if (strcmp(name,"RadarPowerWolff")==0) return RadarPowerWolffEnum;
    760758              else if (strcmp(name,"RheologyBAbsGradient")==0) return RheologyBAbsGradientEnum;
    761759              else if (strcmp(name,"RheologyBInitialguess")==0) return RheologyBInitialguessEnum;
    … …  
    778776              else if (strcmp(name,"BslcRate")==0) return BslcRateEnum;
    779777              else if (strcmp(name,"Gmtslc")==0) return GmtslcEnum;
     778              else if (strcmp(name,"SealevelGrotm1")==0) return SealevelGrotm1Enum;
     779              else if (strcmp(name,"SealevelGrotm2")==0) return SealevelGrotm2Enum;
     780              else if (strcmp(name,"SealevelGrotm3")==0) return SealevelGrotm3Enum;
    780781              else if (strcmp(name,"SealevelRSLBarystatic")==0) return SealevelRSLBarystaticEnum;
    781782              else if (strcmp(name,"SealevelRSLRate")==0) return SealevelRSLRateEnum;
    … …  
    874875              else if (strcmp(name,"SmbReini")==0) return SmbReiniEnum;
    875876              else if (strcmp(name,"SmbRunoff")==0) return SmbRunoffEnum;
    876               else if (strcmp(name,"SmbRunoffSubstep")==0) return SmbRunoffSubstepEnum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"SmbRunoffTransient")==0) return SmbRunoffTransientEnum;
     880              if (strcmp(name,"SmbRunoffSubstep")==0) return SmbRunoffSubstepEnum;
     881              else if (strcmp(name,"SmbRunoffTransient")==0) return SmbRunoffTransientEnum;
    881882              else if (strcmp(name,"SmbS0gcm")==0) return SmbS0gcmEnum;
    882883              else if (strcmp(name,"SmbS0p")==0) return SmbS0pEnum;
    … …  
    997998              else if (strcmp(name,"Outputdefinition24")==0) return Outputdefinition24Enum;
    998999              else if (strcmp(name,"Outputdefinition25")==0) return Outputdefinition25Enum;
    999               else if (strcmp(name,"Outputdefinition26")==0) return Outputdefinition26Enum;
    10001000         else stage=9;
    10011001   }
    10021002   if(stage==9){
    1003               if (strcmp(name,"Outputdefinition27")==0) return Outputdefinition27Enum;
     1003              if (strcmp(name,"Outputdefinition26")==0) return Outputdefinition26Enum;
     1004              else if (strcmp(name,"Outputdefinition27")==0) return Outputdefinition27Enum;
    10041005              else if (strcmp(name,"Outputdefinition28")==0) return Outputdefinition28Enum;
    10051006              else if (strcmp(name,"Outputdefinition29")==0) return Outputdefinition29Enum;
    … …  
    11201121              else if (strcmp(name,"CalvingVonmises")==0) return CalvingVonmisesEnum;
    11211122              else if (strcmp(name,"Cfdragcoeffabsgrad")==0) return CfdragcoeffabsgradEnum;
    1122               else if (strcmp(name,"Cfsurfacelogvel")==0) return CfsurfacelogvelEnum;
    11231123         else stage=10;
    11241124   }
    11251125   if(stage==10){
    1126               if (strcmp(name,"Cfsurfacesquare")==0) return CfsurfacesquareEnum;
     1126              if (strcmp(name,"Cfsurfacelogvel")==0) return CfsurfacelogvelEnum;
     1127              else if (strcmp(name,"Cfsurfacesquare")==0) return CfsurfacesquareEnum;
    11271128              else if (strcmp(name,"Channel")==0) return ChannelEnum;
    11281129              else if (strcmp(name,"ChannelArea")==0) return ChannelAreaEnum;
    … …  
    12431244              else if (strcmp(name,"Inputs")==0) return InputsEnum;
    12441245              else if (strcmp(name,"Internal")==0) return InternalEnum;
    1245               else if (strcmp(name,"Intersect")==0) return IntersectEnum;
    12461246         else stage=11;
    12471247   }
    12481248   if(stage==11){
    1249               if (strcmp(name,"InversionVzObs")==0) return InversionVzObsEnum;
     1249              if (strcmp(name,"Intersect")==0) return IntersectEnum;
     1250              else if (strcmp(name,"InversionVzObs")==0) return InversionVzObsEnum;
    12501251              else if (strcmp(name,"J")==0) return JEnum;
    12511252              else if (strcmp(name,"L1L2Approximation")==0) return L1L2ApproximationEnum;
    … …  
    13661367              else if (strcmp(name,"SMBgcm")==0) return SMBgcmEnum;
    13671368              else if (strcmp(name,"SMBgemb")==0) return SMBgembEnum;
    1368               else if (strcmp(name,"SMBgradients")==0) return SMBgradientsEnum;
    13691369         else stage=12;
    13701370   }
    13711371   if(stage==12){
    1372               if (strcmp(name,"SMBgradientscomponents")==0) return SMBgradientscomponentsEnum;
     1372              if (strcmp(name,"SMBgradients")==0) return SMBgradientsEnum;
     1373              else if (strcmp(name,"SMBgradientscomponents")==0) return SMBgradientscomponentsEnum;
    13731374              else if (strcmp(name,"SMBgradientsela")==0) return SMBgradientselaEnum;
    13741375              else if (strcmp(name,"SMBhenning")==0) return SMBhenningEnum;

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

    r26060 r26126  
    2222                grdmodel               = 0; %grd model (0 by default, 1 for elastic, 2 for Ivins)
    2323                cross_section_shape    = 0; %cross section only used when grd model is Ivins
     24                optim                  = 0; %new optimized version of the GRD code.
    2425        end
    2526        methods
    … …  
    6566                self.cross_section_shape=1; %square as default (see iedge in GiaDeflectionCorex)
    6667
     68                %optim?
     69                self.optim=0;
     70
    6771                %no grd model by default:
    6872                self.grdmodel=0;
    … …  
    8387                        md = checkfield(md,'fieldname','solidearth.settings.grdmodel','values',[1 2]);
    8488                        md = checkfield(md,'fieldname','solidearth.settings.cross_section_shape','numel',[1],'values',[1,2]);
     89                        md = checkfield(md,'fieldname','solidearth.settings.optim','values',[0,1]);
    8590
    8691                        %checks on computational flags
    … …  
    126131                        fielddisplay(self,'grdmodel','type of deformation model, 1 for elastic, 2 for visco-elastic from Ivins');
    127132                        fielddisplay(self,'cross_section_shape','1: square-edged (default). 2: elliptical. See iedge in GiaDeflectionCore');
     133                        fielddisplay(self,'optim','use optimized version of the GRD code? (default 0)');
    128134                end % }}}
    129135                function marshall(self,prefix,md,fid) % {{{
    … …  
    144150                        WriteData(fid,prefix,'object',self,'fieldname','grdmodel','name','md.solidearth.settings.grdmodel','format','Integer');
    145151                        WriteData(fid,prefix,'object',self,'fieldname','cross_section_shape','name','md.solidearth.settings.cross_section_shape','format','Integer');
     152                        WriteData(fid,prefix,'object',self,'fieldname','optim','name','md.solidearth.settings.optim','format','Integer');
    146153                end % }}}
    147154                function savemodeljs(self,fid,modelname) % {{{
    … …  
    158165                        writejsdouble(fid,[modelname '.solidearth.settings.glfraction'],self.glfraction);
    159166                        writejsdouble(fid,[modelname '.solidearth.settings.cross_section_shape'],self.cross_section_shape);
     167                        writejsdouble(fid,[modelname '.solidearth.settings.optim'],self.optim);
    160168                end % }}}
    161169                function self = extrude(self,md) % {{{
Note: See TracChangeset for help on using the changeset viewer.