Changeset 23468


Ignore:
Timestamp:
11/16/18 18:29:43 (6 years ago)
Author:
schlegel
Message:

CHG: add snow density options to GEMB

Location:
issm/trunk-jpl
Files:
13 edited

Legend:

Unmodified
Added
Removed

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

    r23366 r23468  
    5757                        iomodel->FetchDataToInput(elements,"md.smb.zMin",SmbZMinEnum);
    5858                        iomodel->FetchDataToInput(elements,"md.smb.Tmean",SmbTmeanEnum);
     59                        iomodel->FetchDataToInput(elements,"md.smb.Vmean",SmbVmeanEnum);
    5960                        iomodel->FetchDataToInput(elements,"md.smb.C",SmbCEnum);
    6061                        iomodel->FetchDataToInput(elements,"md.smb.Tz",SmbTzEnum);
    … …  
    184185                        parameters->AddObject(iomodel->CopyConstantObject("md.smb.swIdx",SmbSwIdxEnum));
    185186                        parameters->AddObject(iomodel->CopyConstantObject("md.smb.denIdx",SmbDenIdxEnum));
     187                        parameters->AddObject(iomodel->CopyConstantObject("md.smb.dsnowIdx",SmbDsnowIdxEnum));
    186188                        parameters->AddObject(iomodel->CopyConstantObject("md.smb.cldFrac",SmbCldFracEnum));
    187189                        parameters->AddObject(iomodel->CopyConstantObject("md.smb.t0wet",SmbT0wetEnum));

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

    r23366 r23468  
    30263026        IssmDouble dzMin=0.0;
    30273027        IssmDouble Tmean=0.0;
     3028        IssmDouble Vmean=0.0;
    30283029        IssmDouble C=0.0;
    30293030        IssmDouble Tz,Vz=0.0;
    … …  
    30433044        int        aIdx=0;
    30443045        int        denIdx=0;
     3046        int        dsnowIdx=0;
    30453047        int        swIdx=0;
    30463048        IssmDouble cldFrac,t0wet, t0dry, K;
    … …  
    31133115        parameters->FindParam(&denIdx,SmbDenIdxEnum);
    31143116        parameters->FindParam(&swIdx,SmbSwIdxEnum);
     3117        parameters->FindParam(&dsnowIdx,SmbDsnowIdxEnum);
    31153118        parameters->FindParam(&cldFrac,SmbCldFracEnum);
    31163119        parameters->FindParam(&t0wet,SmbT0wetEnum);
    … …  
    31383141        Input* zY_input=this->GetInput(SmbZYEnum); _assert_(zY_input);
    31393142        Input* Tmean_input=this->GetInput(SmbTmeanEnum); _assert_(Tmean_input);
     3143        Input* Vmean_input=this->GetInput(SmbVmeanEnum); _assert_(Vmean_input);
    31403144        Input* C_input=this->GetInput(SmbCEnum); _assert_(C_input);
    31413145        Input* Tz_input=this->GetInput(SmbTzEnum); _assert_(Tz_input);
    … …  
    31613165        zY_input->GetInputValue(&zY,gauss);
    31623166        Tmean_input->GetInputValue(&Tmean,gauss);
     3167        Vmean_input->GetInputValue(&Vmean,gauss);
    31633168        C_input->GetInputValue(&C,gauss);
    31643169        Tz_input->GetInputValue(&Tz,gauss);
    … …  
    33163321
    33173322                /*Add snow/rain to top grid cell adjusting cell depth, temperature and density*/
    3318                 if(isaccumulation)accumulation(&T, &dz, &d, &W, &a, &re, &gdn, &gsp, &m, aIdx, Ta, P, dzMin, aSnow,rho_ice,this->Sid());
     3323                if(isaccumulation)accumulation(&T, &dz, &d, &W, &a, &re, &gdn, &gsp, &m, aIdx, dsnowIdx, Tmean, Ta, P, dzMin, aSnow, C, V, Vmean, rho_ice,this->Sid());
    33193324
    33203325                /*Calculate water production, M [kg m-2] resulting from snow/ice temperature exceeding 273.15 deg K
    … …  
    33273332                /*Calculate upward longwave radiation flux [W m-2] not used in energy balance. Calculated for every
    33283333                 * sub-time step in thermo equations*/
    3329                 ulw = 5.67E-8 * pow(T[0],4.0) * teValue;
     3334                ulw = 5.67E-8 * pow(T[0],4.0) * teValue; // + deltatest here
    33303335
    33313336                /*Calculate net longwave [W m-2]*/

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

    r23243 r23468  
    820820
    821821                // upward longwave contribution
    822                 ulw = - (SB * pow(Ts,4.0)* teValue) * dt ;
     822                ulw = - (SB * pow(Ts,4.0)* teValue) * dt ; //+20
    823823                dT_ulw = ulw / TCs;
    824824
    … …  
    10761076
    10771077} /*}}}*/
    1078 void 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){ /*{{{*/
     1078void accumulation(IssmDouble** pT, IssmDouble** pdz, IssmDouble** pd, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pm, int aIdx, int dsnowIdx, IssmDouble Tmean, IssmDouble T_air, IssmDouble P, IssmDouble dzMin, IssmDouble aSnow, IssmDouble C, IssmDouble V, IssmDouble Vmean, IssmDouble dIce, int sid){ /*{{{*/
    10791079
    10801080        // Adds precipitation and deposition to the model grid
    … …  
    10861086           adjusts the properties of the top grid cell to account for accumulation
    10871087           T_air & T = Air and top grid cell temperatures [K]
     1088                Tmean =  average surface temperature [K]
     1089                V =  wind velocity [m s-1]
     1090                C =  average accumulation rate [kg m-2 yr-1]
    10881091           dz = topgrid cell length [m]
    10891092           d = density of top grid gell [kg m-3]
    … …  
    10951098        // MAIN FUNCTION
    10961099        // specify constants
    1097         const IssmDouble dSnow = 150;    // density of snow [kg m-3]
     1100        IssmDouble dSnow = 150;    // density of snow [kg m-3]
    10981101        const IssmDouble reNew = 0.1;    // new snow grain size [mm]
    10991102        const IssmDouble gdnNew = 1.0;     // new snow dendricity
    … …  
    11301133        gsp=*pgsp;
    11311134        m=*pm;
     1135
     1136        //Density of fresh snow [kg m-3]
     1137        switch (dsnowIdx){
     1138                case 0: // Default value defined above
     1139                        break;
     1140
     1141                case 1: // Density of Antarctica snow
     1142                        dSnow = 350;
     1143                        break;
     1144
     1145                case 2: // Density of Greenland snow, Fausto et al., 2008
     1146                        dSnow = 315;
     1147                        break;
     1148
     1149                case 3: //Surface snow accumulation density from Kaspers et al., 2004, Antarctica
     1150                        //dSnow = alpha1 + beta1*T + delta1*C + epsilon1*W
     1151                        //     7.36x10-2  1.06x10-3  6.69x10-2  4.77x10-3
     1152                        dSnow=(7.362e-2 + 1.06e-3*Tmean + 6.69e-2*C/1000 + 4.77e-3*Vmean)*1000;
     1153                        break;
     1154
     1155                case 4: // Kuipers Munneke and others (2015), Greenland
     1156                        dSnow = 481.0 + 4.834*(Tmean-CtoK);
     1157                        break;
     1158        }
    11321159
    11331160        // determine initial mass
    … …  
    19151942                                c0arth = 0.07 * H;
    19161943                                c1arth = 0.03 * H;
    1917                                 M0 = max(1.435 - (0.151 * log(C)),0.25);
    1918                                 M1 = max(2.366 - (0.293 * log(C)),0.25);
     1944                                M0 = max(1.6599 - (0.1724 * log(C)),0.25);
     1945                                M1 = max(2.0102 - (0.2458 * log(C)),0.25);
     1946                                //M0 = max(1.435 - (0.151 * log(C)),0.25);
     1947                                //M1 = max(2.366 - (0.293 * log(C)),0.25);
    19191948                                c0 = M0*c0arth;
    19201949                                c1 = M1*c1arth;
    … …  
    19261955                                c0arth = 0.07 * H;
    19271956                                c1arth = 0.03 * H;
    1928                                 M0 = max(1.042 - (0.0916 * log(C)),0.25);
    1929                                 M1 = max(1.734 - (0.2039 * log(C)),0.25);
     1957                                //M0 = max(1.042 - (0.0916 * log(C)),0.25);
     1958                                //M1 = max(1.734 - (0.2039 * log(C)),0.25);
     1959                                M0 = max(1.6201 - (0.1450 * log(C)),0.25);
     1960                                M1 = max(2.5577 - (0.2899 * log(C)),0.25);
    19301961                                c0 = M0*c0arth;
    19311962                                c1 = M1*c1arth;

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

    r23366 r23468  
    2929void shortwave(IssmDouble** pswf, int swIdx, int aIdx, IssmDouble dsw, IssmDouble as, IssmDouble* d, IssmDouble* dz, IssmDouble* re, IssmDouble dIce, int m, int sid);
    3030void 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 thermo_scaling, IssmDouble dIce, int sid);
    31 void 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);
     31void accumulation(IssmDouble** pT, IssmDouble** pdz, IssmDouble** pd, IssmDouble** pW, IssmDouble** pa, IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, int* pm, int aIdx, int dsnowIdx, IssmDouble Tmean, IssmDouble Ta, IssmDouble P, IssmDouble dzMin, IssmDouble aSnow, IssmDouble C, IssmDouble V, IssmDouble Vmean, IssmDouble dIce, int sid);
    3232void 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);
    3333void 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

    r23366 r23468  
    279279        SmbAIceEnum,
    280280        SmbAIdxEnum,
     281        SmbDsnowIdxEnum,
    281282        SmbASnowEnum,
    282283        SmbCldFracEnum,
    … …  
    607608        SmbTiniEnum,
    608609        SmbTmeanEnum,
     610        SmbVmeanEnum,
    609611        SmbTzEnum,
    610612        SmbVEnum,

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

    r23366 r23468  
    287287                case SmbAIceEnum : return "SmbAIce";
    288288                case SmbAIdxEnum : return "SmbAIdx";
     289                case SmbDsnowIdxEnum : return "SmbDsnowIdx";
    289290                case SmbASnowEnum : return "SmbASnow";
    290291                case SmbCldFracEnum : return "SmbCldFrac";
    … …  
    613614                case SmbTiniEnum : return "SmbTini";
    614615                case SmbTmeanEnum : return "SmbTmean";
     616                case SmbVmeanEnum : return "SmbVmean";
    615617                case SmbTzEnum : return "SmbTz";
    616618                case SmbVEnum : return "SmbV";

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

    r23366 r23468  
    293293              else if (strcmp(name,"SmbAIce")==0) return SmbAIceEnum;
    294294              else if (strcmp(name,"SmbAIdx")==0) return SmbAIdxEnum;
     295              else if (strcmp(name,"SmbDsnowIdx")==0) return SmbDsnowIdxEnum;
    295296              else if (strcmp(name,"SmbASnow")==0) return SmbASnowEnum;
    296297              else if (strcmp(name,"SmbCldFrac")==0) return SmbCldFracEnum;
    … …  
    382383              else if (strcmp(name,"TransientIsmasstransport")==0) return TransientIsmasstransportEnum;
    383384              else if (strcmp(name,"TransientIsmovingfront")==0) return TransientIsmovingfrontEnum;
    384               else if (strcmp(name,"TransientIsoceancoupling")==0) return TransientIsoceancouplingEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"TransientIsslr")==0) return TransientIsslrEnum;
     388              if (strcmp(name,"TransientIsoceancoupling")==0) return TransientIsoceancouplingEnum;
     389              else if (strcmp(name,"TransientIsslr")==0) return TransientIsslrEnum;
    389390              else if (strcmp(name,"TransientIssmb")==0) return TransientIssmbEnum;
    390391              else if (strcmp(name,"TransientIsstressbalance")==0) return TransientIsstressbalanceEnum;
    … …  
    505506              else if (strcmp(name,"HydrologyDrainageRate")==0) return HydrologyDrainageRateEnum;
    506507              else if (strcmp(name,"Ice")==0) return IceEnum;
    507               else if (strcmp(name,"IceMaskNodeActivation")==0) return IceMaskNodeActivationEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"Input")==0) return InputEnum;
     511              if (strcmp(name,"IceMaskNodeActivation")==0) return IceMaskNodeActivationEnum;
     512              else if (strcmp(name,"Input")==0) return InputEnum;
    512513              else if (strcmp(name,"InversionCostFunctionsCoefficients")==0) return InversionCostFunctionsCoefficientsEnum;
    513514              else if (strcmp(name,"InversionSurfaceObs")==0) return InversionSurfaceObsEnum;
    … …  
    625626              else if (strcmp(name,"SmbTini")==0) return SmbTiniEnum;
    626627              else if (strcmp(name,"SmbTmean")==0) return SmbTmeanEnum;
     628              else if (strcmp(name,"SmbVmean")==0) return SmbVmeanEnum;
    627629              else if (strcmp(name,"SmbTz")==0) return SmbTzEnum;
    628630              else if (strcmp(name,"SmbV")==0) return SmbVEnum;
    629               else if (strcmp(name,"SmbVz")==0) return SmbVzEnum;
    630               else if (strcmp(name,"SmbW")==0) return SmbWEnum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"SmbWini")==0) return SmbWiniEnum;
     634              if (strcmp(name,"SmbVz")==0) return SmbVzEnum;
     635              else if (strcmp(name,"SmbW")==0) return SmbWEnum;
     636              else if (strcmp(name,"SmbWini")==0) return SmbWiniEnum;
    635637              else if (strcmp(name,"SmbZMax")==0) return SmbZMaxEnum;
    636638              else if (strcmp(name,"SmbZMin")==0) return SmbZMinEnum;
    … …  
    750752              else if (strcmp(name,"DamageEvolutionSolution")==0) return DamageEvolutionSolutionEnum;
    751753              else if (strcmp(name,"DataSet")==0) return DataSetEnum;
    752               else if (strcmp(name,"DatasetInput")==0) return DatasetInputEnum;
    753               else if (strcmp(name,"DataSetParam")==0) return DataSetParamEnum;
    754754         else stage=7;
    755755   }
    756756   if(stage==7){
    757               if (strcmp(name,"DefaultAnalysis")==0) return DefaultAnalysisEnum;
     757              if (strcmp(name,"DatasetInput")==0) return DatasetInputEnum;
     758              else if (strcmp(name,"DataSetParam")==0) return DataSetParamEnum;
     759              else if (strcmp(name,"DefaultAnalysis")==0) return DefaultAnalysisEnum;
    758760              else if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum;
    759761              else if (strcmp(name,"DegreeOfChannelization")==0) return DegreeOfChannelizationEnum;
    … …  
    873875              else if (strcmp(name,"L2ProjectionEPLAnalysis")==0) return L2ProjectionEPLAnalysisEnum;
    874876              else if (strcmp(name,"LACrouzeixRaviart")==0) return LACrouzeixRaviartEnum;
    875               else if (strcmp(name,"LambdaS")==0) return LambdaSEnum;
    876               else if (strcmp(name,"LATaylorHood")==0) return LATaylorHoodEnum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"LevelsetAnalysis")==0) return LevelsetAnalysisEnum;
     880              if (strcmp(name,"LambdaS")==0) return LambdaSEnum;
     881              else if (strcmp(name,"LATaylorHood")==0) return LATaylorHoodEnum;
     882              else if (strcmp(name,"LevelsetAnalysis")==0) return LevelsetAnalysisEnum;
    881883              else if (strcmp(name,"LevelsetfunctionPicard")==0) return LevelsetfunctionPicardEnum;
    882884              else if (strcmp(name,"LinearFloatingMeltRate")==0) return LinearFloatingMeltRateEnum;
    … …  
    996998              else if (strcmp(name,"Outputdefinition35")==0) return Outputdefinition35Enum;
    997999              else if (strcmp(name,"Outputdefinition36")==0) return Outputdefinition36Enum;
    998               else if (strcmp(name,"Outputdefinition37")==0) return Outputdefinition37Enum;
    999               else if (strcmp(name,"Outputdefinition38")==0) return Outputdefinition38Enum;
    10001000         else stage=9;
    10011001   }
    10021002   if(stage==9){
    1003               if (strcmp(name,"Outputdefinition39")==0) return Outputdefinition39Enum;
     1003              if (strcmp(name,"Outputdefinition37")==0) return Outputdefinition37Enum;
     1004              else if (strcmp(name,"Outputdefinition38")==0) return Outputdefinition38Enum;
     1005              else if (strcmp(name,"Outputdefinition39")==0) return Outputdefinition39Enum;
    10041006              else if (strcmp(name,"Outputdefinition3")==0) return Outputdefinition3Enum;
    10051007              else if (strcmp(name,"Outputdefinition40")==0) return Outputdefinition40Enum;
    … …  
    11191121              else if (strcmp(name,"SigmaVM")==0) return SigmaVMEnum;
    11201122              else if (strcmp(name,"SmbAnalysis")==0) return SmbAnalysisEnum;
    1121               else if (strcmp(name,"SMBcomponents")==0) return SMBcomponentsEnum;
    1122               else if (strcmp(name,"SMBd18opdd")==0) return SMBd18opddEnum;
    11231123         else stage=10;
    11241124   }
    11251125   if(stage==10){
    1126               if (strcmp(name,"SmbDesfac")==0) return SmbDesfacEnum;
     1126              if (strcmp(name,"SMBcomponents")==0) return SMBcomponentsEnum;
     1127              else if (strcmp(name,"SMBd18opdd")==0) return SMBd18opddEnum;
     1128              else if (strcmp(name,"SmbDesfac")==0) return SmbDesfacEnum;
    11271129              else if (strcmp(name,"SmbDpermil")==0) return SmbDpermilEnum;
    11281130              else if (strcmp(name,"SmbDzAdd")==0) return SmbDzAddEnum;

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

    r22840 r23468  
    3333               
    3434                Tmean = NaN; %mean annual temperature [K]
     35                Vmean = NaN; %mean annual wind velocity [m s-1]
    3536                C     = NaN; %mean annual snow accumulation [kg m-2 yr-1]
    3637                Tz    = NaN; %height above ground at which temperature (T) was sampled [m]
    … …  
    7273                                        % 7 = Greenland semi-emperical model of Kuipers Munneke et al. (2015)
    7374
     75                dsnowIdx = NaN; %model for fresh snow accumulation density (default is 1):
     76                         % 0 = Original GEMB value, 150 kg/m^3
     77                                        % 1 = Antarctica value of fresh snow density, 350 kg/m^3
     78                                        % 2 = Greenland value of fresh snow density, 315 kg/m^3, Fausto et al. (2008)
     79                                        % 3 = Antarctica model of Kaspers et al. (2004)
     80                                        % 4 = Greenland model of Kuipers Munneke et al. (2015)
     81
    7482                zTop  = NaN; % depth over which grid length is constant at the top of the snopack (default 10) [m]
    7583                dzTop = NaN; % initial top vertical grid spacing (default .05) [m]
    … …  
    155163                self.swIdx = 1;
    156164                self.denIdx = 2;
     165                self.dsnowIdx = 1;
    157166                self.zTop=10*ones(mesh.numberofelements,1);
    158167                self.dzTop = .05* ones (mesh.numberofelements,1);
    … …  
    160169                self.InitDensityScaling = 1.0;
    161170                self.ThermoDeltaTScaling = 1/11;
     171
     172                self.Vmean=10.0*ones(mesh.numberofelements,1);
    162173               
    163174                self.zMax=250*ones(mesh.numberofelements,1);
    … …  
    212223                        md = checkfield(md,'fieldname','smb.Tmean','size',[md.mesh.numberofelements 1],'NaN',1,'Inf',1,'>',273-100,'<',273+100); %-100/100 celsius min/max value
    213224                        md = checkfield(md,'fieldname','smb.C','size',[md.mesh.numberofelements 1],'NaN',1,'Inf',1,'>=',0);
     225                        md = checkfield(md,'fieldname','smb.Vmean','size',[md.mesh.numberofelements 1],'NaN',1,'Inf',1,'>=',0);
    214226                        md = checkfield(md,'fieldname','smb.Tz','size',[md.mesh.numberofelements 1],'NaN',1,'Inf',1,'>=',0,'<=',5000);
    215227                        md = checkfield(md,'fieldname','smb.Vz','size',[md.mesh.numberofelements 1],'NaN',1,'Inf',1,'>=',0,'<=',5000);
    … …  
    220232                        md = checkfield(md,'fieldname','smb.swIdx','NaN',1,'Inf',1,'values',[0,1]);
    221233                        md = checkfield(md,'fieldname','smb.denIdx','NaN',1,'Inf',1,'values',[1,2,3,4,5,6,7]);
     234                        md = checkfield(md,'fieldname','smb.dsnowIdx','NaN',1,'Inf',1,'values',[0,1,2,3,4]);
    222235
    223236                        md = checkfield(md,'fieldname','smb.zTop','NaN',1,'Inf',1,'>=',0);
    … …  
    265278                        fielddisplay(self,'isturbulentflux','run turbulant heat fluxes module (default true)');
    266279                        fielddisplay(self,'Ta','2 m air temperature, in Kelvin');
    267                         fielddisplay(self,'V','wind speed (m/s-1)');
     280                        fielddisplay(self,'V','wind speed (m s-1)');
    268281                        fielddisplay(self,'dlwrf','downward shortwave radiation flux [W/m^2]');
    269282                        fielddisplay(self,'dswrf','downward longwave radiation flux [W/m^2]');
    … …  
    273286                        fielddisplay(self,'Tmean','mean annual temperature [K]');
    274287                        fielddisplay(self,'C','mean annual snow accumulation [kg m-2 yr-1]');
     288                        fielddisplay(self,'Vmean','mean annual snow accumulation [m s-1] (default 10 m/s)');
    275289                        fielddisplay(self,'Tz','height above ground at which temperature (T) was sampled [m]');
    276290                        fielddisplay(self,'Vz','height above ground at which wind (V) eas sampled [m]');
    … …  
    331345                                                                        '6 = Antarctica semi-emperical model of Ligtenberg et al. (2011)',...
    332346                                                                        '7 = Greenland semi-emperical model of Kuipers Munneke et al. (2015)'});
     347                        fielddisplay(self,'dsnowIdx',{'model for fresh snow accumulation density (default is 1):',...
     348                                                           '0 = Original GEMB value, 150 kg/m^3',...
     349                                                           '1 = Antarctica value of fresh snow density, 350 kg/m^3',...
     350                                                           '2 = Greenland value of fresh snow density, 315 kg/m^3, Fausto et al. (2008)',...
     351                                                           '3 = Antarctica model of Kaspers et al. (2004), Make sure to set Vmean accurately',...
     352                                                           '4 = Greenland model of Kuipers Munneke et al. (2015)'});
     353
    333354                        fielddisplay(self,'requested_outputs','additional outputs requested');
    334355                                                                       
    … …  
    360381                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','Tmean','format','DoubleMat','mattype',2);
    361382                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','C','format','DoubleMat','mattype',2);
     383                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','Vmean','format','DoubleMat','mattype',2);
    362384                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','Tz','format','DoubleMat','mattype',2);
    363385                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','Vz','format','DoubleMat','mattype',2);
    … …  
    372394                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','swIdx','format','Integer');
    373395                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','denIdx','format','Integer');
     396                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','dsnowIdx','format','Integer');
     397                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','dsnowIdx','format','Integer');
    374398                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','InitDensityScaling','format','Double');
    375399                        WriteData(fid,prefix,'object',self,'class','smb','fieldname','ThermoDeltaTScaling','format','Double');

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

    r22840 r23468  
    3939               
    4040                Tmean = float('NaN')    #mean annual temperature [K]
     41                Vmean = float('NaN')    #mean annual wind velocity [m s-1]
    4142                C     = float('NaN')    #mean annual snow accumulation [kg m-2 yr-1]
    4243                Tz    = float('NaN')    #height above ground at which temperature (T) was sampled [m]
    … …  
    7879                                        # 6 = Antarctica semi-emperical model of Ligtenberg et al. (2011)
    7980                                        # 7 = Greenland semi-emperical model of Kuipers Munneke et al. (2015)
     81                                       
     82                dsnowIdx = float('NaN') #model for fresh snow accumulation density (default is 1):
     83                                        # 0 = Original GEMB value, 150 kg/m^3
     84                                        # 1 = Antarctica value of fresh snow density, 350 kg/m^3
     85                                        # 2 = Greenland value of fresh snow density, 315 kg/m^3, Fausto et al. (2008)
     86                                        # 3 = Antarctica model of Kaspers et al. (2004)
     87                                        # 4 = Greenland model of Kuipers Munneke et al. (2015)
    8088
    8189                zTop  = float('NaN')    # depth over which grid length is constant at the top of the snopack (default 10) [m]
    … …  
    131139                string = "%s\n%s"%(string,fielddisplay(self,'isturbulentflux','run turbulant heat fluxes module (default true)'))
    132140                string = "%s\n%s"%(string,fielddisplay(self,'Ta','2 m air temperature, in Kelvin'))
    133                 string = "%s\n%s"%(string,fielddisplay(self,'V','wind speed (m/s-1)'))
     141                string = "%s\n%s"%(string,fielddisplay(self,'V','wind speed (m s-1)'))
    134142                string = "%s\n%s"%(string,fielddisplay(self,'dlwrf','downward shortwave radiation flux [W/m^2]'))
    135143                string = "%s\n%s"%(string,fielddisplay(self,'dswrf','downward longwave radiation flux [W/m^2]'))
    … …  
    139147                string = "%s\n%s"%(string,fielddisplay(self,'Tmean','mean annual temperature [K]'))
    140148                string = "%s\n%s"%(string,fielddisplay(self,'C','mean annual snow accumulation [kg m-2 yr-1]'))
     149                string = "%s\n%s"%(string,fielddisplay(self,'Vmean','mean annual temperature [m s-1] (default 10 m/s)'))
    141150                string = "%s\n%s"%(string,fielddisplay(self,'Tz','height above ground at which temperature (T) was sampled [m]'))
    142151                string = "%s\n%s"%(string,fielddisplay(self,'Vz','height above ground at which wind (V) eas sampled [m]'))
    … …  
    196205                                                '6 = Antarctica semi-emperical model of Ligtenberg et al. (2011)',
    197206                                                '7 = Greenland semi-emperical model of Kuipers Munneke et al. (2015)']))
     207
     208                string = "%s\n%s"%(string,fielddisplay(self,'dsnowIdx',['model for fresh snow accumulation density (default is 1):',
     209                                                '0 = Original GEMB value, 150 kg/m^3',
     210                                                '1 = Antarctica value of fresh snow density, 350 kg/m^3',
     211                                                '2 = Greenland value of fresh snow density, 315 kg/m^3, Fausto et al. (2008)',
     212                                                '3 = Antarctica model of Kaspers et al. (2004), Make sure to set Vmean accurately',
     213                                                '4 = Greenland model of Kuipers Munneke et al. (2015)']));
     214
    198215                string = "%s\n%s"%(string,fielddisplay(self,'requested_outputs','additional outputs requested'))
    199216                return string
    … …  
    234251                self.swIdx = 1
    235252                self.denIdx = 2
     253                self.dsnowIdx = 1
    236254                self.zTop = 10*np.ones((mesh.numberofelements,))
    237255                self.dzTop = .05* np.ones((mesh.numberofelements,))
    … …  
    239257                self.InitDensityScaling = 1.0
    240258                self.ThermoDeltaTScaling = 1/11.0
     259
     260                self.Vmean = 10*np.ones((mesh.numberofelements,))
    241261               
    242262                self.zMax = 250*np.ones((mesh.numberofelements,))
    … …  
    291311                md = checkfield(md,'fieldname','smb.Tmean','size',[md.mesh.numberofelements],'NaN',1,'Inf',1,'>',273-100,'<',273+100) #-100/100 celsius min/max value
    292312                md = checkfield(md,'fieldname','smb.C','size',[md.mesh.numberofelements],'NaN',1,'Inf',1,'> = ',0)
     313                md = checkfield(md,'fieldname','smb.Vmean','size',[md.mesh.numberofelements],'NaN',1,'Inf',1,'> = ',0)
    293314                md = checkfield(md,'fieldname','smb.Tz','size',[md.mesh.numberofelements],'NaN',1,'Inf',1,'> = ',0,'< = ',5000)
    294315                md = checkfield(md,'fieldname','smb.Vz','size',[md.mesh.numberofelements],'NaN',1,'Inf',1,'> = ',0,'< = ',5000)
    … …  
    299320                md = checkfield(md,'fieldname','smb.swIdx','NaN',1,'Inf',1,'values',[0,1])
    300321                md = checkfield(md,'fieldname','smb.denIdx','NaN',1,'Inf',1,'values',[1,2,3,4,5,6,7])
     322                md = checkfield(md,'fieldname','smb.dsnowIdx','NaN',1,'Inf',1,'values',[0,1,2,3,4])
    301323
    302324                md = checkfield(md,'fieldname','smb.zTop','NaN',1,'Inf',1,'> = ',0)
    … …  
    356378                WriteData(fid,prefix,'object',self,'class','smb','fieldname','Tmean','format','DoubleMat','mattype',2)
    357379                WriteData(fid,prefix,'object',self,'class','smb','fieldname','C','format','DoubleMat','mattype',2)
     380                WriteData(fid,prefix,'object',self,'class','smb','fieldname','Vmean','format','DoubleMat','mattype',2)
    358381                WriteData(fid,prefix,'object',self,'class','smb','fieldname','Tz','format','DoubleMat','mattype',2)
    359382                WriteData(fid,prefix,'object',self,'class','smb','fieldname','Vz','format','DoubleMat','mattype',2)
    … …  
    368391                WriteData(fid,prefix,'object',self,'class','smb','fieldname','swIdx','format','Integer')
    369392                WriteData(fid,prefix,'object',self,'class','smb','fieldname','denIdx','format','Integer')
     393                WriteData(fid,prefix,'object',self,'class','smb','fieldname','dsnowIdx','format','Integer')
    370394                WriteData(fid,prefix,'object',self,'class','smb','fieldname','InitDensityScaling','format','Double')
    371395                WriteData(fid,prefix,'object',self,'class','smb','fieldname','ThermoDeltaTScaling','format','Double')

  • issm/trunk-jpl/test/NightlyRun/test243.m

    r22565 r23468  
    99% Use of Gemb method for SMB computation
    1010md.smb = SMBgemb(md.mesh,md.geometry);
     11md.smb.dsnowIdx = 0;
    1112
    1213%load hourly surface forcing date from 1979 to 2009:

  • issm/trunk-jpl/test/NightlyRun/test243.py

    r22869 r23468  
    2020md.smb = SMBgemb()
    2121md.smb.setdefaultparameters(md.mesh,md.geometry)
     22md.smb.dsnowIdx = 0
    2223
    2324#load hourly surface forcing date from 1979 to 2009:

  • issm/trunk-jpl/test/NightlyRun/test244.m

    r22429 r23468  
    1010% Use of Gemb method for SMB computation
    1111md.smb = SMBgemb(md.mesh,md.geometry);
     12md.smb.dsnowIdx = 0;
    1213
    1314%load hourly surface forcing date from 1979 to 2009:

  • issm/trunk-jpl/test/NightlyRun/test244.py

    r23130 r23468  
    3030md.smb = SMBgemb()
    3131md.smb.setdefaultparameters(md.mesh,md.geometry)
     32md.smb.dsnowIdx = 0
    3233
    3334#load hourly surface forcing date from 1979 to 2009:
Note: See TracChangeset for help on using the changeset viewer.