Changeset 22429


Ignore:
Timestamp:
02/16/18 11:12:12 (7 years ago)
Author:
schlegel
Message:

CHG: make dt thermo much smaller for diffusion stability

Location:
issm/trunk-jpl
Files:
6 edited

Legend:

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

    r22349 r22429  
    27862786
    27872787                /*Snow grain metamorphism:*/
    2788                 if(isgraingrowth)grainGrowth(re, gdn, gsp, T, dz, d, W, smb_dt, m, aIdx,this->Sid());
     2788                if(isgraingrowth)grainGrowth(&re, &gdn, &gsp, T, dz, d, W, smb_dt, m, aIdx,this->Sid());
    27892789
    27902790                /*Snow, firn and ice albedo:*/
    2791                 if(isalbedo)albedo(a,aIdx,re,d,cldFrac,aIce, aSnow,T,W,P,EC,t0wet,t0dry,K,smb_dt,rho_ice,m,this->Sid());
     2791                if(isalbedo)albedo(&a,aIdx,re,d,cldFrac,aIce, aSnow,T,W,P,EC,t0wet,t0dry,K,smb_dt,rho_ice,m,this->Sid());
    27922792
    27932793
     
    27992799
    28002800                /*Thermal profile computation:*/
    2801                 if(isthermal)thermo(&EC, T, dz, d, swf, dlw, Ta, V, eAir, pAir, W[0], smb_dt, m, Vz, Tz,rho_ice,this->Sid());
     2801                if(isthermal)thermo(&EC, &T, dz, d, swf, dlw, Ta, V, eAir, pAir, W[0], smb_dt, m, Vz, Tz,rho_ice,this->Sid());
    28022802
    28032803                /*Change in thickness of top cell due to evaporation/condensation  assuming same density as top cell.
     
    28132813
    28142814                /*Allow non-melt densification and determine compaction [m]*/
    2815                 if(isdensification)densification(d,dz, T, re, denIdx, C, smb_dt, Tmean,rho_ice,m,this->Sid());
     2815                if(isdensification)densification(&d,&dz, T, re, denIdx, C, smb_dt, Tmean,rho_ice,m,this->Sid());
    28162816
    28172817                /*Calculate upward longwave radiation flux [W m-2] not used in energy balance. Calculated for every
  • issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/Gembx.cpp

    r22249 r22429  
    147147
    148148} /*}}}*/
    149 void grainGrowth(IssmDouble* re, IssmDouble* gdn, IssmDouble* gsp, IssmDouble* T,IssmDouble* dz,IssmDouble* d, IssmDouble* W,IssmDouble smb_dt,int m,int aIdx,int sid){ /*{{{*/
     149void grainGrowth(IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, IssmDouble* T,IssmDouble* dz,IssmDouble* d, IssmDouble* W,IssmDouble smb_dt,int m,int aIdx,int sid){ /*{{{*/
    150150
    151151        /*Created by: Alex S. Gardner, University of Alberta
     
    191191        IssmDouble  Q=0.0;
    192192
     193        /*output: */
     194        IssmDouble* re=NULL;
     195        IssmDouble* gdn=NULL;
     196        IssmDouble* gsp=NULL;
     197
    193198        if(VerboseSmb() && sid==0 && IssmComm::GetRank()==0)_printf0_("   grain growth module\n");
     199
     200        /*Recover pointers: */
     201        re=*pre;
     202        gdn=*pgdn;
     203        gsp=*pgsp;
    194204       
    195205        /*only when aIdx = 1 or 2 do we run grainGrowth: */
     
    321331        xDelete<IssmDouble>(lwc);
    322332
     333        /*Assign output pointers:*/
     334        *pre=re;
     335        *pgdn=gdn;
     336        *pgsp=gsp;
     337
    323338}  /*}}}*/
    324 void albedo(IssmDouble* a, 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* TK, IssmDouble* W, IssmDouble P, IssmDouble EC, IssmDouble t0wet, IssmDouble t0dry, IssmDouble K, IssmDouble dt, IssmDouble dIce, int m,int sid) { /*{{{*/
    325340
    326341        //// Calculates Snow, firn and ice albedo as a function of:
     
    366381        //   [273 272.5 ... 265], [0 0.001 ... 0], 0, 0.01, 15, 15, 7, 3600)
    367382
     383        /*output: */
     384        IssmDouble* a=NULL;
     385
    368386        if(VerboseSmb() && sid==0 && IssmComm::GetRank()==0)_printf0_("   albedo module\n");
     387
     388        /*Recover pointers: */
     389        a=*pa;
    369390       
    370391        //some constants:
     
    471492        else if (a[0] < 0) _printf_("albedo is negative\n");
    472493        else if (xIsNan(a[0])) _error_("albedo == NAN\n");
     494
     495        /*Assign output pointers:*/
     496        *pa=a;
     497
    473498}  /*}}}*/
    474 void thermo(IssmDouble* pEC, IssmDouble* T, 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) { /*{{{*/
     499void 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) { /*{{{*/
    475500
    476501        /* ENGLACIAL THERMODYNAMICS*/
     
    544569        /*outputs:*/
    545570        IssmDouble EC=0.0;
     571        IssmDouble* T=NULL;
    546572
    547573        if(VerboseSmb() && sid==0 && IssmComm::GetRank()==0)_printf0_("   thermal module\n");
     574
     575        /*Recover pointers: */
     576        T=*pT;
    548577
    549578        ds = d[0];      // density of top grid cell
     
    626655
    627656        // determine minimum acceptable delta t (diffusion number > 1/2) [s]
     657        // NS: 2.16.18 divided dt by 9 for stability (changed 3 to 27)
    628658        dt=1e12;
    629         for(int i=0;i<m;i++)dt = fmin(dt,CI * pow(dz[i],2) * d[i]  / (3 * K[i]));
     659        for(int i=0;i<m;i++)dt = fmin(dt,CI * pow(dz[i],2) * d[i]  / (27 * K[i]));
    630660
    631661        // smallest possible even integer of 60 min where diffusion number > 1/2
     
    833863        /*Assign output pointers:*/
    834864        *pEC=EC;
     865        *pT=T;
    835866
    836867}  /*}}}*/
     
    17541785
    17551786} /*}}}*/
    1756 void densification(IssmDouble* d,IssmDouble* dz, IssmDouble* T, IssmDouble* re, int denIdx, IssmDouble C, IssmDouble dt, IssmDouble Tmean, IssmDouble dIce, int m, int sid){ /*{{{*/
     1787void densification(IssmDouble** pd,IssmDouble** pdz, IssmDouble* T, IssmDouble* re, int denIdx, IssmDouble C, IssmDouble dt, IssmDouble Tmean, IssmDouble dIce, int m, int sid){ /*{{{*/
    17571788
    17581789        //// 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???]
     
    18061837        IssmDouble c1=0.0;
    18071838        IssmDouble H=0.0;
     1839
     1840        /*output: */
     1841        IssmDouble* dz=NULL;
     1842        IssmDouble* d=NULL;
    18081843       
    18091844        if(VerboseSmb() && sid==0 && IssmComm::GetRank()==0)_printf0_("   densification module\n");
     1845
     1846        /*Recover pointers: */
     1847        dz=*pdz;
     1848        d=*pd;
    18101849
    18111850        // initial mass
     
    18781917        xDelete<IssmDouble>(obp);
    18791918
     1919        /*Assign output pointers:*/
     1920        *pdz=dz;
     1921        *pd=d;
     1922
    18801923} /*}}}*/
    18811924void turbulentFlux(IssmDouble* pshf, IssmDouble* plhf, IssmDouble* pEC, IssmDouble Ta, IssmDouble Ts, IssmDouble V, IssmDouble eAir, IssmDouble pAir, IssmDouble ds, IssmDouble Ws, IssmDouble Vz, IssmDouble Tz, IssmDouble dIce, int sid){ /*{{{*/
  • issm/trunk-jpl/src/c/modules/SurfaceMassBalancex/SurfaceMassBalancex.h

    r22249 r22429  
    2424void       GembgridInitialize(IssmDouble** pdz, int* psize, IssmDouble zTop, IssmDouble dzTop, IssmDouble zMax, IssmDouble zY);
    2525IssmDouble Marbouty(IssmDouble T, IssmDouble d, IssmDouble dT);
    26 void 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);
     26void grainGrowth(IssmDouble** pre, IssmDouble** pgdn, IssmDouble** pgsp, IssmDouble* T,IssmDouble* dz,IssmDouble* d, IssmDouble* W,IssmDouble smb_dt,int m,int aIdx, int sid);
     27void 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);
    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);
     29void 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);
    3030void 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);
    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);
    32 void densification(IssmDouble* d,IssmDouble* dz, IssmDouble* T, IssmDouble* re, int denIdx, IssmDouble C, IssmDouble dt, IssmDouble Tmean, IssmDouble dIce, int m, int sid);
     32void densification(IssmDouble** pd,IssmDouble** pdz, IssmDouble* T, IssmDouble* re, int denIdx, IssmDouble C, IssmDouble dt, IssmDouble Tmean, IssmDouble dIce, int m, int sid);
    3333void turbulentFlux(IssmDouble* pshf, IssmDouble* plhf, IssmDouble* pEC, IssmDouble Ta, IssmDouble Ts, IssmDouble V, IssmDouble eAir, IssmDouble pAir, IssmDouble ds, IssmDouble Ws, IssmDouble Vz, IssmDouble Tz, IssmDouble dIce, int sid);
    3434#endif  /* _SurfaceMassBalancex_H*/
  • issm/trunk-jpl/test/NightlyRun/test244.m

    r22250 r22429  
    7373%  nond_sampling study
    7474md.qmu.method=dakota_method('nond_samp');
    75 md.qmu.method(end)=dmeth_params_set(md.qmu.method(end),'seed',1234,'samples',4,'sample_type','lhs');
     75md.qmu.method(end)=dmeth_params_set(md.qmu.method(end),'seed',1234,'samples',3,'sample_type','lhs');
    7676dver=textscan(IssmConfig('_DAKOTA_VERSION_'),'%[0123456789].%[0123456789].%[0123456789]');
    7777if ((str2num(dver{1}{1})==4 && str2num(dver{2}{1})>2) || str2num(dver{1}{1})>4)
Note: See TracChangeset for help on using the changeset viewer.