Changeset 27936


Ignore:
Timestamp:
10/03/23 09:26:58 (19 months ago)
Author:
koitf1
Message:

CHG:islaminar flag renamed to istransition, hydrology sheet alpha/beta and omega changed from input type to parameter type

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

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.cpp

    r27927 r27936  
    146146        iomodel->FetchDataToInput(inputs,elements,"md.hydrology.bump_height",HydrologyBumpHeightEnum);
    147147        iomodel->FetchDataToInput(inputs,elements,"md.hydrology.sheet_conductivity",HydrologySheetConductivityEnum);
    148         iomodel->FetchDataToInput(inputs,elements,"md.hydrology.omega",HydrologyOmegaEnum);
    149         iomodel->FetchDataToInput(inputs,elements,"md.hydrology.sheet_alpha",HydrologySheetAlphaEnum);
    150         iomodel->FetchDataToInput(inputs,elements,"md.hydrology.sheet_beta",HydrologySheetBetaEnum);
    151148        iomodel->FetchDataToInput(inputs,elements,"md.hydrology.channel_conductivity",HydrologyChannelConductivityEnum);
    152149        iomodel->FetchDataToInput(inputs,elements,"md.hydrology.neumannflux",HydrologyNeumannfluxEnum);
     
    188185        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.channel_alpha",HydrologyChannelAlphaEnum));
    189186        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.channel_beta",HydrologyChannelBetaEnum));
    190         parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.islaminar",HydrologyIsLaminarEnum));
     187        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.sheet_alpha",HydrologySheetAlphaEnum));
     188        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.sheet_beta",HydrologySheetBetaEnum));
     189        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.omega",HydrologyOmegaEnum));
     190        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.istransition",HydrologyIsTransitionEnum));
    191191        parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.englacial_void_ratio",HydrologyEnglacialVoidRatioEnum));
    192192
     
    222222        /*Intermediaries */
    223223        IssmDouble  Jdet,dphi[3],h,k;
    224         IssmDouble  alpha,beta,omega,h_r;
     224        IssmDouble  h_r;
    225225        IssmDouble  A,B,n,phi_old,phi,phi_0,H,b,v1;
    226226        IssmDouble* xyz_list = NULL;
    227 
    228         /*Hard coded coefficients*/
    229         const IssmPDouble ALPHA = 5./4.;
    230         const IssmPDouble BETA  = 3./2.;
    231227
    232228        /*Fetch number of nodes and dof for this finite element*/
     
    242238
    243239        /*Get all inputs and parameters*/
    244         bool islaminar;
    245         element->FindParam(&islaminar,HydrologyIsLaminarEnum);
     240        bool istransition;
     241        element->FindParam(&istransition,HydrologyIsTransitionEnum);
     242        IssmDouble alpha     = element->FindParam(HydrologySheetAlphaEnum);
     243        IssmDouble beta      = element->FindParam(HydrologySheetBetaEnum);
     244        IssmDouble omega     = element->FindParam(HydrologyOmegaEnum);
    246245        IssmDouble dt        = element->FindParam(TimesteppingTimeStepEnum);
    247246        IssmDouble rho_water = element->FindParam(MaterialsRhoFreshwaterEnum);
     
    252251        Input* hr_input      = element->GetInput(HydrologyBumpHeightEnum);_assert_(hr_input);
    253252        Input* k_input   = element->GetInput(HydrologySheetConductivityEnum);_assert_(k_input);
    254         Input* alpha_input   = element->GetInput(HydrologySheetAlphaEnum);_assert_(alpha_input);
    255         Input* beta_input    = element->GetInput(HydrologySheetBetaEnum);_assert_(beta_input);
    256         Input* omega_input   = element->GetInput(HydrologyOmegaEnum);_assert_(omega_input);
    257253        Input* phi_input = element->GetInput(HydraulicPotentialEnum);      _assert_(phi_input);
    258254        Input* h_input   = element->GetInput(HydrologySheetThicknessEnum); _assert_(h_input);
     
    274270                h_input->GetInputValue(&h,gauss);
    275271                k_input->GetInputValue(&k,gauss);
    276                 omega_input->GetInputValue(&omega,gauss);
    277                 alpha_input->GetInputValue(&alpha,gauss);
    278                 beta_input->GetInputValue(&beta,gauss);
    279272                B_input->GetInputValue(&B,gauss);
    280273                n_input->GetInputValue(&n,gauss);
     
    290283                if(normgradphi < AEPS) normgradphi = AEPS;
    291284               
    292                 /*Use Laminar flow if specified*/
     285                /*Use transition model if specified*/
    293286                IssmDouble nu = mu_water/rho_water;
    294287                IssmDouble coeff;
    295                 if(islaminar==1 && omega>=AEPS){
     288                if(istransition==1 && omega>=AEPS){
    296289                        IssmDouble hratio = fabs(h/h_r);
    297290                        IssmDouble coarg = 1. + 4.*pow(hratio,3-2*alpha)*omega*k*pow(h,3)*normgradphi/nu;
     
    300293                else {
    301294                        /*If omega is zero, use standard model, otherwise transition model*/
    302                         coeff = k*pow(h,ALPHA)*pow(normgradphi,BETA-2.);
     295                        coeff = k*pow(h,alpha)*pow(normgradphi,beta-2.);
    303296                }
    304297
     
    468461        /*Intermediaries*/
    469462   IssmDouble  dphi[3],h,k,phi;
    470         IssmDouble  alpha,beta,omega,h_r;
     463        IssmDouble  h_r;
    471464        IssmDouble  oceanLS,iceLS;
    472465        IssmDouble* xyz_list = NULL;
    473 
    474         /*Hard coded coefficients*/
    475         const IssmPDouble ALPHA = 5./4.;
    476         const IssmPDouble BETA  = 3./2.;
    477466
    478467        /*Fetch number vertices for this element*/
     
    495484
    496485        /*Retrieve all inputs and parameters*/
    497         bool islaminar;
    498         element->FindParam(&islaminar,HydrologyIsLaminarEnum);
     486        bool istransition;
     487        element->FindParam(&istransition,HydrologyIsTransitionEnum);
     488        IssmDouble alpha     = element->FindParam(HydrologySheetAlphaEnum);
     489        IssmDouble beta      = element->FindParam(HydrologySheetBetaEnum);
     490        IssmDouble omega     = element->FindParam(HydrologyOmegaEnum);
    499491        element->GetVerticesCoordinates(&xyz_list);
    500492        IssmDouble rho_water = element->FindParam(MaterialsRhoFreshwaterEnum);
    501493        IssmDouble mu_water  = element->FindParam(MaterialsMuWaterEnum);
    502494        Input *k_input       = element->GetInput(HydrologySheetConductivityEnum); _assert_(k_input);
    503         Input *omega_input   = element->GetInput(HydrologyOmegaEnum);             _assert_(omega_input);
    504         Input *alpha_input   = element->GetInput(HydrologySheetAlphaEnum);        _assert_(alpha_input);
    505         Input *beta_input    = element->GetInput(HydrologySheetBetaEnum);         _assert_(beta_input);
    506495        Input *phi_input     = element->GetInput(HydraulicPotentialEnum);         _assert_(phi_input);
    507496        Input *hr_input      = element->GetInput(HydrologyBumpHeightEnum);        _assert_(hr_input);
     
    521510        hr_input->GetInputValue(&h_r,gauss);
    522511      k_input->GetInputValue(&k,gauss);
    523         omega_input->GetInputValue(&omega,gauss);
    524         alpha_input->GetInputValue(&alpha,gauss);
    525         beta_input->GetInputValue(&beta,gauss);
    526512                oceanLS_input->GetInputValue(&oceanLS,gauss);
    527513                iceLS_input->GetInputValue(&iceLS,gauss);
     
    541527         IssmDouble nu = mu_water/rho_water;
    542528        IssmDouble coeff;
    543         if(islaminar==1 && omega>=AEPS){
     529        if(istransition==1 && omega>=AEPS){
    544530                IssmDouble hratio = fabs(h/h_r);
    545531                IssmDouble coarg = 1. + 4.*pow(hratio,3-2*alpha)*omega*k*pow(h,3)*normgradphi/nu;
     
    548534        }
    549535        else {
    550                 coeff = k*pow(h,ALPHA)*pow(normgradphi,BETA-2.)/max(AEPS,h);  // divide by h to get speed instead of discharge
     536                coeff = k*pow(h,alpha)*pow(normgradphi,beta-2.)/max(AEPS,h);  // divide by h to get speed instead of discharge
    551537        }
    552538
  • issm/trunk-jpl/src/c/classes/Loads/Channel.cpp

    r27926 r27936  
    2020
    2121#define C_W         4.22e3   /*specific heat capacity of water (J/kg/K)*/
    22 #define ALPHA_C     5./4.
    23 #define BETA_C      3./2.
    24 /*Make sure these are the same as in HydrologyGlaDSAnalysis::CreateKMatrix*/
    25 #define ALPHA_S     5./4.
    26 #define BETA_S      3./2.
    2722#define AEPS        2.2204460492503131E-015
    2823
     
    367362        IssmDouble  Jdet,v1,qc,fFactor,Afactor,Bfactor,Xifactor;
    368363        IssmDouble  A,B,n,phi_old,phi,phi_0,dPw,ks,kc,Ngrad;
    369         IssmDouble  alpha_s,beta_s,omega,h_r;
     364        IssmDouble  h_r;
    370365        IssmDouble  H,h,b,dphi[2],dphids,dphimds,db[2],dbds;
    371366        IssmDouble  xyz_list[NUMVERTICES][3];
     
    383378        GetVerticesCoordinates(&xyz_list_tria[0][0],tria->vertices,3);
    384379
    385         bool islaminar;
    386         element->FindParam(&islaminar,HydrologyIsLaminarEnum);
     380        bool istransition;
     381        element->FindParam(&istransition,HydrologyIsTransitionEnum);
    387382        IssmDouble L         = element->FindParam(MaterialsLatentheatEnum);
    388383        IssmDouble mu_water  = element->FindParam(MaterialsMuWaterEnum);
     
    394389        IssmDouble alpha_c   = element->FindParam(HydrologyChannelAlphaEnum);
    395390        IssmDouble beta_c    = element->FindParam(HydrologyChannelBetaEnum);
     391        IssmDouble alpha_s   = element->FindParam(HydrologySheetAlphaEnum);
     392        IssmDouble beta_s    = element->FindParam(HydrologySheetBetaEnum);
     393        IssmDouble omega     = element->FindParam(HydrologyOmegaEnum);
    396394
    397395        Input* h_input      = element->GetInput(HydrologySheetThicknessEnum);      _assert_(h_input);
     
    402400        Input* ks_input     = element->GetInput(HydrologySheetConductivityEnum);   _assert_(ks_input);
    403401        Input* kc_input     = element->GetInput(HydrologyChannelConductivityEnum); _assert_(kc_input);
    404         Input* omega_input  = element->GetInput(HydrologyOmegaEnum); _assert_(ks_input);
    405         Input* alpha_s_input= element->GetInput(HydrologySheetAlphaEnum); _assert_(alpha_s_input);
    406         Input* beta_s_input = element->GetInput(HydrologySheetBetaEnum); _assert_(beta_s_input);
    407402        Input* hr_input     = element->GetInput(HydrologyBumpHeightEnum);_assert_(hr_input);
    408403        Input* phi_input    = element->GetInput(HydraulicPotentialEnum);           _assert_(phi_input);
     
    432427                ks_input->GetInputValue(&ks,gauss);
    433428                kc_input->GetInputValue(&kc,gauss);
    434                 alpha_s_input->GetInputValue(&alpha_s,gauss);
    435                 beta_s_input->GetInputValue(&beta_s,gauss);
    436                 omega_input->GetInputValue(&omega,gauss);
    437429                hr_input->GetInputValue(&h_r,gauss);
    438430                B_input->GetInputValue(&B,gauss);
     
    451443                if(Ngrad<AEPS) Ngrad = AEPS;
    452444               
    453                 /*Compute the effective conductivity Kc = k h^alpha |grad Phi|^{beta-2} (same for sheet) and use laminar if specified*/
     445                /*Compute the effective conductivity Kc = k h^alpha |grad Phi|^{beta-2} (same for sheet) and use transition model if specified*/
    454446                IssmDouble Kc;
    455447                IssmDouble Ks;
    456448                IssmDouble nu = mu_water/rho_water;
    457                 if(islaminar==1 && omega>=AEPS){
     449                if(istransition==1 && omega>=AEPS){
    458450                        IssmDouble hratio = h/h_r;
    459451                        IssmDouble coarg = 1. + 4.*omega*pow(hratio,3-2*alpha_s)*ks*pow(h,3)*Ngrad/nu;
     
    462454                }
    463455                else {
    464                         Ks = ks*pow(h,ALPHA_S)*pow(Ngrad,BETA_S-2.);
    465                         Kc = kc * pow(this->S,ALPHA_C) * pow(Ngrad,BETA_C-2.);
     456                        Ks = ks*pow(h,alpha_s)*pow(Ngrad,beta_s-2.);
     457                        Kc = kc * pow(this->S,alpha_c) * pow(Ngrad,beta_c-2.);
    466458                }
    467459
     
    529521        IssmDouble  A,B,n,phi_old,phi,phi_0,dphimds,dphi[2];
    530522        IssmDouble  H,h,b,db[2],dphids,qc,dPw,ks,kc,Ngrad;
    531         IssmDouble  alpha_s,beta_s,omega,h_r;
     523        IssmDouble  h_r;
    532524        IssmDouble  xyz_list[NUMVERTICES][3];
    533525        IssmDouble  xyz_list_tria[3][3];
     
    542534        GetVerticesCoordinates(&xyz_list_tria[0][0],tria->vertices,3);
    543535
    544         bool islaminar;
    545         element->FindParam(&islaminar,HydrologyIsLaminarEnum);
     536        bool istransition;
     537        element->FindParam(&istransition,HydrologyIsTransitionEnum);
    546538        IssmDouble L         = element->FindParam(MaterialsLatentheatEnum);
    547539        IssmDouble mu_water  = element->FindParam(MaterialsMuWaterEnum);
     
    551543        IssmDouble lc        = element->FindParam(HydrologyChannelSheetWidthEnum);
    552544        IssmDouble c_t       = element->FindParam(HydrologyPressureMeltCoefficientEnum);
     545        IssmDouble alpha_s   = element->FindParam(HydrologySheetAlphaEnum);
     546        IssmDouble beta_s    = element->FindParam(HydrologySheetBetaEnum);
     547        IssmDouble omega     = element->FindParam(HydrologyOmegaEnum);
    553548
    554549        Input* h_input      = element->GetInput(HydrologySheetThicknessEnum);      _assert_(h_input);
     
    559554        Input* ks_input     = element->GetInput(HydrologySheetConductivityEnum);   _assert_(ks_input);
    560555        Input* kc_input     = element->GetInput(HydrologyChannelConductivityEnum); _assert_(kc_input);
    561         Input* omega_input  = element->GetInput(HydrologyOmegaEnum); _assert_(omega_input);
    562556        Input* phi_input    = element->GetInput(HydraulicPotentialEnum);           _assert_(phi_input);
    563         Input* alpha_s_input= element->GetInput(HydrologySheetAlphaEnum);_assert_(alpha_s_input);
    564         Input* beta_s_input = element->GetInput(HydrologySheetBetaEnum);_assert_(beta_s_input);
    565557        Input* hr_input     = element->GetInput(HydrologyBumpHeightEnum);_assert_(hr_input);
    566558
     
    585577                ks_input->GetInputValue(&ks,gauss);
    586578                kc_input->GetInputValue(&kc,gauss);
    587                 omega_input->GetInputValue(&omega,gauss);
    588579                B_input->GetInputValue(&B,gauss);
    589580                n_input->GetInputValue(&n,gauss);
     
    591582                b_input->GetInputValue(&b,gauss);
    592583                H_input->GetInputValue(&H,gauss);
    593                 alpha_s_input->GetInputValue(&alpha_s,gauss);
    594                 beta_s_input->GetInputValue(&beta_s,gauss);
    595584                hr_input->GetInputValue(&h_r,gauss);
    596585
     
    606595
    607596
    608                 /*Approx. discharge in the sheet flowing folwing in the direction of the channel ofver a width lc, use laminar if specified*/
     597                /*Approx. discharge in the sheet flowing folwing in the direction of the channel ofver a width lc, use transition model if specified*/
    609598                IssmDouble Ks;
    610                 if (islaminar==1 && omega>=AEPS){
     599                if (istransition==1 && omega>=AEPS){
    611600                IssmDouble hratio = h/h_r;
    612601                        IssmDouble nu = mu_water/rho_water;
     
    615604                }
    616605                else {
    617                         Ks = ks * pow(h,ALPHA_S) * pow(Ngrad,BETA_S-2.);
     606                        Ks = ks * pow(h,alpha_s) * pow(Ngrad,beta_s-2.);
    618607                }
    619608
     
    683672        /*Intermediaries */
    684673        IssmDouble  A,B,n,phi,phi_0,ks,kc,Ngrad;
    685         IssmDouble  alpha_s,beta_s,omega,h_r;
     674        IssmDouble  h_r;
    686675        IssmDouble  H,h,b,dphi[2],dphids,dphimds,db[2],dbds;
    687676        IssmDouble  xyz_list[NUMVERTICES][3];
     
    692681        GetVerticesCoordinates(&xyz_list_tria[0][0],tria->vertices,3);
    693682
    694         bool islaminar;
    695         element->FindParam(&islaminar,HydrologyIsLaminarEnum);
     683        bool istransition;
     684        element->FindParam(&istransition,HydrologyIsTransitionEnum);
    696685        IssmDouble L         = element->FindParam(MaterialsLatentheatEnum);
    697686        IssmDouble rho_ice   = element->FindParam(MaterialsRhoIceEnum);
     
    704693        IssmDouble alpha_c   = element->FindParam(HydrologyChannelAlphaEnum);
    705694        IssmDouble beta_c    = element->FindParam(HydrologyChannelBetaEnum);
     695        IssmDouble alpha_s   = element->FindParam(HydrologySheetAlphaEnum);
     696        IssmDouble beta_s    = element->FindParam(HydrologySheetBetaEnum);
     697        IssmDouble omega     = element->FindParam(HydrologyOmegaEnum);
    706698
    707699        Input* h_input      = element->GetInput(HydrologySheetThicknessEnum);      _assert_(h_input);
     
    712704        Input* ks_input     = element->GetInput(HydrologySheetConductivityEnum);   _assert_(ks_input);
    713705        Input* kc_input     = element->GetInput(HydrologyChannelConductivityEnum); _assert_(kc_input);
    714         Input* omega_input  = element->GetInput(HydrologyOmegaEnum); _assert_(omega_input);
    715706        Input* phi_input    = element->GetInput(HydraulicPotentialEnum);           _assert_(phi_input);
    716         Input* alpha_s_input= element->GetInput(HydrologySheetAlphaEnum);_assert_(alpha_s_input);
    717         Input* beta_s_input = element->GetInput(HydrologySheetBetaEnum); _assert_(beta_s_input);
    718707        Input* hr_input     = element->GetInput(HydrologyBumpHeightEnum);_assert_(hr_input);
    719708
     
    731720        ks_input->GetInputValue(&ks,gauss);
    732721        kc_input->GetInputValue(&kc,gauss);
    733         omega_input->GetInputValue(&omega,gauss);
    734722        B_input->GetInputValue(&B,gauss);
    735723        n_input->GetInputValue(&n,gauss);
     
    738726        H_input->GetInputValue(&H,gauss);
    739727        hr_input->GetInputValue(&h_r,gauss);
    740         alpha_s_input->GetInputValue(&alpha_s,gauss);
    741         beta_s_input->GetInputValue(&beta_s,gauss);
    742728
    743729
     
    755741        IssmDouble dPw = dphids - dphimds;
    756742
    757         /*Approx. discharge in the sheet flowing folwing in the direction of the channel ofver a width lc, use laminar if necessary*/
     743        /*Approx. discharge in the sheet flowing folwing in the direction of the channel ofver a width lc, use transition model if necessary*/
    758744        IssmDouble qc;
    759         if (islaminar==1 && omega>=AEPS){
     745        if (istransition==1 && omega>=AEPS){
    760746        IssmDouble hratio = h/h_r;
    761747                IssmDouble nu = mu_water/rho_water;
     
    764750        }
    765751        else {
    766                 qc = - ks * pow(h,ALPHA_S) * pow(Ngrad,BETA_S-2.) * dphids;
     752                qc = - ks * pow(h,alpha_s) * pow(Ngrad,beta_s-2.) * dphids;
    767753        }
    768754
     
    771757
    772758        IssmDouble C = C_W*c_t*rho_water;
    773         IssmDouble Qprime = -kc * pow(Ngrad,BETA_C-2.)*dphids;
     759        IssmDouble Qprime = -kc * pow(Ngrad,beta_c-2.)*dphids;
    774760        IssmDouble N = phi_0 - phi;
    775761
     
    788774
    789775                IssmDouble alpha = 1./(rho_ice*L)*(
    790                                         fabs(Qprime*pow(Snew,ALPHA_C-1.)*dphids)
    791                                         + C*Qprime*pow(Snew,ALPHA_C-1.)*dPw
     776                                        fabs(Qprime*pow(Snew,alpha_c-1.)*dphids)
     777                                        + C*Qprime*pow(Snew,alpha_c-1.)*dPw
    792778                                        ) - 2./pow(n,n)*A*pow(fabs(N),n-1.)*N;
    793779
     
    811797
    812798        /*Compute new channel discharge for output only*/
    813         IssmDouble Kc = kc * pow(this->S,ALPHA_C) * pow(Ngrad,BETA_C-2.);
     799        IssmDouble Kc = kc * pow(this->S,alpha_c) * pow(Ngrad,beta_c-2.);
    814800        this->discharge = -Kc*dphids;
    815801
  • issm/trunk-jpl/src/c/shared/Enum/Enum.vim

    r27925 r27936  
    267267syn keyword cConstant HydrologyEnglacialVoidRatioEnum
    268268syn keyword cConstant HydrologyIschannelsEnum
    269 syn keyword cConstant HydrologyIsLaminarEnum
     269syn keyword cConstant HydrologyIsTransitionEnum
    270270syn keyword cConstant HydrologyIsWaterPressureArmaEnum
    271271syn keyword cConstant HydrologyMeltFlagEnum
     
    273273syn keyword cConstant HydrologyNumBasinsEnum
    274274syn keyword cConstant HydrologyNumRequestedOutputsEnum
     275syn keyword cConstant HydrologyOmegaEnum
    275276syn keyword cConstant HydrologyPressureMeltCoefficientEnum
    276277syn keyword cConstant HydrologyRelaxationEnum
    277278syn keyword cConstant HydrologyRequestedOutputsEnum
    278279syn keyword cConstant HydrologySedimentKmaxEnum
     280syn keyword cConstant HydrologySheetAlphaEnum
     281syn keyword cConstant HydrologySheetBetaEnum
    279282syn keyword cConstant HydrologyStepsPerStepEnum
    280283syn keyword cConstant HydrologyStorageEnum
     
    898901syn keyword cConstant HydrologyMoulinInputEnum
    899902syn keyword cConstant HydrologyNeumannfluxEnum
    900 syn keyword cConstant HydrologyOmegaEnum
    901903syn keyword cConstant HydrologyReynoldsEnum
    902 syn keyword cConstant HydrologySheetAlphaEnum
    903 syn keyword cConstant HydrologySheetBetaEnum
    904904syn keyword cConstant HydrologySheetConductivityEnum
    905905syn keyword cConstant HydrologySheetThicknessEnum
  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r27925 r27936  
    261261        HydrologyEnglacialVoidRatioEnum,
    262262        HydrologyIschannelsEnum,
    263         HydrologyIsLaminarEnum,
     263        HydrologyIsTransitionEnum,
    264264        HydrologyIsWaterPressureArmaEnum,
    265265        HydrologyMeltFlagEnum,
     
    267267        HydrologyNumBasinsEnum,
    268268        HydrologyNumRequestedOutputsEnum,
     269        HydrologyOmegaEnum,
    269270        HydrologyPressureMeltCoefficientEnum,
    270271        HydrologyRelaxationEnum,
    271272        HydrologyRequestedOutputsEnum,
    272273        HydrologySedimentKmaxEnum,
     274        HydrologySheetAlphaEnum,
     275        HydrologySheetBetaEnum,
    273276        HydrologyStepsPerStepEnum,
    274277        HydrologyStorageEnum,
     
    894897        HydrologyMoulinInputEnum,
    895898        HydrologyNeumannfluxEnum,
    896         HydrologyOmegaEnum,
    897899        HydrologyReynoldsEnum,
    898         HydrologySheetAlphaEnum,
    899         HydrologySheetBetaEnum,
    900900        HydrologySheetConductivityEnum,
    901901        HydrologySheetThicknessEnum,
  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r27925 r27936  
    269269                case HydrologyEnglacialVoidRatioEnum : return "HydrologyEnglacialVoidRatio";
    270270                case HydrologyIschannelsEnum : return "HydrologyIschannels";
    271                 case HydrologyIsLaminarEnum : return "HydrologyIsLaminar";
     271                case HydrologyIsTransitionEnum : return "HydrologyIsTransition";
    272272                case HydrologyIsWaterPressureArmaEnum : return "HydrologyIsWaterPressureArma";
    273273                case HydrologyMeltFlagEnum : return "HydrologyMeltFlag";
     
    275275                case HydrologyNumBasinsEnum : return "HydrologyNumBasins";
    276276                case HydrologyNumRequestedOutputsEnum : return "HydrologyNumRequestedOutputs";
     277                case HydrologyOmegaEnum : return "HydrologyOmega";
    277278                case HydrologyPressureMeltCoefficientEnum : return "HydrologyPressureMeltCoefficient";
    278279                case HydrologyRelaxationEnum : return "HydrologyRelaxation";
    279280                case HydrologyRequestedOutputsEnum : return "HydrologyRequestedOutputs";
    280281                case HydrologySedimentKmaxEnum : return "HydrologySedimentKmax";
     282                case HydrologySheetAlphaEnum : return "HydrologySheetAlpha";
     283                case HydrologySheetBetaEnum : return "HydrologySheetBeta";
    281284                case HydrologyStepsPerStepEnum : return "HydrologyStepsPerStep";
    282285                case HydrologyStorageEnum : return "HydrologyStorage";
     
    900903                case HydrologyMoulinInputEnum : return "HydrologyMoulinInput";
    901904                case HydrologyNeumannfluxEnum : return "HydrologyNeumannflux";
    902                 case HydrologyOmegaEnum : return "HydrologyOmega";
    903905                case HydrologyReynoldsEnum : return "HydrologyReynolds";
    904                 case HydrologySheetAlphaEnum : return "HydrologySheetAlpha";
    905                 case HydrologySheetBetaEnum : return "HydrologySheetBeta";
    906906                case HydrologySheetConductivityEnum : return "HydrologySheetConductivity";
    907907                case HydrologySheetThicknessEnum : return "HydrologySheetThickness";
  • issm/trunk-jpl/src/c/shared/Enum/Enumjl.vim

    r27925 r27936  
    260260syn keyword juliaConstC HydrologyEnglacialVoidRatioEnum
    261261syn keyword juliaConstC HydrologyIschannelsEnum
    262 syn keyword juliaConstC HydrologyIsLaminarEnum
     262syn keyword juliaConstC HydrologyIsTransitionEnum
    263263syn keyword juliaConstC HydrologyIsWaterPressureArmaEnum
    264264syn keyword juliaConstC HydrologyMeltFlagEnum
     
    266266syn keyword juliaConstC HydrologyNumBasinsEnum
    267267syn keyword juliaConstC HydrologyNumRequestedOutputsEnum
     268syn keyword juliaConstC HydrologyOmegaEnum
    268269syn keyword juliaConstC HydrologyPressureMeltCoefficientEnum
    269270syn keyword juliaConstC HydrologyRelaxationEnum
    270271syn keyword juliaConstC HydrologyRequestedOutputsEnum
    271272syn keyword juliaConstC HydrologySedimentKmaxEnum
     273syn keyword juliaConstC HydrologySheetAlphaEnum
     274syn keyword juliaConstC HydrologySheetBetaEnum
    272275syn keyword juliaConstC HydrologyStepsPerStepEnum
    273276syn keyword juliaConstC HydrologyStorageEnum
     
    891894syn keyword juliaConstC HydrologyMoulinInputEnum
    892895syn keyword juliaConstC HydrologyNeumannfluxEnum
    893 syn keyword juliaConstC HydrologyOmegaEnum
    894896syn keyword juliaConstC HydrologyReynoldsEnum
    895 syn keyword juliaConstC HydrologySheetAlphaEnum
    896 syn keyword juliaConstC HydrologySheetBetaEnum
    897897syn keyword juliaConstC HydrologySheetConductivityEnum
    898898syn keyword juliaConstC HydrologySheetThicknessEnum
  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r27925 r27936  
    275275              else if (strcmp(name,"HydrologyEnglacialVoidRatio")==0) return HydrologyEnglacialVoidRatioEnum;
    276276              else if (strcmp(name,"HydrologyIschannels")==0) return HydrologyIschannelsEnum;
    277               else if (strcmp(name,"HydrologyIsLaminar")==0) return HydrologyIsLaminarEnum;
     277              else if (strcmp(name,"HydrologyIsTransition")==0) return HydrologyIsTransitionEnum;
    278278              else if (strcmp(name,"HydrologyIsWaterPressureArma")==0) return HydrologyIsWaterPressureArmaEnum;
    279279              else if (strcmp(name,"HydrologyMeltFlag")==0) return HydrologyMeltFlagEnum;
     
    281281              else if (strcmp(name,"HydrologyNumBasins")==0) return HydrologyNumBasinsEnum;
    282282              else if (strcmp(name,"HydrologyNumRequestedOutputs")==0) return HydrologyNumRequestedOutputsEnum;
     283              else if (strcmp(name,"HydrologyOmega")==0) return HydrologyOmegaEnum;
    283284              else if (strcmp(name,"HydrologyPressureMeltCoefficient")==0) return HydrologyPressureMeltCoefficientEnum;
    284285              else if (strcmp(name,"HydrologyRelaxation")==0) return HydrologyRelaxationEnum;
    285286              else if (strcmp(name,"HydrologyRequestedOutputs")==0) return HydrologyRequestedOutputsEnum;
    286287              else if (strcmp(name,"HydrologySedimentKmax")==0) return HydrologySedimentKmaxEnum;
     288              else if (strcmp(name,"HydrologySheetAlpha")==0) return HydrologySheetAlphaEnum;
     289              else if (strcmp(name,"HydrologySheetBeta")==0) return HydrologySheetBetaEnum;
    287290              else if (strcmp(name,"HydrologyStepsPerStep")==0) return HydrologyStepsPerStepEnum;
    288291              else if (strcmp(name,"HydrologyStorage")==0) return HydrologyStorageEnum;
     
    380383              else if (strcmp(name,"MassFluxSegments")==0) return MassFluxSegmentsEnum;
    381384              else if (strcmp(name,"MassFluxSegmentsPresent")==0) return MassFluxSegmentsPresentEnum;
    382               else if (strcmp(name,"MasstransportHydrostaticAdjustment")==0) return MasstransportHydrostaticAdjustmentEnum;
    383               else if (strcmp(name,"MasstransportIsfreesurface")==0) return MasstransportIsfreesurfaceEnum;
    384               else if (strcmp(name,"MasstransportMinThickness")==0) return MasstransportMinThicknessEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"MasstransportNumRequestedOutputs")==0) return MasstransportNumRequestedOutputsEnum;
     388              if (strcmp(name,"MasstransportHydrostaticAdjustment")==0) return MasstransportHydrostaticAdjustmentEnum;
     389              else if (strcmp(name,"MasstransportIsfreesurface")==0) return MasstransportIsfreesurfaceEnum;
     390              else if (strcmp(name,"MasstransportMinThickness")==0) return MasstransportMinThicknessEnum;
     391              else if (strcmp(name,"MasstransportNumRequestedOutputs")==0) return MasstransportNumRequestedOutputsEnum;
    389392              else if (strcmp(name,"MasstransportPenaltyFactor")==0) return MasstransportPenaltyFactorEnum;
    390393              else if (strcmp(name,"MasstransportRequestedOutputs")==0) return MasstransportRequestedOutputsEnum;
     
    503506              else if (strcmp(name,"SealevelchangeTidalK2")==0) return SealevelchangeTidalK2Enum;
    504507              else if (strcmp(name,"SealevelchangeTidalH2")==0) return SealevelchangeTidalH2Enum;
    505               else if (strcmp(name,"SealevelchangeTidalL2")==0) return SealevelchangeTidalL2Enum;
    506               else if (strcmp(name,"SolidearthSettingsSealevelLoading")==0) return SolidearthSettingsSealevelLoadingEnum;
    507               else if (strcmp(name,"SolidearthSettingsGRD")==0) return SolidearthSettingsGRDEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"SolidearthSettingsRunFrequency")==0) return SolidearthSettingsRunFrequencyEnum;
     511              if (strcmp(name,"SealevelchangeTidalL2")==0) return SealevelchangeTidalL2Enum;
     512              else if (strcmp(name,"SolidearthSettingsSealevelLoading")==0) return SolidearthSettingsSealevelLoadingEnum;
     513              else if (strcmp(name,"SolidearthSettingsGRD")==0) return SolidearthSettingsGRDEnum;
     514              else if (strcmp(name,"SolidearthSettingsRunFrequency")==0) return SolidearthSettingsRunFrequencyEnum;
    512515              else if (strcmp(name,"SolidearthSettingsTimeAcc")==0) return SolidearthSettingsTimeAccEnum;
    513516              else if (strcmp(name,"SolidearthSettingsHoriz")==0) return SolidearthSettingsHorizEnum;
     
    626629              else if (strcmp(name,"SmbSemicTmin")==0) return SmbSemicTminEnum;
    627630              else if (strcmp(name,"SmbSemicTmid")==0) return SmbSemicTmidEnum;
    628               else if (strcmp(name,"SmbSemicTmax")==0) return SmbSemicTmaxEnum;
    629               else if (strcmp(name,"SmbStepsPerStep")==0) return SmbStepsPerStepEnum;
    630               else if (strcmp(name,"SmbSwIdx")==0) return SmbSwIdxEnum;
    631631         else stage=6;
    632632   }
    633633   if(stage==6){
    634               if (strcmp(name,"SmbSWgrad")==0) return SmbSWgradEnum;
     634              if (strcmp(name,"SmbSemicTmax")==0) return SmbSemicTmaxEnum;
     635              else if (strcmp(name,"SmbStepsPerStep")==0) return SmbStepsPerStepEnum;
     636              else if (strcmp(name,"SmbSwIdx")==0) return SmbSwIdxEnum;
     637              else if (strcmp(name,"SmbSWgrad")==0) return SmbSWgradEnum;
    635638              else if (strcmp(name,"SmbT0dry")==0) return SmbT0dryEnum;
    636639              else if (strcmp(name,"SmbT0wet")==0) return SmbT0wetEnum;
     
    749752              else if (strcmp(name,"BasalforcingsGroundediceMeltingRateObs")==0) return BasalforcingsGroundediceMeltingRateObsEnum;
    750753              else if (strcmp(name,"BasalforcingsLinearBasinId")==0) return BasalforcingsLinearBasinIdEnum;
    751               else if (strcmp(name,"BasalforcingsPerturbationMeltingRate")==0) return BasalforcingsPerturbationMeltingRateEnum;
    752               else if (strcmp(name,"BasalforcingsSpatialDeepwaterElevation")==0) return BasalforcingsSpatialDeepwaterElevationEnum;
    753               else if (strcmp(name,"BasalforcingsSpatialDeepwaterMeltingRate")==0) return BasalforcingsSpatialDeepwaterMeltingRateEnum;
    754754         else stage=7;
    755755   }
    756756   if(stage==7){
    757               if (strcmp(name,"BasalforcingsSpatialUpperwaterElevation")==0) return BasalforcingsSpatialUpperwaterElevationEnum;
     757              if (strcmp(name,"BasalforcingsPerturbationMeltingRate")==0) return BasalforcingsPerturbationMeltingRateEnum;
     758              else if (strcmp(name,"BasalforcingsSpatialDeepwaterElevation")==0) return BasalforcingsSpatialDeepwaterElevationEnum;
     759              else if (strcmp(name,"BasalforcingsSpatialDeepwaterMeltingRate")==0) return BasalforcingsSpatialDeepwaterMeltingRateEnum;
     760              else if (strcmp(name,"BasalforcingsSpatialUpperwaterElevation")==0) return BasalforcingsSpatialUpperwaterElevationEnum;
    758761              else if (strcmp(name,"BasalforcingsSpatialUpperwaterMeltingRate")==0) return BasalforcingsSpatialUpperwaterMeltingRateEnum;
    759762              else if (strcmp(name,"BasalforcingsIsmip6BasinId")==0) return BasalforcingsIsmip6BasinIdEnum;
     
    872875              else if (strcmp(name,"FrictionM")==0) return FrictionMEnum;
    873876              else if (strcmp(name,"FrictionP")==0) return FrictionPEnum;
    874               else if (strcmp(name,"FrictionPressureAdjustedTemperature")==0) return FrictionPressureAdjustedTemperatureEnum;
    875               else if (strcmp(name,"FrictionQ")==0) return FrictionQEnum;
    876               else if (strcmp(name,"FrictionSedimentCompressibilityCoefficient")==0) return FrictionSedimentCompressibilityCoefficientEnum;
    877877         else stage=8;
    878878   }
    879879   if(stage==8){
    880               if (strcmp(name,"FrictionTillFrictionAngle")==0) return FrictionTillFrictionAngleEnum;
     880              if (strcmp(name,"FrictionPressureAdjustedTemperature")==0) return FrictionPressureAdjustedTemperatureEnum;
     881              else if (strcmp(name,"FrictionQ")==0) return FrictionQEnum;
     882              else if (strcmp(name,"FrictionSedimentCompressibilityCoefficient")==0) return FrictionSedimentCompressibilityCoefficientEnum;
     883              else if (strcmp(name,"FrictionTillFrictionAngle")==0) return FrictionTillFrictionAngleEnum;
    881884              else if (strcmp(name,"FrictionWaterLayer")==0) return FrictionWaterLayerEnum;
    882885              else if (strcmp(name,"FrictionWaterPressure")==0) return FrictionWaterPressureEnum;
     
    921924              else if (strcmp(name,"HydrologyMoulinInput")==0) return HydrologyMoulinInputEnum;
    922925              else if (strcmp(name,"HydrologyNeumannflux")==0) return HydrologyNeumannfluxEnum;
    923               else if (strcmp(name,"HydrologyOmega")==0) return HydrologyOmegaEnum;
    924926              else if (strcmp(name,"HydrologyReynolds")==0) return HydrologyReynoldsEnum;
    925               else if (strcmp(name,"HydrologySheetAlpha")==0) return HydrologySheetAlphaEnum;
    926               else if (strcmp(name,"HydrologySheetBeta")==0) return HydrologySheetBetaEnum;
    927927              else if (strcmp(name,"HydrologySheetConductivity")==0) return HydrologySheetConductivityEnum;
    928928              else if (strcmp(name,"HydrologySheetThickness")==0) return HydrologySheetThicknessEnum;
  • issm/trunk-jpl/src/m/classes/hydrologyglads.m

    r27926 r27936  
    1111                cavity_spacing            = 0.;
    1212                bump_height               = NaN;
    13                 omega                     = 0; % TH
    14                 sheet_alpha               = 5.0/4.0; % TH
    15                 sheet_beta                = 3.0/2.0; % TH
     13                omega                     = 0;
     14                sheet_alpha               = NaN;
     15                sheet_beta                = NaN;
    1616
    1717                %Channels
     
    1919                channel_conductivity = NaN;
    2020                channel_sheet_width  = 0.;
    21                 channel_alpha        = 5.0/4.0; % TH
    22                 channel_beta         = 3.0/2.0; % TH
     21                channel_alpha        = NaN;
     22                channel_beta         = NaN;
    2323
    2424                %Other
     
    2929                requested_outputs    = {};
    3030                melt_flag            = 0;
    31                 islaminar            = 0;
     31                istransition         = 0;
    3232        end
    3333        methods
     
    5151                        self.pressure_melt_coefficient = 7.5e-8; %K/Pa (See table 1 in Erder et al. 2013)
    5252                        self.cavity_spacing = 2.; %m
    53                         self.omega = 1./2000.; % TH
     53                        self.sheet_alpha = 5.0/4.0;
     54                        self.sheet_beta = 3.0/2.0;
     55                        self.omega = 1./2000.;
    5456
    5557                        %Channel parameters
     
    5759                        self.channel_conductivity = 5.e-2; %Dow's default, Table uses 0.1
    5860                        self.channel_sheet_width = 2.; %m
     61                        self.channel_alpha = 5.0/4.0;
     62                        self.channel_beta = 3.0/2.0;
    5963
    60                         %Otherself.omega = 1./2000.; % TH
     64                        %Otherself.omega = 1./2000.;
    6165                        self.englacial_void_ratio = 1.e-5;% Dow's default, Table from Werder et al. uses 1e-3;
    6266                        self.requested_outputs={'default'};
    6367                        self.melt_flag=0;
    64                         self.islaminar = 0; %by default use GlaDS default turbulent code
     68                        self.istransition = 0; %by default use GlaDS default turbulent code
    6569                end % }}}
    6670                function md = checkconsistency(self,md,solution,analyses) % {{{
     
    7680                        md = checkfield(md,'fieldname','hydrology.cavity_spacing','numel',[1],'>',0);
    7781                        md = checkfield(md,'fieldname','hydrology.bump_height','size',[md.mesh.numberofvertices 1],'>=',0,'NaN',1,'Inf',1);
    78                         md = checkfield(md,'fieldname','hydrology.omega', 'numel', [1], '>=', 0); % TH
    79                         md = checkfield(md,'fieldname','hydrology.sheet_alpha', 'numel', [1], '>', 0); % TH
    80                         md = checkfield(md,'fieldname','hydrology.sheet_beta', 'numel', [1], '>', 0); % TH
     82                        md = checkfield(md,'fieldname','hydrology.omega', 'numel', [1], '>=', 0);
     83                        md = checkfield(md,'fieldname','hydrology.sheet_alpha', 'numel', [1], '>', 0);
     84                        md = checkfield(md,'fieldname','hydrology.sheet_beta', 'numel', [1], '>', 0);
    8185
    8286                        %Channels
     
    8488                        md = checkfield(md,'fieldname','hydrology.channel_conductivity','size',[md.mesh.numberofvertices 1],'>=',0,'NaN',1,'Inf',1);
    8589                        md = checkfield(md,'fieldname','hydrology.channel_sheet_width','numel',[1],'>=',0);
    86                         md = checkfield(md,'fieldname','hydrology.channel_alpha', 'numel', [1], '>', 0); % TH
    87                         md = checkfield(md,'fieldname','hydrology.channel_beta', 'numel', [1], '>', 0); % TH
     90                        md = checkfield(md,'fieldname','hydrology.channel_alpha', 'numel', [1], '>', 0);
     91                        md = checkfield(md,'fieldname','hydrology.channel_beta', 'numel', [1], '>', 0);
    8892
    8993                        %Other
     
    9498                        md = checkfield(md,'fieldname','hydrology.requested_outputs','stringrow',1);
    9599                        md = checkfield(md,'fieldname','hydrology.melt_flag','numel',[1],'values',[0 1 2]);
    96                         md = checkfield(md,'fieldname','hydrology.islaminar','numel',[1],'values',[0 1]);
     100                        md = checkfield(md,'fieldname','hydrology.istransition','numel',[1],'values',[0 1]);
    97101                        if self.melt_flag==1 || self.melt_flag==2
    98102                                md = checkfield(md,'fieldname','basalforcings.groundedice_melting_rate','NaN',1,'Inf',1,'timeseries',1);
     
    104108                        fielddisplay(self,'pressure_melt_coefficient','Pressure melt coefficient (c_t) [K Pa^-1]');
    105109                        fielddisplay(self,'sheet_conductivity','sheet conductivity (k) [m^(7/4) kg^(-1/2)]');
    106                         fielddisplay(self,'sheet_alpha','First sheet-flow exponent (alpha_s) []'); % TH
    107                         fielddisplay(self,'sheet_beta','Second sheet-flow exponent (beta_s) []'); % TH
     110                        fielddisplay(self,'sheet_alpha','First sheet-flow exponent (alpha_s) []');
     111                        fielddisplay(self,'sheet_beta','Second sheet-flow exponent (beta_s) []');
    108112                        fielddisplay(self,'cavity_spacing','cavity spacing (l_r) [m]');
    109113                        fielddisplay(self,'bump_height','typical bump height (h_r) [m]');
    110                         fielddisplay(self,'omega','transition parameter (omega) []'); % TH
     114                        fielddisplay(self,'omega','transition parameter (omega) []');
    111115                        disp(sprintf('      CHANNELS'));
    112116                        fielddisplay(self,'ischannels','Do we allow for channels? 1: yes, 0: no');
    113117                        fielddisplay(self,'channel_conductivity','channel conductivity (k_c) [m^(3/2) kg^(-1/2)]');
    114                         fielddisplay(self,'channel_alpha','First channel-flow exponent (alpha_s) []'); % TH
    115                         fielddisplay(self,'channel_beta','Second channel-flow exponent (beta_s) []'); % TH
     118                        fielddisplay(self,'channel_alpha','First channel-flow exponent (alpha_s) []');
     119                        fielddisplay(self,'channel_beta','Second channel-flow exponent (beta_s) []');
    116120                        fielddisplay(self,'channel_sheet_width','channel sheet width [m]');
    117121                        disp(sprintf('      OTHER'));
     
    122126                        fielddisplay(self,'requested_outputs','additional outputs requested');
    123127                        fielddisplay(self,'melt_flag','User specified basal melt? 0: no (default), 1: use md.basalforcings.groundedice_melting_rate');
    124                         fielddisplay(self,'islaminar','do we use laminar [1] or turbulent physics [0, default]');
     128                        fielddisplay(self,'istransition','do we use transition [1] or turbulent physics [0, default]');
    125129                end % }}}
    126130                function marshall(self,prefix,md,fid) % {{{
     
    131135                        WriteData(fid,prefix,'name','md.hydrology.model','data',5,'format','Integer');
    132136
    133                         %Sheet
     137                                                %Sheet
    134138                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','pressure_melt_coefficient','format','Double');
    135139                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_conductivity','format','DoubleMat','mattype',1);
    136140                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','cavity_spacing','format','Double');
    137141                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','bump_height','format','DoubleMat','mattype',1);
    138                         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','omega','format','Double'); % TH
    139                         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_alpha','format','Double'); % TH
    140                         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_beta','format','Double'); % TH
     142                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','omega','format','Double');
     143                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_alpha','format','Double');
     144                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_beta','format','Double');
    141145
    142146                        %Channels
     
    144148                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_conductivity','format','DoubleMat','mattype',1);
    145149                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_sheet_width','format','Double');
    146                         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_alpha','format','Double'); % TH
    147                         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_beta','format','Double'); % TH
     150                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_alpha','format','Double');
     151                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_beta','format','Double');
    148152
    149153                        %Others
     
    153157                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','englacial_void_ratio','format','Double');
    154158                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','melt_flag','format','Integer');
    155                         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','islaminar','format','Boolean');
     159                        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','istransition','format','Boolean');
    156160                        outputs = self.requested_outputs;
    157161                        pos  = find(ismember(outputs,'default'));
  • issm/trunk-jpl/src/m/classes/hydrologyglads.py

    r27926 r27936  
    1414
    1515    def __init__(self, *args):  # {{{
    16         # Sheet
     16       # Sheet
    1717        self.pressure_melt_coefficient = 0.
    1818        self.sheet_conductivity = np.nan
    1919        self.cavity_spacing = 0.
    2020        self.bump_height = np.nan
    21         self.omega = 0.; #TH
    22         self.sheet_alpha = 5.0/4.0; # TH
    23         self.sheet_beta = 3.0/2.0; # TH
     21        self.omega = 0.;
     22        self.sheet_alpha = np.nan;
     23        self.sheet_beta = np.nan;
    2424
    2525        # Channels
     
    2727        self.channel_conductivity = np.nan
    2828        self.channel_sheet_width = 0.
    29         self.channel_alpha = 5.0/4.0; # TH
    30         self.channel_beta = 3.0/2.0; # TH
     29        self.channel_alpha = np.nan;
     30        self.channel_beta = np.nan;
    3131
    3232        # Other
     
    3737        self.requested_outputs = []
    3838        self.melt_flag = 0
    39         self.islaminar = 0
     39        self.istransition = 0
    4040
    4141        nargs = len(args)
     
    7272        s += '{}\n'.format(fielddisplay(self, 'requested_outputs', 'additional outputs requested'))
    7373        s += '{}\n'.format(fielddisplay(self, 'melt_flag', 'User specified basal melt? 0: no (default), 1: use md.basalforcings.groundedice_melting_rate'))
    74         s += '{}\n'.format(fielddisplay(self, 'islaminar','do we use laminar [1] or turbulent physics [0, default]'))
     74        s += '{}\n'.format(fielddisplay(self, 'istransition','do we use transition [1] or turbulent physics [0, default]'))
    7575        return string
    7676    # }}}
     
    9696        self.pressure_melt_coefficient = 7.5e-8  #K / Pa (See table 1 in Erder et al. 2013)
    9797        self.cavity_spacing = 2.  #m
    98         self.omega = 1./2000.; #TH
     98        self.sheet_alpha = 5.0/4.0;
     99        self.sheet_beta = 3.0/2.0;
     100        self.omega = 1./2000.;
    99101
    100102        # Channel parameters
     
    102104        self.channel_conductivity = 5.e-2  #Dow's default, Table uses 0.1
    103105        self.channel_sheet_width = 2.  #m
    104         self.islaminar = 0  #by default use turbulent physics
     106        self.channel_alpha = 5.0/4.0;
     107        self.channel_beta = 3.0/2.0;
    105108
    106109        # Other
     
    108111        self.requested_outputs = ['default']
    109112        self.melt_flag = 0
    110         self.islaminar = 0  #by default use turbulent physics
     113        self.istransition = 0  #by default use turbulent physics
    111114
    112115        return self
     
    123126        md = checkfield(md, 'fieldname', 'hydrology.cavity_spacing', 'numel', [1], '>', 0)
    124127        md = checkfield(md, 'fieldname', 'hydrology.bump_height', 'size', [md.mesh.numberofvertices], '>=', 0, 'np.nan', 1, 'Inf', 1)
    125         md = checkfield(md,'fieldname','hydrology.omega', 'numel', [1], '>=', 0); # TH
    126         md = checkfield(md,'fieldname','hydrology.sheet_alpha', 'numel', [1], '>', 0); # TH
    127         md = checkfield(md,'fieldname','hydrology.sheet_beta', 'numel', [1], '>', 0); # TH
     128        md = checkfield(md,'fieldname','hydrology.omega', 'numel', [1], '>=', 0);
     129        md = checkfield(md,'fieldname','hydrology.sheet_alpha', 'numel', [1], '>', 0);
     130        md = checkfield(md,'fieldname','hydrology.sheet_beta', 'numel', [1], '>', 0);
    128131
    129132        # Channels
     
    131134        md = checkfield(md, 'fieldname', 'hydrology.channel_conductivity', 'size', [md.mesh.numberofvertices], '>', 0)
    132135        md = checkfield(md, 'fieldname', 'hydrology.channel_sheet_width', 'numel', [1], '>=', 0)
    133         md = checkfield(md,'fieldname','hydrology.channel_alpha', 'numel', [1], '>', 0); # TH
    134         md = checkfield(md,'fieldname','hydrology.channel_beta', 'numel', [1], '>', 0); # TH
     136        md = checkfield(md,'fieldname','hydrology.channel_alpha', 'numel', [1], '>', 0);
     137        md = checkfield(md,'fieldname','hydrology.channel_beta', 'numel', [1], '>', 0);
    135138
    136139        # Other
     
    141144        md = checkfield(md, 'fieldname', 'hydrology.requested_outputs', 'stringrow', 1)
    142145        md = checkfield(md, 'fieldname', 'hydrology.melt_flag', 'numel', [1], 'values', [0, 1])
    143         md = checkfield(md, 'fieldname', 'hydrology.islaminar', 'numel', [1], 'values', [0, 1])
     146        md = checkfield(md, 'fieldname', 'hydrology.istransition', 'numel', [1], 'values', [0, 1])
    144147        if self.melt_flag == 1 or self.melt_flag == 2:
    145148            md = checkfield(md, 'fieldname', 'basalforcings.groundedice_melting_rate', 'NaN', 1, 'Inf', 1, 'timeseries', 1)
     
    151154        WriteData(fid, prefix, 'name', 'md.hydrology.model', 'data', 5, 'format', 'Integer')
    152155
    153         # Sheet
     156       # Sheet
    154157        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'pressure_melt_coefficient', 'format', 'Double')
    155158        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'sheet_conductivity', 'format', 'DoubleMat', 'mattype', 1)
    156159        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'cavity_spacing', 'format', 'Double')
    157160        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'bump_height', 'format', 'DoubleMat', 'mattype', 1)
    158         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','omega','format','Double'); # TH
    159         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_alpha','format','Double'); # TH
    160         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_beta','format','Double'); # TH
     161        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','omega','format','Double');
     162        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_alpha','format','Double');
     163        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','sheet_beta','format','Double');
    161164
    162165        # Channels
     
    164167        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'channel_conductivity', 'format', 'DoubleMat', 'mattype', 1)
    165168        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'channel_sheet_width', 'format', 'Double')
    166         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_alpha','format','Double'); # TH
    167         WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_beta','format','Double'); # TH
     169        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_alpha','format','Double');
     170        WriteData(fid,prefix,'object',self,'class','hydrology','fieldname','channel_beta','format','Double');
    168171
    169172        # Others
     
    173176        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'englacial_void_ratio', 'format', 'Double')
    174177        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'melt_flag', 'format', 'Integer')
    175         WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'islaminar', 'format', 'Boolean')
     178        WriteData(fid, prefix, 'object', self, 'class', 'hydrology', 'fieldname', 'istransition', 'format', 'Boolean')
    176179
    177180        outputs = self.requested_outputs
Note: See TracChangeset for help on using the changeset viewer.