Changeset 21773


Ignore:
Timestamp:
06/30/17 09:33:47 (8 years ago)
Author:
Mathieu Morlighem
Message:

CHG: added class to handle constrained adaptive time stepping

Location:
issm/trunk-jpl/src/c
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/classes/FemModel.cpp

    r21762 r21773  
    20712071
    20722072        /*intermediary: */
    2073         Element *element     = NULL;
    20742073        IssmDouble   min_dt      = 0;
    20752074        IssmDouble   node_min_dt = 0;
    20762075
    20772076        /*Go through elements, and figure out the minimum of the time steps for each element (using CFL criterion): */
    2078         element=(Element*)elements->GetObjectByOffset(0); min_dt=element->TimeAdapt();
     2077        Element* element=(Element*)elements->GetObjectByOffset(0); min_dt=element->TimeAdapt();
    20792078
    20802079        for (i=1;i<elements->Size();i++){
     
    20882087        ISSM_MPI_Bcast(&node_min_dt,1,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
    20892088        min_dt=node_min_dt;
     2089
     2090        /*Constrain dt if necessary*/
     2091        int timestepping_type;
     2092        IssmDouble dt_low,dt_high;
     2093        parameters->FindParam(&timestepping_type,TimesteppingTypeEnum);
     2094        if(timestepping_type==TimeAdaptationEnum){
     2095                parameters->FindParam(&dt_low,TimesteppingTimeStepMinEnum);
     2096                parameters->FindParam(&dt_high,TimesteppingTimeStepMaxEnum);
     2097                if(min_dt<dt_low)  min_dt = dt_low;
     2098                if(min_dt>dt_high) min_dt = dt_high;
     2099        }
    20902100
    20912101        /*Assign output pointers:*/
  • issm/trunk-jpl/src/c/classes/IoModel.cpp

    r21750 r21773  
    603603                                        if(strcmp(record_name,"md.hydrology.model")==0) integer = IoCodeToEnumHydrology(integer);
    604604                                        if(strcmp(record_name,"md.materials.type")==0) integer = IoCodeToEnumMaterials(integer);
     605                                        if(strcmp(record_name,"md.timestepping.type")==0) integer = IoCodeToEnumTimestepping(integer);
    605606
    606607                                        /*Broadcast to other cpus*/
  • issm/trunk-jpl/src/c/cores/transient_core.cpp

    r21680 r21773  
    4444        femmodel->parameters->FindParam(&time,TimeEnum);
    4545        femmodel->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
    46         femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum);
    4746        femmodel->parameters->FindParam(&yts,ConstantsYtsEnum);
    4847        femmodel->parameters->FindParam(&dakota_analysis,QmuIsdakotaEnum);
     
    7776                        if(time+dt>finaltime) dt=finaltime-time;
    7877                        femmodel->parameters->SetParam(dt,TimesteppingTimeStepEnum);
     78                }
     79                else{
     80                        femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum);
    7981                }
    8082                step+=1;
  • issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp

    r21749 r21773  
    1919
    2020        int         i,j,m,k;
    21         int         numoutputs,materialtype,smb_model,basalforcing_model;
     21        int         numoutputs,materialtype,smb_model,basalforcing_model,timestepping_type;
    2222        char**      requestedoutputs = NULL;
    2323        IssmDouble  time;
     
    4646        parameters->AddObject(iomodel->CopyConstantObject("md.settings.recording_frequency",SettingsRecordingFrequencyEnum));
    4747        parameters->AddObject(iomodel->CopyConstantObject("md.constants.yts",ConstantsYtsEnum));
    48         parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.start_time",TimesteppingStartTimeEnum));
    49         parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.final_time",TimesteppingFinalTimeEnum));
    50         parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.time_adapt",TimesteppingTimeAdaptEnum));
    51         parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.time_step",TimesteppingTimeStepEnum));
    52         parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.cfl_coefficient",TimesteppingCflCoefficientEnum));
    53         parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.interp_forcings",TimesteppingInterpForcingsEnum));
    5448        parameters->AddObject(iomodel->CopyConstantObject("md.settings.lowmem",SettingsLowmemEnum));
    5549        parameters->AddObject(iomodel->CopyConstantObject("md.debug.profiling",DebugProfilingEnum));
     
    142136        parameters->AddObject(new IntParam(SolutionTypeEnum,solution_type));
    143137
     138        /*Time stepping*/
     139        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.type",TimesteppingTypeEnum));
     140        iomodel->FindConstant(&timestepping_type,"md.timestepping.type");
     141        switch(timestepping_type){
     142                case FixedTimesteppingEnum:
     143                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.start_time",TimesteppingStartTimeEnum));
     144                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.final_time",TimesteppingFinalTimeEnum));
     145                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.time_adapt",TimesteppingTimeAdaptEnum));
     146                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.time_step",TimesteppingTimeStepEnum));
     147                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.cfl_coefficient",TimesteppingCflCoefficientEnum));
     148                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.interp_forcings",TimesteppingInterpForcingsEnum));
     149                        break;
     150                case TimeAdaptationEnum:
     151                        parameters->AddObject(new BoolParam(TimesteppingTimeAdaptEnum,true));
     152                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.start_time",TimesteppingStartTimeEnum));
     153                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.final_time",TimesteppingFinalTimeEnum));
     154                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.time_step_min",TimesteppingTimeStepMinEnum));
     155                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.time_step_max",TimesteppingTimeStepMaxEnum));
     156                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.cfl_coefficient",TimesteppingCflCoefficientEnum));
     157                        parameters->AddObject(iomodel->CopyConstantObject("md.timestepping.interp_forcings",TimesteppingInterpForcingsEnum));
     158                        break;
     159                default:
     160                        _error_(EnumToStringx(timestepping_type) <<" not supported yet");
     161        }
    144162        iomodel->FindConstant(&time,"md.timestepping.start_time");
    145163        parameters->AddObject(new DoubleParam(TimeEnum,time)); 
  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r21752 r21773  
    3030        SmoothThicknessMultiplierEnum,
    3131        LevelsetStabilizationEnum,
     32        TimesteppingTypeEnum,
     33        FixedTimesteppingEnum,
     34        TimeAdaptationEnum,
     35        TimesteppingTimeStepMinEnum,
     36        TimesteppingTimeStepMaxEnum,
    3237        /*}}}*/
    3338        /*Model fields {{{*/
  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r21752 r21773  
    3838                case SmoothThicknessMultiplierEnum : return "SmoothThicknessMultiplier";
    3939                case LevelsetStabilizationEnum : return "LevelsetStabilization";
     40                case TimesteppingTypeEnum : return "TimesteppingType";
     41                case FixedTimesteppingEnum : return "FixedTimestepping";
     42                case TimeAdaptationEnum : return "TimeAdaptation";
     43                case TimesteppingTimeStepMinEnum : return "TimesteppingTimeStepMin";
     44                case TimesteppingTimeStepMaxEnum : return "TimesteppingTimeStepMax";
    4045                case AutodiffIsautodiffEnum : return "AutodiffIsautodiff";
    4146                case AutodiffNumDependentsEnum : return "AutodiffNumDependents";
  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r21752 r21773  
    3838              else if (strcmp(name,"SmoothThicknessMultiplier")==0) return SmoothThicknessMultiplierEnum;
    3939              else if (strcmp(name,"LevelsetStabilization")==0) return LevelsetStabilizationEnum;
     40              else if (strcmp(name,"TimesteppingType")==0) return TimesteppingTypeEnum;
     41              else if (strcmp(name,"FixedTimestepping")==0) return FixedTimesteppingEnum;
     42              else if (strcmp(name,"TimeAdaptation")==0) return TimeAdaptationEnum;
     43              else if (strcmp(name,"TimesteppingTimeStepMin")==0) return TimesteppingTimeStepMinEnum;
     44              else if (strcmp(name,"TimesteppingTimeStepMax")==0) return TimesteppingTimeStepMaxEnum;
    4045              else if (strcmp(name,"AutodiffIsautodiff")==0) return AutodiffIsautodiffEnum;
    4146              else if (strcmp(name,"AutodiffNumDependents")==0) return AutodiffNumDependentsEnum;
     
    132137              else if (strcmp(name,"Hydrologydc")==0) return HydrologydcEnum;
    133138              else if (strcmp(name,"SedimentHead")==0) return SedimentHeadEnum;
    134               else if (strcmp(name,"SedimentHeadOld")==0) return SedimentHeadOldEnum;
     139         else stage=2;
     140   }
     141   if(stage==2){
     142              if (strcmp(name,"SedimentHeadOld")==0) return SedimentHeadOldEnum;
    135143              else if (strcmp(name,"SedimentHeadResidual")==0) return SedimentHeadResidualEnum;
    136144              else if (strcmp(name,"EffectivePressure")==0) return EffectivePressureEnum;
    137145              else if (strcmp(name,"EplHead")==0) return EplHeadEnum;
    138146              else if (strcmp(name,"EplHeadOld")==0) return EplHeadOldEnum;
    139          else stage=2;
    140    }
    141    if(stage==2){
    142               if (strcmp(name,"EplHeadSlopeX")==0) return EplHeadSlopeXEnum;
     147              else if (strcmp(name,"EplHeadSlopeX")==0) return EplHeadSlopeXEnum;
    143148              else if (strcmp(name,"EplHeadSlopeY")==0) return EplHeadSlopeYEnum;
    144149              else if (strcmp(name,"EplZigZagCounter")==0) return EplZigZagCounterEnum;
     
    255260              else if (strcmp(name,"StressIntensityFactor")==0) return StressIntensityFactorEnum;
    256261              else if (strcmp(name,"CalvingLaw")==0) return CalvingLawEnum;
    257               else if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
     262         else stage=3;
     263   }
     264   if(stage==3){
     265              if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
    258266              else if (strcmp(name,"CalvingMeltingrate")==0) return CalvingMeltingrateEnum;
    259267              else if (strcmp(name,"CalvingLevermann")==0) return CalvingLevermannEnum;
    260268              else if (strcmp(name,"CalvingDev")==0) return CalvingDevEnum;
    261269              else if (strcmp(name,"CalvingMinthickness")==0) return CalvingMinthicknessEnum;
    262          else stage=3;
    263    }
    264    if(stage==3){
    265               if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum;
     270              else if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum;
    266271              else if (strcmp(name,"CalvinglevermannCoeff")==0) return CalvinglevermannCoeffEnum;
    267272              else if (strcmp(name,"CalvinglevermannMeltingrate")==0) return CalvinglevermannMeltingrateEnum;
     
    378383              else if (strcmp(name,"BalancethicknessD0")==0) return BalancethicknessD0Enum;
    379384              else if (strcmp(name,"Smb")==0) return SmbEnum;
    380               else if (strcmp(name,"SmbAnalysis")==0) return SmbAnalysisEnum;
     385         else stage=4;
     386   }
     387   if(stage==4){
     388              if (strcmp(name,"SmbAnalysis")==0) return SmbAnalysisEnum;
    381389              else if (strcmp(name,"SmbSolution")==0) return SmbSolutionEnum;
    382390              else if (strcmp(name,"SmbNumRequestedOutputs")==0) return SmbNumRequestedOutputsEnum;
    383391              else if (strcmp(name,"SmbRequestedOutputs")==0) return SmbRequestedOutputsEnum;
    384392              else if (strcmp(name,"SmbIsInitialized")==0) return SmbIsInitializedEnum;
    385          else stage=4;
    386    }
    387    if(stage==4){
    388               if (strcmp(name,"SmbDzini")==0) return SmbDziniEnum;
     393              else if (strcmp(name,"SmbDzini")==0) return SmbDziniEnum;
    389394              else if (strcmp(name,"SmbDini")==0) return SmbDiniEnum;
    390395              else if (strcmp(name,"SmbReini")==0) return SmbReiniEnum;
     
    501506              else if (strcmp(name,"MeltingOffset")==0) return MeltingOffsetEnum;
    502507              else if (strcmp(name,"Misfit")==0) return MisfitEnum;
    503               else if (strcmp(name,"Pressure")==0) return PressureEnum;
     508         else stage=5;
     509   }
     510   if(stage==5){
     511              if (strcmp(name,"Pressure")==0) return PressureEnum;
    504512              else if (strcmp(name,"PressurePicard")==0) return PressurePicardEnum;
    505513              else if (strcmp(name,"AndroidFrictionCoefficient")==0) return AndroidFrictionCoefficientEnum;
    506514              else if (strcmp(name,"ResetPenalties")==0) return ResetPenaltiesEnum;
    507515              else if (strcmp(name,"SurfaceAbsVelMisfit")==0) return SurfaceAbsVelMisfitEnum;
    508          else stage=5;
    509    }
    510    if(stage==5){
    511               if (strcmp(name,"SurfaceArea")==0) return SurfaceAreaEnum;
     516              else if (strcmp(name,"SurfaceArea")==0) return SurfaceAreaEnum;
    512517              else if (strcmp(name,"SurfaceAverageVelMisfit")==0) return SurfaceAverageVelMisfitEnum;
    513518              else if (strcmp(name,"SurfaceLogVelMisfit")==0) return SurfaceLogVelMisfitEnum;
     
    624629              else if (strcmp(name,"Outputdefinition22")==0) return Outputdefinition22Enum;
    625630              else if (strcmp(name,"Outputdefinition23")==0) return Outputdefinition23Enum;
    626               else if (strcmp(name,"Outputdefinition24")==0) return Outputdefinition24Enum;
     631         else stage=6;
     632   }
     633   if(stage==6){
     634              if (strcmp(name,"Outputdefinition24")==0) return Outputdefinition24Enum;
    627635              else if (strcmp(name,"Outputdefinition25")==0) return Outputdefinition25Enum;
    628636              else if (strcmp(name,"Outputdefinition26")==0) return Outputdefinition26Enum;
    629637              else if (strcmp(name,"Outputdefinition27")==0) return Outputdefinition27Enum;
    630638              else if (strcmp(name,"Outputdefinition28")==0) return Outputdefinition28Enum;
    631          else stage=6;
    632    }
    633    if(stage==6){
    634               if (strcmp(name,"Outputdefinition29")==0) return Outputdefinition29Enum;
     639              else if (strcmp(name,"Outputdefinition29")==0) return Outputdefinition29Enum;
    635640              else if (strcmp(name,"Outputdefinition30")==0) return Outputdefinition30Enum;
    636641              else if (strcmp(name,"Outputdefinition31")==0) return Outputdefinition31Enum;
     
    747752              else if (strcmp(name,"InputFileName")==0) return InputFileNameEnum;
    748753              else if (strcmp(name,"LockFileName")==0) return LockFileNameEnum;
    749               else if (strcmp(name,"RestartFileName")==0) return RestartFileNameEnum;
     754         else stage=7;
     755   }
     756   if(stage==7){
     757              if (strcmp(name,"RestartFileName")==0) return RestartFileNameEnum;
    750758              else if (strcmp(name,"ToolkitsOptionsAnalyses")==0) return ToolkitsOptionsAnalysesEnum;
    751759              else if (strcmp(name,"ToolkitsOptionsStrings")==0) return ToolkitsOptionsStringsEnum;
    752760              else if (strcmp(name,"QmuErrName")==0) return QmuErrNameEnum;
    753761              else if (strcmp(name,"QmuInName")==0) return QmuInNameEnum;
    754          else stage=7;
    755    }
    756    if(stage==7){
    757               if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
     762              else if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
    758763              else if (strcmp(name,"Regular")==0) return RegularEnum;
    759764              else if (strcmp(name,"Scaled")==0) return ScaledEnum;
     
    870875              else if (strcmp(name,"ElementSId")==0) return ElementSIdEnum;
    871876              else if (strcmp(name,"VectorParam")==0) return VectorParamEnum;
    872               else if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
     877         else stage=8;
     878   }
     879   if(stage==8){
     880              if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
    873881              else if (strcmp(name,"Segment")==0) return SegmentEnum;
    874882              else if (strcmp(name,"SegmentRiftfront")==0) return SegmentRiftfrontEnum;
    875883              else if (strcmp(name,"SpcDynamic")==0) return SpcDynamicEnum;
    876884              else if (strcmp(name,"SpcStatic")==0) return SpcStaticEnum;
    877          else stage=8;
    878    }
    879    if(stage==8){
    880               if (strcmp(name,"SpcTransient")==0) return SpcTransientEnum;
     885              else if (strcmp(name,"SpcTransient")==0) return SpcTransientEnum;
    881886              else if (strcmp(name,"StringArrayParam")==0) return StringArrayParamEnum;
    882887              else if (strcmp(name,"StringParam")==0) return StringParamEnum;
     
    993998              else if (strcmp(name,"P2")==0) return P2Enum;
    994999              else if (strcmp(name,"P2bubble")==0) return P2bubbleEnum;
    995               else if (strcmp(name,"P2bubblecondensed")==0) return P2bubblecondensedEnum;
     1000         else stage=9;
     1001   }
     1002   if(stage==9){
     1003              if (strcmp(name,"P2bubblecondensed")==0) return P2bubblecondensedEnum;
    9961004              else if (strcmp(name,"P2xP1")==0) return P2xP1Enum;
    9971005              else if (strcmp(name,"P1xP2")==0) return P1xP2Enum;
    9981006              else if (strcmp(name,"P1xP3")==0) return P1xP3Enum;
    9991007              else if (strcmp(name,"P2xP4")==0) return P2xP4Enum;
    1000          else stage=9;
    1001    }
    1002    if(stage==9){
    1003               if (strcmp(name,"P1P1")==0) return P1P1Enum;
     1008              else if (strcmp(name,"P1P1")==0) return P1P1Enum;
    10041009              else if (strcmp(name,"P1P1GLS")==0) return P1P1GLSEnum;
    10051010              else if (strcmp(name,"MINI")==0) return MINIEnum;
  • issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.cpp

    r21654 r21773  
    121121        }
    122122}/*}}}*/
     123int IoCodeToEnumTimestepping(int enum_in){/*{{{*/
     124        switch(enum_in){
     125                case 1: return FixedTimesteppingEnum;
     126                case 2: return TimeAdaptationEnum;
     127                default: _error_("Marshalled materials code \""<<enum_in<<"\" not supported yet");
     128        }
     129}/*}}}*/
    123130
    124131int IoCodeToEnumVertexEquation(int enum_in){/*{{{*/
  • issm/trunk-jpl/src/c/shared/io/Marshalling/IoCodeConversions.h

    r21049 r21773  
    99int IoCodeToEnumHydrology(int enum_in);
    1010int IoCodeToEnumMaterials(int enum_in);
     11int IoCodeToEnumTimestepping(int enum_in);
    1112
    1213int IoCodeToEnumVertexEquation(int enum_in);
Note: See TracChangeset for help on using the changeset viewer.