Changeset 22450


Ignore:
Timestamp:
02/22/18 13:30:15 (7 years ago)
Author:
schlegel
Message:

CHG: add albedo and emissivity input to SMBgemb class

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

Legend:

Unmodified
Added
Removed

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

    r22448 r22450  
    6969                        iomodel->FetchDataToInput(elements,"md.smb.Tini",SmbTiniEnum);
    7070                        iomodel->FetchDataToInput(elements,"md.smb.Sizeini",SmbSizeiniEnum);
     71                        iomodel->FetchDataToInput(elements,"md.smb.aValue",SmbAValueEnum);
     72                        iomodel->FetchDataToInput(elements,"md.smb.teValue",SmbTeValueEnum);
    7173                        break;
    7274                case SMBpddEnum:

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

    r22448 r22450  
    25562556        IssmDouble eAir=0.0;
    25572557        IssmDouble pAir=0.0;
     2558        IssmDouble teValue=1.0;
     2559        IssmDouble aValue=0.0;
    25582560        int        aIdx=0;
    25592561        int        denIdx=0;
    … …  
    26592661        Input* eAir_input=this->GetInput(SmbEAirEnum); _assert_(eAir_input);
    26602662        Input* pAir_input=this->GetInput(SmbPAirEnum); _assert_(pAir_input);
     2663        Input* teValue_input=this->GetInput(SmbTeValueEnum); _assert_(teValue_input);
     2664        Input* aValue_input=this->GetInput(SmbAValueEnum); _assert_(aValue_input);
    26612665        Input* isinitialized_input=this->GetInput(SmbIsInitializedEnum); _assert_(isinitialized_input);
    26622666        /*Retrieve input values:*/
    26632667        Gauss* gauss=this->NewGauss(1); gauss->GaussPoint(0);
     2668
     2669        if (aIdx == 0) aValue_input->GetInputValue(&aValue,gauss);
    26642670
    26652671        zTop_input->GetInputValue(&zTop,gauss);
    … …  
    26732679        Tz_input->GetInputValue(&Tz,gauss);
    26742680        Vz_input->GetInputValue(&Vz,gauss);
     2681        teValue_input->GetInputValue(&teValue,gauss);
    26752682        isinitialized_input->GetInputValue(&isinitialized);
    26762683        /*}}}*/
    … …  
    27192726                        a = xNewZeroInit<IssmDouble>(m); for(int i=0;i<m;i++)a[i]=aini[0];         //set albedo equal to fresh snow [fraction]
    27202727                        T = xNewZeroInit<IssmDouble>(m); for(int i=0;i<m;i++)T[i]=Tmean;         //set initial grid cell temperature to the annual mean temperature [K]
    2721                         /*            /!\ Default value of T can not be retrived from SMBgemb.m (like other snow properties) because don't know Tmean yet when set default values.
    2722                                                           Default value of 0C given in SMBgemb.m is overwritten here with value of Tmean*/
     2728                        /*/!\ Default value of T can not be retrived from SMBgemb.m (like other snow properties)
     2729                         *    because don't know Tmean yet when set default values.
     2730                         *    Default value of 0C given in SMBgemb.m is overwritten here with value of Tmean*/
    27232731
    27242732                        //fixed lower temperature bounday condition - T is fixed
    … …  
    27972805                eAir_input->GetInputValue(&eAir,gauss,t);  //screen level vapor pressure [Pa]
    27982806                pAir_input->GetInputValue(&pAir,gauss,t);  // screen level air pressure [Pa]
     2807                teValue_input->GetInputValue(&teValue,gauss);  // screen level air pressure [Pa]
     2808                if(aIdx == 0) aValue_input->GetInputValue(&aValue,gauss);  // screen level air pressure [Pa]
    27992809                //_printf_("Time: " << t << " Ta: " << Ta << " V: " << V << " dlw: " << dlw << " dsw: " << dsw << " P: " << P << " eAir: " << eAir << " pAir: " << pAir << "\n");
    28002810                /*}}}*/
    … …  
    28042814
    28052815                /*Snow, firn and ice albedo:*/
    2806                 if(isalbedo)albedo(&a,aIdx,re,d,cldFrac,aIce, aSnow,T,W,P,EC,t0wet,t0dry,K,smb_dt,rho_ice,m,this->Sid());
     2816                if(isalbedo)albedo(&a,aIdx,re,d,cldFrac,aIce, aSnow,aValue,T,W,P,EC,t0wet,t0dry,K,smb_dt,rho_ice,m,this->Sid());
    28072817
    28082818
    … …  
    28142824
    28152825                /*Thermal profile computation:*/
    2816                 if(isthermal)thermo(&EC, &T, dz, d, swf, dlw, Ta, V, eAir, pAir, W[0], smb_dt, m, Vz, Tz,rho_ice,this->Sid());
     2826                if(isthermal)thermo(&EC, &T, dz, d, swf, dlw, Ta, V, eAir, pAir, teValue, W[0], smb_dt, m, Vz, Tz,rho_ice,this->Sid());
    28172827
    28182828                /*Change in thickness of top cell due to evaporation/condensation  assuming same density as top cell.
    … …  
    28212831
    28222832                /*Add snow/rain to top grid cell adjusting cell depth, temperature and density*/
    2823                 if(isaccumulation)accumulation(&T, &dz, &d, &W, &a, &re, &gdn, &gsp, &m, Ta, P, dzMin, aSnow,rho_ice,this->Sid());
     2833                if(isaccumulation)accumulation(&T, &dz, &d, &W, &a, &re, &gdn, &gsp, &m, aIdx, Ta, P, dzMin, aSnow,rho_ice,this->Sid());
    28242834
    28252835                /*Calculate water production, M [kg m-2] resulting from snow/ice temperature exceeding 273.15 deg K
    … …  
    28322842                /*Calculate upward longwave radiation flux [W m-2] not used in energy balance. Calculated for every
    28332843                 * sub-time step in thermo equations*/
    2834                 ulw = 5.67E-8 * pow(T[0],4.0);
     2844                ulw = 5.67E-8 * pow(T[0],4.0) * teValue;
    28352845
    28362846                /*Calculate net longwave [W m-2]*/
    … …  
    28402850                if(isturbulentflux)turbulentFlux(&shf, &lhf, &dayEC, Ta, T[0], V, eAir, pAir, d[0], W[0], Vz, Tz,rho_ice,this->Sid());
    28412851
    2842                 /*Verbose some resuls in debug mode: {{{*/
     2852                /*Verbose some results in debug mode: {{{*/
    28432853                if(VerboseSmb() && 0){
    28442854                        _printf_("smb log: count[" << count << "] m[" << m << "] "
    … …  
    28522862                                                << "gsp[" << cellsum(gsp,m)  << "] "
    28532863                                                << "swf[" << netSW << "] "
     2864                                                << "a[" << a << "] "
     2865                                                << "te[" << teValue << "] "
    28542866                                                << "\n");
    28552867                } /*}}}*/

  • issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/Gembx.cpp

    r22432 r22450  
    337337
    338338}  /*}}}*/
    339 void albedo(IssmDouble** pa, int aIdx, IssmDouble* re, IssmDouble* d, IssmDouble cldFrac, IssmDouble aIce, IssmDouble aSnow, IssmDouble* TK, IssmDouble* W, IssmDouble P, IssmDouble EC, IssmDouble t0wet, IssmDouble t0dry, IssmDouble K, IssmDouble dt, IssmDouble dIce, int m,int sid) { /*{{{*/
     339void albedo(IssmDouble** pa, int aIdx, IssmDouble* re, IssmDouble* d, IssmDouble cldFrac, IssmDouble aIce, IssmDouble aSnow, IssmDouble aValue, IssmDouble* TK, IssmDouble* W, IssmDouble P, IssmDouble EC, IssmDouble t0wet, IssmDouble t0dry, IssmDouble K, IssmDouble dt, IssmDouble dIce, int m,int sid) { /*{{{*/
    340340
    341341        //// Calculates Snow, firn and ice albedo as a function of:
     342        //   0 : direct input from aValue parameter
    342343        //   1 : effective grain radius (Gardner & Sharp, 2009)
    343344        //   2 : effective grain radius (Brun et al., 2009)
    … …  
    347348        //// Inputs
    348349        // aIdx      = albedo method to use
     350
     351        // Method 0
     352        //  aValue   = direct input value for albedo
    349353
    350354        // Methods 1 & 2
    … …  
    392396        const IssmDouble dSnow = 300;   // density of fresh snow [kg m-3]       
    393397
    394         if(aIdx==1){
     398        if (aIdx==0){
     399                  for(int i=0;i<m;i++)a[i] = aValue;
     400        }
     401        else if(aIdx==1){
    395402        //function of effective grain radius
    396403       
    … …  
    400407        //determine broadband albedo
    401408        a[0]= 1.48 - pow(S,-0.07);
     409                  for(int i=1;i<m;i++)a[i]=a[0];
    402410        }
    403411        else if(aIdx==2){
    … …  
    421429        // broadband surface albedo
    422430        a[0] = sF[0]*a0 + sF[1]*a1 + sF[2]*a2;
    423 
     431                  for(int i=1;i<m;i++)a[i]=a[0];
    424432        }
    425433        else if(aIdx==3){
    … …  
    429437        // calculate albedo
    430438        a[0] = aIce + (d[0] - dIce)*(aSnow - aIce) / (dSnow - dIce) + (0.05 * (cldFrac - 0.5));
     439                  for(int i=1;i<m;i++)a[i]=a[0];
    431440        }
    432441        else if(aIdx==4){
    … …  
    486495
    487496        }
    488         else _error_("albedo method switch should range from 1 to 4!");
     497        else _error_("albedo method switch should range from 0 to 4!");
    489498       
    490499        // Check for erroneous values
    … …  
    497506
    498507}  /*}}}*/
    499 void thermo(IssmDouble* pEC, IssmDouble** pT, IssmDouble* dz, IssmDouble* d, IssmDouble* swf, IssmDouble dlwrf, IssmDouble Ta, IssmDouble V, IssmDouble eAir, IssmDouble pAir, IssmDouble Ws, IssmDouble dt0, int m, IssmDouble Vz, IssmDouble Tz, IssmDouble dIce, int sid) { /*{{{*/
     508void thermo(IssmDouble* pEC, IssmDouble** pT, IssmDouble* dz, IssmDouble* d, IssmDouble* swf, IssmDouble dlwrf, IssmDouble Ta, IssmDouble V, IssmDouble eAir, IssmDouble pAir, IssmDouble teValue, IssmDouble Ws, IssmDouble dt0, int m, IssmDouble Vz, IssmDouble Tz, IssmDouble dIce, int sid) { /*{{{*/
    500509
    501510        /* ENGLACIAL THERMODYNAMICS*/
    … …  
    808817
    809818                // upward longwave contribution
    810                 ulw = - SB * pow(Ts,4.0) * dt;
     819                ulw = - SB * pow(Ts,4.0)* teValue * dt;
    811820                dT_ulw = ulw / TCs;
    812821               
    … …  
    10691078
    10701079} /*}}}*/
    1071 void accumulation(IssmDouble** pT, IssmDouble** pdz, IssmDouble** pd, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pm, IssmDouble T_air, IssmDouble P, IssmDouble dzMin, IssmDouble aSnow, IssmDouble dIce, int sid){ /*{{{*/
     1080void accumulation(IssmDouble** pT, IssmDouble** pdz, IssmDouble** pd, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pm, int aIdx, IssmDouble T_air, IssmDouble P, IssmDouble dzMin, IssmDouble aSnow, IssmDouble dIce, int sid){ /*{{{*/
    10721081
    10731082        // Adds precipitation and deposition to the model grid
    … …  
    11621171
    11631172                                // adjust a, re, gdn & gsp
    1164                                 a[0] = (aSnow * P + a[0] * mInit[0])/mass;
     1173                                if(aIdx>0)a[0] = (aSnow * P + a[0] * mInit[0])/mass;
    11651174                                re[0] = (reNew * P + re[0] * mInit[0])/mass;
    11661175                                gdn[0] = (gdnNew * P + gdn[0] * mInit[0])/mass;
    … …  
    17871796void densification(IssmDouble** pd,IssmDouble** pdz, IssmDouble* T, IssmDouble* re, int denIdx, IssmDouble C, IssmDouble dt, IssmDouble Tmean, IssmDouble dIce, int m, int sid){ /*{{{*/
    17881797
    1789         //// THIS NEEDS TO BE DOUBLE CHECKED AS THERE SEAMS TO BE LITTLE DENSIFICATION IN THE MODEL OUTOUT [MAYBE COMPATION IS COMPNSATED FOR BY TRACES OF SNOW???]
     1798        //// THIS NEEDS TO BE DOUBLE CHECKED AS THERE SEAMS TO BE LITTLE DENSIFICATION IN THE MODEL OUTOUT [MAYBE COMPACTION IS COMPENSATED FOR BY TRACES OF SNOW???]
    17901799
    17911800        //// FUNCTION INFO

  • issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.h

    r22429 r22450  
    2525IssmDouble Marbouty(IssmDouble T, IssmDouble d, IssmDouble dT);
    2626void grainGrowth(IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, IssmDouble* T,IssmDouble* dz,IssmDouble* d, IssmDouble* W,IssmDouble smb_dt,int m,int aIdx, int sid);
    27 void albedo(IssmDouble** a,int aIdx, IssmDouble* re, IssmDouble* d, IssmDouble cldFrac, IssmDouble aIce, IssmDouble aSnow, IssmDouble* T, IssmDouble* W, IssmDouble P, IssmDouble EC, IssmDouble t0wet, IssmDouble t0dry, IssmDouble K, IssmDouble dt, IssmDouble dIce, int m, int sid);
     27void albedo(IssmDouble** a,int aIdx, IssmDouble* re, IssmDouble* d, IssmDouble cldFrac, IssmDouble aIce, IssmDouble aSnow, IssmDouble aValue, IssmDouble* T, IssmDouble* W, IssmDouble P, IssmDouble EC, IssmDouble t0wet, IssmDouble t0dry, IssmDouble K, IssmDouble dt, IssmDouble dIce, int m, int sid);
    2828void shortwave(IssmDouble** pswf, int swIdx, int aIdx, IssmDouble dsw, IssmDouble as, IssmDouble* d, IssmDouble* dz, IssmDouble* re, IssmDouble dIce, int m, int sid);
    29 void thermo(IssmDouble* pEC, IssmDouble** T, IssmDouble* dz, IssmDouble* d, IssmDouble* swf, IssmDouble dlw, IssmDouble Ta, IssmDouble V, IssmDouble eAir, IssmDouble pAir, IssmDouble Ws, IssmDouble dt0, int m, IssmDouble Vz, IssmDouble Tz, IssmDouble dIce, int sid);
    30 void accumulation(IssmDouble** pT, IssmDouble** pdz, IssmDouble** pd, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pm, IssmDouble Ta, IssmDouble P, IssmDouble dzMin, IssmDouble aSnow, IssmDouble dIce, int sid);
     29void thermo(IssmDouble* pEC, IssmDouble** T, IssmDouble* dz, IssmDouble* d, IssmDouble* swf, IssmDouble dlw, IssmDouble Ta, IssmDouble V, IssmDouble eAir, IssmDouble pAir, IssmDouble teValue, IssmDouble Ws, IssmDouble dt0, int m, IssmDouble Vz, IssmDouble Tz, IssmDouble dIce, int sid);
     30void accumulation(IssmDouble** pT, IssmDouble** pdz, IssmDouble** pd, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pm, int aIdx,IssmDouble Ta, IssmDouble P, IssmDouble dzMin, IssmDouble aSnow, IssmDouble dIce, int sid);
    3131void melt(IssmDouble* pM, IssmDouble* pR, IssmDouble* pmAdd, IssmDouble* pdz_add, IssmDouble** pT, IssmDouble** pd, IssmDouble** pdz, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pn, IssmDouble dzMin, IssmDouble zMax, IssmDouble zMin, IssmDouble zTop, IssmDouble dIce, int sid);
    3232void densification(IssmDouble** pd,IssmDouble** pdz, IssmDouble* T, IssmDouble* re, int denIdx, IssmDouble C, IssmDouble dt, IssmDouble Tmean, IssmDouble dIce, int m, int sid);

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

    r22448 r22450  
    459459        SmbAEnum,
    460460        SmbTEnum,
     461        SmbAValueEnum,
     462        SmbTeValueEnum,
    461463        SmbIsgraingrowthEnum,
    462464        SmbIsalbedoEnum,

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

    r22448 r22450  
    461461                case SmbAEnum : return "SmbA";
    462462                case SmbTEnum : return "SmbT";
     463                case SmbAValueEnum : return "SmbAValue";
     464                case SmbTeValueEnum : return "SmbTeValue";
    463465                case SmbIsgraingrowthEnum : return "SmbIsgraingrowth";
    464466                case SmbIsalbedoEnum : return "SmbIsalbedo";

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

    r22448 r22450  
    470470              else if (strcmp(name,"SmbA")==0) return SmbAEnum;
    471471              else if (strcmp(name,"SmbT")==0) return SmbTEnum;
     472              else if (strcmp(name,"SmbAValue")==0) return SmbAValueEnum;
     473              else if (strcmp(name,"SmbTeValue")==0) return SmbTeValueEnum;
    472474              else if (strcmp(name,"SmbIsgraingrowth")==0) return SmbIsgraingrowthEnum;
    473475              else if (strcmp(name,"SmbIsalbedo")==0) return SmbIsalbedoEnum;
    … …  
    504506              else if (strcmp(name,"SMBd18opdd")==0) return SMBd18opddEnum;
    505507              else if (strcmp(name,"SmbDpermil")==0) return SmbDpermilEnum;
    506               else if (strcmp(name,"SmbF")==0) return SmbFEnum;
    507               else if (strcmp(name,"SMBgradients")==0) return SMBgradientsEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"SmbMonthlytemperatures")==0) return SmbMonthlytemperaturesEnum;
     511              if (strcmp(name,"SmbF")==0) return SmbFEnum;
     512              else if (strcmp(name,"SMBgradients")==0) return SMBgradientsEnum;
     513              else if (strcmp(name,"SmbMonthlytemperatures")==0) return SmbMonthlytemperaturesEnum;
    512514              else if (strcmp(name,"SmbHref")==0) return SmbHrefEnum;
    513515              else if (strcmp(name,"SmbSmbref")==0) return SmbSmbrefEnum;
    … …  
    627629              else if (strcmp(name,"GiaW")==0) return GiaWEnum;
    628630              else if (strcmp(name,"SaveResults")==0) return SaveResultsEnum;
    629               else if (strcmp(name,"BoolExternalResult")==0) return BoolExternalResultEnum;
    630               else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
     634              if (strcmp(name,"BoolExternalResult")==0) return BoolExternalResultEnum;
     635              else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
     636              else if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
    635637              else if (strcmp(name,"IntExternalResult")==0) return IntExternalResultEnum;
    636638              else if (strcmp(name,"IntMatExternalResult")==0) return IntMatExternalResultEnum;
    … …  
    750752              else if (strcmp(name,"VyObs")==0) return VyObsEnum;
    751753              else if (strcmp(name,"Numberedcostfunction")==0) return NumberedcostfunctionEnum;
    752               else if (strcmp(name,"Absolute")==0) return AbsoluteEnum;
    753               else if (strcmp(name,"Incremental")==0) return IncrementalEnum;
    754754         else stage=7;
    755755   }
    756756   if(stage==7){
    757               if (strcmp(name,"AugmentedLagrangianR")==0) return AugmentedLagrangianREnum;
     757              if (strcmp(name,"Absolute")==0) return AbsoluteEnum;
     758              else if (strcmp(name,"Incremental")==0) return IncrementalEnum;
     759              else if (strcmp(name,"AugmentedLagrangianR")==0) return AugmentedLagrangianREnum;
    758760              else if (strcmp(name,"AugmentedLagrangianRhop")==0) return AugmentedLagrangianRhopEnum;
    759761              else if (strcmp(name,"AugmentedLagrangianRlambda")==0) return AugmentedLagrangianRlambdaEnum;
    … …  
    873875              else if (strcmp(name,"LoveKernelsImag")==0) return LoveKernelsImagEnum;
    874876              else if (strcmp(name,"EsaUmotion")==0) return EsaUmotionEnum;
    875               else if (strcmp(name,"EsaNmotion")==0) return EsaNmotionEnum;
    876               else if (strcmp(name,"EsaEmotion")==0) return EsaEmotionEnum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"EsaXmotion")==0) return EsaXmotionEnum;
     880              if (strcmp(name,"EsaNmotion")==0) return EsaNmotionEnum;
     881              else if (strcmp(name,"EsaEmotion")==0) return EsaEmotionEnum;
     882              else if (strcmp(name,"EsaXmotion")==0) return EsaXmotionEnum;
    881883              else if (strcmp(name,"EsaYmotion")==0) return EsaYmotionEnum;
    882884              else if (strcmp(name,"EsaHemisphere")==0) return EsaHemisphereEnum;
    … …  
    996998              else if (strcmp(name,"BalancethicknessSoftSolution")==0) return BalancethicknessSoftSolutionEnum;
    997999              else if (strcmp(name,"BalancevelocityAnalysis")==0) return BalancevelocityAnalysisEnum;
    998               else if (strcmp(name,"BalancevelocitySolution")==0) return BalancevelocitySolutionEnum;
    999               else if (strcmp(name,"L2ProjectionEPLAnalysis")==0) return L2ProjectionEPLAnalysisEnum;
    10001000         else stage=9;
    10011001   }
    10021002   if(stage==9){
    1003               if (strcmp(name,"L2ProjectionBaseAnalysis")==0) return L2ProjectionBaseAnalysisEnum;
     1003              if (strcmp(name,"BalancevelocitySolution")==0) return BalancevelocitySolutionEnum;
     1004              else if (strcmp(name,"L2ProjectionEPLAnalysis")==0) return L2ProjectionEPLAnalysisEnum;
     1005              else if (strcmp(name,"L2ProjectionBaseAnalysis")==0) return L2ProjectionBaseAnalysisEnum;
    10041006              else if (strcmp(name,"BedSlopeSolution")==0) return BedSlopeSolutionEnum;
    10051007              else if (strcmp(name,"DamageEvolutionSolution")==0) return DamageEvolutionSolutionEnum;
    … …  
    11191121              else if (strcmp(name,"DataSet")==0) return DataSetEnum;
    11201122              else if (strcmp(name,"Constraints")==0) return ConstraintsEnum;
    1121               else if (strcmp(name,"Loads")==0) return LoadsEnum;
    1122               else if (strcmp(name,"Materials")==0) return MaterialsEnum;
    11231123         else stage=10;
    11241124   }
    11251125   if(stage==10){
    1126               if (strcmp(name,"Nodes")==0) return NodesEnum;
     1126              if (strcmp(name,"Loads")==0) return LoadsEnum;
     1127              else if (strcmp(name,"Materials")==0) return MaterialsEnum;
     1128              else if (strcmp(name,"Nodes")==0) return NodesEnum;
    11271129              else if (strcmp(name,"Contours")==0) return ContoursEnum;
    11281130              else if (strcmp(name,"Parameters")==0) return ParametersEnum;

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

    r22267 r22450  
    3636                Tz    = NaN; %height above ground at which temperature (T) was sampled [m]
    3737                Vz    = NaN; %height above ground at which wind (V) eas sampled [m]
     38
     39                %optional inputs:
     40                aValue = NaN; %Albedo forcing at every element.  Used only if aIdx == 0.
     41                teValue = NaN; %Outward longwave radiation thermal emissivity forcing at every element (default in code is 1)
    3842       
    3943                % Initialization of snow properties
    … …  
    5155                %settings:
    5256                aIdx   = NaN; %method for calculating albedo and subsurface absorption (default is 1)
    53                               % 1: effective grain radius [Gardner & Sharp, 2009]
     57                           % 0: direct input from aValue parameter
     58                           % 1: effective grain radius [Gardner & Sharp, 2009]
    5459                                          % 2: effective grain radius [Brun et al., 2009]
    5560                                          % 3: density and cloud amount [Greuell & Konzelmann, 1994]
    … …  
    114119                        self.eAir=project3d(md,'vector',self.eAir,'type','node');
    115120                        self.pAir=project3d(md,'vector',self.pAir,'type','node');
     121
     122                        if aIdx == 0 & ~isnan(self.aValue)
     123                                self.aValue=project3d(md,'vector',self.aValue,'type','node');
     124                        end
     125                        if ~isnan(self.teValue)
     126                                self.teValue=project3d(md,'vector',self.teValue,'type','node');
     127                        end
     128
    116129
    117130                end % }}}
    … …  
    150163                self.t0dry = 30;
    151164                self.K = 7;
     165
     166                self.teValue = ones(mesh.numberofelements,1);
     167                self.aValue = self.aSnow*ones(mesh.numberofelements,1);
    152168       
    153169                self.Dzini=0.05*ones(mesh.numberofelements,2);
    … …  
    167183                function md = checkconsistency(self,md,solution,analyses) % {{{
    168184
    169 
    170185                        md = checkfield(md,'fieldname','smb.isgraingrowth','values',[0 1]);
    171186                        md = checkfield(md,'fieldname','smb.isalbedo','values',[0 1]);
    … …  
    189204                        md = checkfield(md,'fieldname','smb.Vz','size',[md.mesh.numberofelements 1],'NaN',1,'Inf',1,'>=',0,'<=',5000);
    190205
    191                         md = checkfield(md,'fieldname','smb.aIdx','NaN',1,'Inf',1,'values',[1,2,3,4]);
     206                        md = checkfield(md,'fieldname','smb.teValue','timeseries',1,'NaN',1,'Inf',1,'>=',0,'<=',1);
     207
     208                        md = checkfield(md,'fieldname','smb.aIdx','NaN',1,'Inf',1,'values',[0,1,2,3,4]);
    192209                        md = checkfield(md,'fieldname','smb.swIdx','NaN',1,'Inf',1,'values',[0,1]);
    193210                        md = checkfield(md,'fieldname','smb.denIdx','NaN',1,'Inf',1,'values',[1,2,3,4,5]);
    … …  
    201218
    202219                        switch self.aIdx,
     220                                case 0
     221                                        md = checkfield(md,'fieldname','smb.aValue','timeseries',1,'NaN',1,'Inf',1,'>=',0,'<=',1);
    203222                                case {1 2}
    204223                                        md = checkfield(md,'fieldname','smb.aSnow','NaN',1,'Inf',1,'>=',.64,'<=',.89);
    … …  
    252271                        fielddisplay(self,'outputFreq','output frequency in days (default is monthly, 30)');
    253272                        fielddisplay(self,'aIdx',{'method for calculating albedo and subsurface absorption (default is 1)',...
     273                                               '0: direct input from aValue parameter',...
    254274                                                                        '1: effective grain radius [Gardner & Sharp, 2009]',...
    255275                                                                        '2: effective grain radius [Brun et al., 2009]',...
    256276                                                                        '3: density and cloud amount [Greuell & Konzelmann, 1994]',...
    257277                                                                        '4: exponential time decay & wetness [Bougamont & Bamber, 2005]'});
     278
     279                        fielddisplay(self,'teValue','Outward longwave radiation thermal emissivity forcing at every element (default in code is 1)');
    258280                               
    259281                        %snow properties init
    260                         fielddisplay(self,'Dzini','Initial cel depth when restart [m]');
     282                        fielddisplay(self,'Dzini','Initial cell depth when restart [m]');
    261283                        fielddisplay(self,'Dini','Initial snow density when restart [kg m-3]');
    262284                        fielddisplay(self,'Reini','Initial grain size when restart [mm]');
    … …  
    272294                        %additional albedo parameters:
    273295                        switch self.aIdx
     296                        case 0
     297                                fielddisplay(self,'aValue','Albedo forcing at every element.  Used only if aIdx == 0.');
    274298                        case {1 2}
    275299                                fielddisplay(self,'aSnow','new snow albedo (0.64 - 0.89)');
    … …  
    339363                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','t0dry','format','Double');
    340364                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','K','format','Double');
     365
     366                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','aValue','format','DoubleMat','mattype',2,'timeserieslength',md.mesh.numberofelements+1,'yts',md.constants.yts);
     367                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','teValue','format','DoubleMat','mattype',2,'timeserieslength',md.mesh.numberofelements+1,'yts',md.constants.yts);
    341368           
    342369                        %snow properties init

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

    r22267 r22450  
    4242                Tz    = float('NaN')    #height above ground at which temperature (T) was sampled [m]
    4343                Vz    = float('NaN')    #height above ground at which wind (V) eas sampled [m]
     44
     45                #optional inputs:
     46                aValue  = float('NaN') #Albedo forcing at every element.  Used only if aIdx == 0.
     47                teValue = float('NaN') #Outward longwave radiation thermal emissivity forcing at every element (default in code is 1)
    4448       
    4549                # Initialization of snow properties
    … …  
    5761                #settings:
    5862                aIdx   = float('NaN')   #method for calculating albedo and subsurface absorption (default is 1)
    59                                         # 1: effective grain radius [Gardner & Sharp, 2009]
     63                           # 0: direct input from aValue parameter
     64                                          # 1: effective grain radius [Gardner & Sharp, 2009]
    6065                                          # 2: effective grain radius [Brun et al., 2009]
    6166                                          # 3: density and cloud amount [Greuell & Konzelmann, 1994]
    … …  
    135140                string = "%s\n%s"%(string,fielddisplay(self,'outputFreq','output frequency in days (default is monthly, 30)'))
    136141                string = "%s\n%s"%(string,fielddisplay(self,'aIdx',['method for calculating albedo and subsurface absorption (default is 1)',
     142                                 '0: direct input from aValue parameter',
    137143                                                '1: effective grain radius [Gardner & Sharp, 2009]',
    138144                                                '2: effective grain radius [Brun et al., 2009]',
    139145                                                '3: density and cloud amount [Greuell & Konzelmann, 1994]',
    140146                                                '4: exponential time decay & wetness [Bougamont & Bamber, 2005]']))
     147
     148                string = "%s\n%s"%(string,fielddisplay(self,'teValue','Outward longwave radiation thermal emissivity forcing at every element (default in code is 1)'))
    141149                               
    142150                #snow properties init
    143                 string = "%s\n%s"%(string,fielddisplay(self,'Dzini','Initial cel depth when restart [m]'))
     151                string = "%s\n%s"%(string,fielddisplay(self,'Dzini','Initial cell depth when restart [m]'))
    144152                string = "%s\n%s"%(string,fielddisplay(self,'Dini','Initial snow density when restart [kg m-3]'))
    145153                string = "%s\n%s"%(string,fielddisplay(self,'Reini','Initial grain size when restart [mm]'))
    … …  
    156164                        string = "%s\n%s"%(string,fielddisplay(self,'aSnow','new snow albedo (0.64 - 0.89)'))
    157165                        string = "%s\n%s"%(string,fielddisplay(self,'aIce','albedo of ice (0.27-0.58)'))
     166                elif self.aIdx == 0:
     167                        string = "%s\n%s"%(string,fielddisplay(self,'aValue','Albedo forcing at every element.  Used only if aIdx == 0.'))
    158168                elif self.aIdx == 3:
    159169                        string = "%s\n%s"%(string,fielddisplay(self,'cldFrac','average cloud amount'))
    … …  
    183193                self.eAir = project3d(md,'vector',self.eAir,'type','node')
    184194                self.pAir = project3d(md,'vector',self.pAir,'type','node')
     195
     196                if aIdx == 0 and np.isnan(self.aValue):
     197                        self.aValue=project3d(md,'vector',self.aValue,'type','node');
     198                if np.isnan(self.teValue):
     199                        self.teValue=project3d(md,'vector',self.teValue,'type','node');
     200
    185201                return self
    186202        #}}}
    … …  
    211227                self.zY = 1.10*np.ones((mesh.numberofelements,))
    212228                self.outputFreq = 30
    213                
     229
    214230                #additional albedo parameters
    215231                self.aSnow = 0.85
    … …  
    219235                self.t0dry = 30
    220236                self.K = 7
     237
     238                self.teValue = np.ones((mesh.numberofelements,));
     239                self.aValue = self.aSnow*np.ones(mesh.numberofelements,);
    221240       
    222241                self.Dzini = 0.05*np.ones((mesh.numberofelements,2))
    … …  
    257276                md = checkfield(md,'fieldname','smb.Vz','size',[md.mesh.numberofelements],'NaN',1,'Inf',1,'> = ',0,'< = ',5000)
    258277
    259                 md = checkfield(md,'fieldname','smb.aIdx','NaN',1,'Inf',1,'values',[1,2,3,4])
     278                md = checkfield(md,'fieldname','smb.teValue','timeseries',1,'NaN',1,'Inf',1,'>=',0,'<=',1);
     279
     280                md = checkfield(md,'fieldname','smb.aIdx','NaN',1,'Inf',1,'values',[0,1,2,3,4])
    260281                md = checkfield(md,'fieldname','smb.swIdx','NaN',1,'Inf',1,'values',[0,1])
    261282                md = checkfield(md,'fieldname','smb.denIdx','NaN',1,'Inf',1,'values',[1,2,3,4,5])
    … …  
    271292                        md = checkfield(md,'fieldname','smb.aSnow','NaN',1,'Inf',1,'> = ',.64,'< = ',.89)
    272293                        md = checkfield(md,'fieldname','smb.aIce','NaN',1,'Inf',1,'> = ',.27,'< = ',.58)
     294                elif self.aIdx == 0:
     295                        md = checkfield(md,'fieldname','smb.aValue','timeseries',1,'NaN',1,'Inf',1,'>=',0,'<=',1);
    273296                elif self.aIdx == 3:
    274297                        md = checkfield(md,'fieldname','smb.cldFrac','NaN',1,'Inf',1,'> = ',0,'< = ',1)
    … …  
    333356                WriteData(fid,prefix,'object',self,'class','smb','fieldname','t0dry','format','Double')
    334357                WriteData(fid,prefix,'object',self,'class','smb','fieldname','K','format','Double')
     358
     359                WriteData(fid,prefix,'object',self,'class','smb','fieldname','aValue','format','DoubleMat','mattype',2,'timeserieslength',md.mesh.numberofelements+1,'yts',md.constants.yts);
     360                WriteData(fid,prefix,'object',self,'class','smb','fieldname','teValue','format','DoubleMat','mattype',2,'timeserieslength',md.mesh.numberofelements+1,'yts',md.constants.yts);
    335361           
    336362                #snow properties init
Note: See TracChangeset for help on using the changeset viewer.