Changeset 26750


Ignore:
Timestamp:
12/27/21 08:23:23 (3 years ago)
Author:
vverjans
Message:

NEW: added SMBforcing to stochasticforcingfields

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

Legend:

Unmodified
Added
Removed

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

    r26550 r26750  
    2424
    2525        int    smb_model;
    26         bool   isdelta18o,ismungsm,isd18opd,issetpddfac,isprecipscaled,istemperaturescaled,isfirnwarming;
     26        bool   isdelta18o,ismungsm,isd18opd,issetpddfac,isprecipscaled,istemperaturescaled,isfirnwarming,isstochastic;
    2727
    2828        /*Update elements: */
    … …  
    3838        /*Figure out smb model: */
    3939        iomodel->FindConstant(&smb_model,"md.smb.model");
     40        iomodel->FindConstant(&isstochastic,"md.stochasticforcing.isstochasticforcing");
    4041        switch(smb_model){
    4142                case SMBforcingEnum:
    4243                        iomodel->FetchDataToInput(inputs,elements,"md.smb.mass_balance",SmbMassBalanceEnum,0.);
     44                        if(isstochastic){
     45                                iomodel->FetchDataToInput(inputs,elements,"md.stochasticforcing.default_id",StochasticForcingDefaultIdEnum);
     46                                iomodel->FetchDataToInput(inputs,elements,"md.smb.mass_balance",BaselineSmbMassBalanceEnum,0.);
     47                        }
    4348                        break;
    4449                case SMBgembEnum:

  • issm/trunk-jpl/src/c/modules/StochasticForcingx/StochasticForcingx.cpp

    r26676 r26750  
    5858      multivariateNormal(&temparray,dimtot,0.0,covariance,fixedseed);
    5959      for(int i=0;i<dimtot;i++) noiseterms[i]=temparray[i];
    60       xDelete<IssmDouble>(temparray);
     60                xDelete<IssmDouble>(temparray);
    6161   }
    6262   ISSM_MPI_Bcast(noiseterms,dimtot,ISSM_MPI_DOUBLE,0,IssmComm::GetComm());
    … …  
    142142                                        }
    143143                                        break;
     144                                case SMBforcingEnum:
     145                                        /*Delete SmbMassBalanceEnum at previous time step (required if it is transient)*/
     146                                        femmodel->inputs->DeleteInput(SmbMassBalanceEnum);
     147                                        for(Object* &object:femmodel->elements->objects){
     148                                                Element* element = xDynamicCast<Element*>(object);
     149                                                int numvertices  = element->GetNumberOfVertices();
     150                                                IssmDouble baselinesmb;
     151                                                IssmDouble smb_tot[numvertices];
     152                                                Input* baselinesmb_input  = NULL;
     153                                                baselinesmb_input = element->GetInput(BaselineSmbMassBalanceEnum); _assert_(baselinesmb_input);
     154                                                element->GetInputValue(&dimensionid,StochasticForcingDefaultIdEnum);
     155                                                Gauss* gauss = element->NewGauss();
     156                                                for(int i=0;i<numvertices;i++){
     157                                                        gauss->GaussVertex(i);
     158                                                        baselinesmb_input->GetInputValue(&baselinesmb,gauss);
     159                                                        smb_tot[i] = baselinesmb+noisefield[dimensionid];
     160                                                }
     161                                                element->AddInput(SmbMassBalanceEnum,&smb_tot[0],P1DGEnum);
     162                                                delete gauss;
     163                                        }
     164                                        break;
    144165                                //VV working(29Nov2021)
    145166                                case FrictionWaterPressureEnum:

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

    r26747 r26750  
    604604syn keyword cConstant BaselineCalvingCalvingrateEnum
    605605syn keyword cConstant BaselineFrictionEffectivePressureEnum
     606syn keyword cConstant BaselineSmbMassBalanceEnum
    606607syn keyword cConstant BedEnum
    607608syn keyword cConstant BedGRDEnum
    … …  
    15821583syn keyword cType Cfsurfacesquare
    15831584syn keyword cType Channel
     1585syn keyword cType classes
    15841586syn keyword cType Constraint
    15851587syn keyword cType Constraints
    … …  
    15881590syn keyword cType ControlInput
    15891591syn keyword cType Covertree
     1592syn keyword cType DatasetInput
    15901593syn keyword cType DataSetParam
    1591 syn keyword cType DatasetInput
    15921594syn keyword cType Definition
    15931595syn keyword cType DependentObject
    … …  
    16021604syn keyword cType ElementInput
    16031605syn keyword cType ElementMatrix
     1606syn keyword cType Elements
    16041607syn keyword cType ElementVector
    1605 syn keyword cType Elements
    16061608syn keyword cType ExponentialVariogram
    16071609syn keyword cType ExternalResult
    … …  
    16101612syn keyword cType Friction
    16111613syn keyword cType Gauss
     1614syn keyword cType GaussianVariogram
     1615syn keyword cType gaussobjects
    16121616syn keyword cType GaussPenta
    16131617syn keyword cType GaussSeg
    16141618syn keyword cType GaussTetra
    16151619syn keyword cType GaussTria
    1616 syn keyword cType GaussianVariogram
    16171620syn keyword cType GenericExternalResult
    16181621syn keyword cType GenericOption
    … …  
    16301633syn keyword cType IssmDirectApplicInterface
    16311634syn keyword cType IssmParallelDirectApplicInterface
     1635syn keyword cType krigingobjects
    16321636syn keyword cType Load
    16331637syn keyword cType Loads
    … …  
    16401644syn keyword cType Matice
    16411645syn keyword cType Matlitho
     1646syn keyword cType matrixobjects
    16421647syn keyword cType MatrixParam
    16431648syn keyword cType Misfit
    … …  
    16521657syn keyword cType Observations
    16531658syn keyword cType Option
     1659syn keyword cType Options
    16541660syn keyword cType OptionUtilities
    1655 syn keyword cType Options
    16561661syn keyword cType Param
    16571662syn keyword cType Parameters
    … …  
    16671672syn keyword cType Regionaloutput
    16681673syn keyword cType Results
     1674syn keyword cType Riftfront
    16691675syn keyword cType RiftStruct
    1670 syn keyword cType Riftfront
    16711676syn keyword cType SealevelGeometry
    16721677syn keyword cType Seg
    16731678syn keyword cType SegInput
     1679syn keyword cType Segment
    16741680syn keyword cType SegRef
    1675 syn keyword cType Segment
    16761681syn keyword cType SpcDynamic
    16771682syn keyword cType SpcStatic
    … …  
    16921697syn keyword cType Vertex
    16931698syn keyword cType Vertices
    1694 syn keyword cType classes
    1695 syn keyword cType gaussobjects
    1696 syn keyword cType krigingobjects
    1697 syn keyword cType matrixobjects
    16981699syn keyword cType AdjointBalancethickness2Analysis
    16991700syn keyword cType AdjointBalancethicknessAnalysis

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

    r26747 r26750  
    600600        BaselineCalvingCalvingrateEnum,
    601601        BaselineFrictionEffectivePressureEnum,
     602        BaselineSmbMassBalanceEnum,
    602603        BedEnum,
    603604        BedGRDEnum,

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

    r26747 r26750  
    606606                case BaselineCalvingCalvingrateEnum : return "BaselineCalvingCalvingrate";
    607607                case BaselineFrictionEffectivePressureEnum : return "BaselineFrictionEffectivePressure";
     608                case BaselineSmbMassBalanceEnum : return "BaselineSmbMassBalance";
    608609                case BedEnum : return "Bed";
    609610                case BedGRDEnum : return "BedGRD";

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

    r26747 r26750  
    597597syn keyword juliaConstC BaselineCalvingCalvingrateEnum
    598598syn keyword juliaConstC BaselineFrictionEffectivePressureEnum
     599syn keyword juliaConstC BaselineSmbMassBalanceEnum
    599600syn keyword juliaConstC BedEnum
    600601syn keyword juliaConstC BedGRDEnum

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

    r26747 r26750  
    618618              else if (strcmp(name,"BaselineCalvingCalvingrate")==0) return BaselineCalvingCalvingrateEnum;
    619619              else if (strcmp(name,"BaselineFrictionEffectivePressure")==0) return BaselineFrictionEffectivePressureEnum;
     620              else if (strcmp(name,"BaselineSmbMassBalance")==0) return BaselineSmbMassBalanceEnum;
    620621              else if (strcmp(name,"Bed")==0) return BedEnum;
    621622              else if (strcmp(name,"BedGRD")==0) return BedGRDEnum;
    … …  
    628629              else if (strcmp(name,"BottomPressure")==0) return BottomPressureEnum;
    629630              else if (strcmp(name,"BottomPressureOld")==0) return BottomPressureOldEnum;
    630               else if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"CalvingHabFraction")==0) return CalvingHabFractionEnum;
     634              if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
     635              else if (strcmp(name,"CalvingHabFraction")==0) return CalvingHabFractionEnum;
    635636              else if (strcmp(name,"CalvingMeltingrate")==0) return CalvingMeltingrateEnum;
    636637              else if (strcmp(name,"CalvingStressThresholdFloatingice")==0) return CalvingStressThresholdFloatingiceEnum;
    … …  
    751752              else if (strcmp(name,"HydrologyMoulinInput")==0) return HydrologyMoulinInputEnum;
    752753              else if (strcmp(name,"HydrologyNeumannflux")==0) return HydrologyNeumannfluxEnum;
    753               else if (strcmp(name,"HydrologyReynolds")==0) return HydrologyReynoldsEnum;
    754754         else stage=7;
    755755   }
    756756   if(stage==7){
    757               if (strcmp(name,"HydrologySheetConductivity")==0) return HydrologySheetConductivityEnum;
     757              if (strcmp(name,"HydrologyReynolds")==0) return HydrologyReynoldsEnum;
     758              else if (strcmp(name,"HydrologySheetConductivity")==0) return HydrologySheetConductivityEnum;
    758759              else if (strcmp(name,"HydrologySheetThickness")==0) return HydrologySheetThicknessEnum;
    759760              else if (strcmp(name,"HydrologySheetThicknessOld")==0) return HydrologySheetThicknessOldEnum;
    … …  
    874875              else if (strcmp(name,"SealevelGErotm3")==0) return SealevelGErotm3Enum;
    875876              else if (strcmp(name,"SealevelRSLBarystatic")==0) return SealevelRSLBarystaticEnum;
    876               else if (strcmp(name,"SealevelRSLRate")==0) return SealevelRSLRateEnum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"SealevelUGrd")==0) return SealevelUGrdEnum;
     880              if (strcmp(name,"SealevelRSLRate")==0) return SealevelRSLRateEnum;
     881              else if (strcmp(name,"SealevelUGrd")==0) return SealevelUGrdEnum;
    881882              else if (strcmp(name,"SealevelNGrd")==0) return SealevelNGrdEnum;
    882883              else if (strcmp(name,"SealevelUEastEsa")==0) return SealevelUEastEsaEnum;
    … …  
    997998              else if (strcmp(name,"SmbS0gcm")==0) return SmbS0gcmEnum;
    998999              else if (strcmp(name,"SmbS0p")==0) return SmbS0pEnum;
    999               else if (strcmp(name,"SmbS0t")==0) return SmbS0tEnum;
    10001000         else stage=9;
    10011001   }
    10021002   if(stage==9){
    1003               if (strcmp(name,"SmbSizeini")==0) return SmbSizeiniEnum;
     1003              if (strcmp(name,"SmbS0t")==0) return SmbS0tEnum;
     1004              else if (strcmp(name,"SmbSizeini")==0) return SmbSizeiniEnum;
    10041005              else if (strcmp(name,"SmbSmbCorr")==0) return SmbSmbCorrEnum;
    10051006              else if (strcmp(name,"SmbSmbref")==0) return SmbSmbrefEnum;
    … …  
    11201121              else if (strcmp(name,"Outputdefinition16")==0) return Outputdefinition16Enum;
    11211122              else if (strcmp(name,"Outputdefinition17")==0) return Outputdefinition17Enum;
    1122               else if (strcmp(name,"Outputdefinition18")==0) return Outputdefinition18Enum;
    11231123         else stage=10;
    11241124   }
    11251125   if(stage==10){
    1126               if (strcmp(name,"Outputdefinition19")==0) return Outputdefinition19Enum;
     1126              if (strcmp(name,"Outputdefinition18")==0) return Outputdefinition18Enum;
     1127              else if (strcmp(name,"Outputdefinition19")==0) return Outputdefinition19Enum;
    11271128              else if (strcmp(name,"Outputdefinition20")==0) return Outputdefinition20Enum;
    11281129              else if (strcmp(name,"Outputdefinition21")==0) return Outputdefinition21Enum;
    … …  
    12431244              else if (strcmp(name,"IntInput")==0) return IntInputEnum;
    12441245              else if (strcmp(name,"DoubleInput")==0) return DoubleInputEnum;
    1245               else if (strcmp(name,"BoolParam")==0) return BoolParamEnum;
    12461246         else stage=11;
    12471247   }
    12481248   if(stage==11){
    1249               if (strcmp(name,"Boundary")==0) return BoundaryEnum;
     1249              if (strcmp(name,"BoolParam")==0) return BoolParamEnum;
     1250              else if (strcmp(name,"Boundary")==0) return BoundaryEnum;
    12501251              else if (strcmp(name,"BuddJacka")==0) return BuddJackaEnum;
    12511252              else if (strcmp(name,"CalvingDev2")==0) return CalvingDev2Enum;
    … …  
    13661367              else if (strcmp(name,"IcefrontMassFlux")==0) return IcefrontMassFluxEnum;
    13671368              else if (strcmp(name,"IcefrontMassFluxLevelset")==0) return IcefrontMassFluxLevelsetEnum;
    1368               else if (strcmp(name,"Incremental")==0) return IncrementalEnum;
    13691369         else stage=12;
    13701370   }
    13711371   if(stage==12){
    1372               if (strcmp(name,"Indexed")==0) return IndexedEnum;
     1372              if (strcmp(name,"Incremental")==0) return IncrementalEnum;
     1373              else if (strcmp(name,"Indexed")==0) return IndexedEnum;
    13731374              else if (strcmp(name,"IntExternalResult")==0) return IntExternalResultEnum;
    13741375              else if (strcmp(name,"ElementInput")==0) return ElementInputEnum;
    … …  
    14891490              else if (strcmp(name,"Regionaloutput")==0) return RegionaloutputEnum;
    14901491              else if (strcmp(name,"Regular")==0) return RegularEnum;
    1491               else if (strcmp(name,"RecoveryAnalysis")==0) return RecoveryAnalysisEnum;
    14921492         else stage=13;
    14931493   }
    14941494   if(stage==13){
    1495               if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
     1495              if (strcmp(name,"RecoveryAnalysis")==0) return RecoveryAnalysisEnum;
     1496              else if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
    14961497              else if (strcmp(name,"SamplingAnalysis")==0) return SamplingAnalysisEnum;
    14971498              else if (strcmp(name,"SamplingSolution")==0) return SamplingSolutionEnum;

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

    r26676 r26750  
    4848                        for field=self.fields
    4949                                %Checking agreement of classes
    50                                 if(contains(field,'SMB'))
     50                                if(contains(field,'SMBautoregression'))
     51                                        mdname = structstoch.mdnames(find(strcmp(structstoch.fields,char(field))));
     52                                        if~(isequal(class(md.smb),char(mdname)))
     53                                                error('md.smb does not agree with stochasticforcing field %s', char(field));
     54                                        end
     55                                end
     56                                if(contains(field,'SMBforcing'))
    5157                                        mdname = structstoch.mdnames(find(strcmp(structstoch.fields,char(field))));
    5258                                        if~(isequal(class(md.smb),char(mdname)))
    … …  
    156162
    157163                                %Scaling covariance matrix (scale column-by-column and row-by-row)
    158                                 scaledfields = {'DefaultCalving','FloatingMeltRate','SMBautoregression'}; %list of fields that need scaling *1/yts
     164                                scaledfields = {'DefaultCalving','FloatingMeltRate','SMBautoregression','SMBforcing'}; %list of fields that need scaling *1/yts
    159165                                tempcovariance = self.covariance; %copy of covariance to avoid writing back in member variable
    160166                                for i=1:num_fields
    … …  
    199205                'FrictionWaterPressure',...
    200206                'FrontalForcingsRignotAutoregression',...
    201                 'SMBautoregression'
     207                'SMBautoregression',...
     208                'SMBforcing'
    202209                };
    203210        structure.mdnames = {...
    … …  
    206213                'friction',...
    207214                'frontalforcingsrignotautoregression',...
    208                 'SMBautoregression'
     215                'SMBautoregression',...
     216                'SMBforcing'
    209217        };
    210218end % }}}

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

    r26677 r26750  
    6767        for field in self.fields:
    6868            # Checking agreement of classes
    69             if 'SMB' in field:
     69            if 'SMBautoregression' in field:
     70                mdname = structstoch[field]
     71                if (type(md.smb).__name__ != mdname):
     72                    raise TypeError('md.smb does not agree with stochasticforcing field {}'.format(field))
     73            if 'SMBforcing' in field:
    7074                mdname = structstoch[field]
    7175                if (type(md.smb).__name__ != mdname):
    … …  
    144148
    145149            # Scaling covariance matrix (scale column-by-column and row-by-row)
    146             scaledfields = ['DefaultCalving','FloatingMeltRate','SMBautoregression'] # list of fields that need scaling * 1/yts
     150            scaledfields = ['DefaultCalving','FloatingMeltRate','SMBautoregression','SMBforcing'] # list of fields that need scaling * 1/yts
    147151            tempcovariance = np.copy(self.covariance)
    148152            for i in range(num_fields):
    … …  
    184188                     'FrictionWaterPressure': 'friction',
    185189                     'FrontalForcingsRignotAutoregression': 'frontalforcingsrignotautoregression',
    186                      'SMBautoregression': 'SMBautoregression'}
     190                     'SMBautoregression': 'SMBautoregression',
     191                     'SMBforcing': 'SMBforcing'}
    187192        return structure
    188193    # }}}
Note: See TracChangeset for help on using the changeset viewer.