Changeset 19373

Timestamp:

05/21/15 11:43:14 (10 years ago)

Author:

seroussi

Message:

CHG: updated friction law water layer for Derrick

Location:

issm/trunk-jpl

Files:

• src/c/analyses/StressbalanceAnalysis.cpp (modified) (1 diff)
• src/c/classes/Loads/Friction.cpp (modified) (3 diffs)
• src/c/shared/Enum/EnumDefinitions.h (modified) (1 diff)
• src/c/shared/Enum/EnumToStringx.cpp (modified) (2 diffs)
• src/c/shared/Enum/StringToEnumx.cpp (modified) (9 diffs)
• src/m/classes/frictionwaterlayer.m (modified) (4 diffs)
• src/m/enum/EnumDefinitions.py (modified) (1 diff)
• src/m/enum/FrictionFEnum.m (added)
• test/Archives/Archive438.nc (modified) ( previous)
• test/Archives/Archive439.nc (modified) ( previous)
• test/NightlyRun/test438.m (modified) (1 diff)
• test/NightlyRun/test439.m (modified) (1 diff)

issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

@@ -891 +891 @@
         if(frictionlaw==4 || frictionlaw==6) parameters->AddObject(iomodel->CopyConstantObject(FrictionGammaEnum));
         if(frictionlaw==3) parameters->AddObject(iomodel->CopyConstantObject(FrictionCouplingEnum));
+        if(frictionlaw==5) parameters->AddObject(iomodel->CopyConstantObject(FrictionFEnum));
 
 }/*}}}*/

issm/trunk-jpl/src/c/classes/Loads/Friction.cpp

@@ -390 +390 @@
         IssmDouble  r,s;
         IssmDouble  drag_p, drag_q;
-        IssmDouble  Neff;
+        IssmDouble  Neff,F;
         IssmDouble  thickness,bed;
         IssmDouble  vx,vy,vz,vmag;
@@ -397 +397 @@
 
         /*Recover parameters: */
+        element->parameters->FindParam(&F,FrictionFEnum);
         element->GetInputValue(&drag_p,FrictionPEnum);
         element->GetInputValue(&drag_q,FrictionQEnum);
@@ -413 +414 @@
         //From bed and thickness, compute effective pressure when drag is viscous:
         if(bed>0) bed=0;
-        Neff=gravity*(rho_ice*thickness+rho_water*(bed-water_layer));
+        if(water_layer==0) Neff=gravity*rho_ice*thickness+gravity*rho_water*bed;
+        else if(water_layer>0) Neff=gravity*rho_ice*thickness*F;
+        else _error_("negative water layer thickness");
         if(Neff<0) Neff=0;
 

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

@@ -99 +99 @@
         FrictionLawEnum,
         FrictionGammaEnum,
+        FrictionFEnum,
         FrictionWaterLayerEnum,
         FrictionEffectivePressureEnum,

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

@@ -105 +105 @@
                 case FrictionLawEnum : return "FrictionLaw";
                 case FrictionGammaEnum : return "FrictionGamma";
+                case FrictionFEnum : return "FrictionF";
                 case FrictionWaterLayerEnum : return "FrictionWaterLayer";
                 case FrictionEffectivePressureEnum : return "FrictionEffectivePressure";
@@ -576 +577 @@
                 case TemperatureEnum : return "Temperature";
                 case TemperaturePicardEnum : return "TemperaturePicard";
+                case TemperaturePDDEnum : return "TemperaturePDD";
                 case ThicknessAbsMisfitEnum : return "ThicknessAbsMisfit";
                 case SurfaceAbsMisfitEnum : return "SurfaceAbsMisfit";

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

@@ -105 +105 @@
               else if (strcmp(name,"FrictionLaw")==0) return FrictionLawEnum;
               else if (strcmp(name,"FrictionGamma")==0) return FrictionGammaEnum;
+              else if (strcmp(name,"FrictionF")==0) return FrictionFEnum;
               else if (strcmp(name,"FrictionWaterLayer")==0) return FrictionWaterLayerEnum;
               else if (strcmp(name,"FrictionEffectivePressure")==0) return FrictionEffectivePressureEnum;
@@ -136 +137 @@
               else if (strcmp(name,"HydrologydcEplCompressibility")==0) return HydrologydcEplCompressibilityEnum;
               else if (strcmp(name,"HydrologydcEplPorosity")==0) return HydrologydcEplPorosityEnum;
-              else if (strcmp(name,"HydrologydcEplInitialThickness")==0) return HydrologydcEplInitialThicknessEnum;
          else stage=2;
    }
    if(stage==2){
-              if (strcmp(name,"HydrologydcEplColapseThickness")==0) return HydrologydcEplColapseThicknessEnum;
+              if (strcmp(name,"HydrologydcEplInitialThickness")==0) return HydrologydcEplInitialThicknessEnum;
+              else if (strcmp(name,"HydrologydcEplColapseThickness")==0) return HydrologydcEplColapseThicknessEnum;
               else if (strcmp(name,"HydrologydcEplMaxThickness")==0) return HydrologydcEplMaxThicknessEnum;
               else if (strcmp(name,"HydrologydcEplThickness")==0) return HydrologydcEplThicknessEnum;
@@ -259 +260 @@
               else if (strcmp(name,"MeshElements")==0) return MeshElementsEnum;
               else if (strcmp(name,"MeshLowerelements")==0) return MeshLowerelementsEnum;
-              else if (strcmp(name,"MeshNumberofelements2d")==0) return MeshNumberofelements2dEnum;
          else stage=3;
    }
    if(stage==3){
-              if (strcmp(name,"MeshNumberofelements")==0) return MeshNumberofelementsEnum;
+              if (strcmp(name,"MeshNumberofelements2d")==0) return MeshNumberofelements2dEnum;
+              else if (strcmp(name,"MeshNumberofelements")==0) return MeshNumberofelementsEnum;
               else if (strcmp(name,"MeshNumberoflayers")==0) return MeshNumberoflayersEnum;
               else if (strcmp(name,"MeshNumberofvertices2d")==0) return MeshNumberofvertices2dEnum;
@@ -382 +383 @@
               else if (strcmp(name,"SurfaceforcingsHref")==0) return SurfaceforcingsHrefEnum;
               else if (strcmp(name,"SurfaceforcingsSmbref")==0) return SurfaceforcingsSmbrefEnum;
-              else if (strcmp(name,"SurfaceforcingsBPos")==0) return SurfaceforcingsBPosEnum;
          else stage=4;
    }
    if(stage==4){
-              if (strcmp(name,"SurfaceforcingsBNeg")==0) return SurfaceforcingsBNegEnum;
+              if (strcmp(name,"SurfaceforcingsBPos")==0) return SurfaceforcingsBPosEnum;
+              else if (strcmp(name,"SurfaceforcingsBNeg")==0) return SurfaceforcingsBNegEnum;
               else if (strcmp(name,"SMBhenning")==0) return SMBhenningEnum;
               else if (strcmp(name,"SMBcomponents")==0) return SMBcomponentsEnum;
@@ -505 +506 @@
               else if (strcmp(name,"Matice")==0) return MaticeEnum;
               else if (strcmp(name,"Matdamageice")==0) return MatdamageiceEnum;
-              else if (strcmp(name,"Matpar")==0) return MatparEnum;
          else stage=5;
    }
    if(stage==5){
-              if (strcmp(name,"Node")==0) return NodeEnum;
+              if (strcmp(name,"Matpar")==0) return MatparEnum;
+              else if (strcmp(name,"Node")==0) return NodeEnum;
               else if (strcmp(name,"Numericalflux")==0) return NumericalfluxEnum;
               else if (strcmp(name,"NumericalfluxType")==0) return NumericalfluxTypeEnum;
@@ -588 +589 @@
               else if (strcmp(name,"Temperature")==0) return TemperatureEnum;
               else if (strcmp(name,"TemperaturePicard")==0) return TemperaturePicardEnum;
+              else if (strcmp(name,"TemperaturePDD")==0) return TemperaturePDDEnum;
               else if (strcmp(name,"ThicknessAbsMisfit")==0) return ThicknessAbsMisfitEnum;
               else if (strcmp(name,"SurfaceAbsMisfit")==0) return SurfaceAbsMisfitEnum;
@@ -627 +629 @@
               else if (strcmp(name,"StressTensor")==0) return StressTensorEnum;
               else if (strcmp(name,"StressTensorxx")==0) return StressTensorxxEnum;
-              else if (strcmp(name,"StressTensorxy")==0) return StressTensorxyEnum;
-              else if (strcmp(name,"StressTensorxz")==0) return StressTensorxzEnum;
          else stage=6;
    }
    if(stage==6){
-              if (strcmp(name,"StressTensoryy")==0) return StressTensoryyEnum;
+              if (strcmp(name,"StressTensorxy")==0) return StressTensorxyEnum;
+              else if (strcmp(name,"StressTensorxz")==0) return StressTensorxzEnum;
+              else if (strcmp(name,"StressTensoryy")==0) return StressTensoryyEnum;
               else if (strcmp(name,"StressTensoryz")==0) return StressTensoryzEnum;
               else if (strcmp(name,"StressTensorzz")==0) return StressTensorzzEnum;
@@ -750 +752 @@
               else if (strcmp(name,"Outputdefinition61")==0) return Outputdefinition61Enum;
               else if (strcmp(name,"Outputdefinition62")==0) return Outputdefinition62Enum;
-              else if (strcmp(name,"Outputdefinition63")==0) return Outputdefinition63Enum;
-              else if (strcmp(name,"Outputdefinition64")==0) return Outputdefinition64Enum;
          else stage=7;
    }
    if(stage==7){
-              if (strcmp(name,"Outputdefinition65")==0) return Outputdefinition65Enum;
+              if (strcmp(name,"Outputdefinition63")==0) return Outputdefinition63Enum;
+              else if (strcmp(name,"Outputdefinition64")==0) return Outputdefinition64Enum;
+              else if (strcmp(name,"Outputdefinition65")==0) return Outputdefinition65Enum;
               else if (strcmp(name,"Outputdefinition66")==0) return Outputdefinition66Enum;
               else if (strcmp(name,"Outputdefinition67")==0) return Outputdefinition67Enum;
@@ -873 +875 @@
               else if (strcmp(name,"QmuInName")==0) return QmuInNameEnum;
               else if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
-              else if (strcmp(name,"Regular")==0) return RegularEnum;
-              else if (strcmp(name,"Scaled")==0) return ScaledEnum;
          else stage=8;
    }
    if(stage==8){
-              if (strcmp(name,"Separate")==0) return SeparateEnum;
+              if (strcmp(name,"Regular")==0) return RegularEnum;
+              else if (strcmp(name,"Scaled")==0) return ScaledEnum;
+              else if (strcmp(name,"Separate")==0) return SeparateEnum;
               else if (strcmp(name,"Sset")==0) return SsetEnum;
               else if (strcmp(name,"Verbose")==0) return VerboseEnum;

issm/trunk-jpl/src/m/classes/frictionwaterlayer.m

@@ -7 +7 @@
         properties (SetAccess=public) 
                 coefficient = NaN;
+                f           = NaN;
                 p           = NaN;
                 q           = NaN;
@@ -31 +32 @@
 
                         md = checkfield(md,'fieldname','friction.coefficient','timeseries',1,'NaN',1);
+                        md = checkfield(md,'fieldname','friction.f','size',[1 1],'NaN',1);
                         md = checkfield(md,'fieldname','friction.q','NaN',1,'size',[md.mesh.numberofelements 1]);
                         md = checkfield(md,'fieldname','friction.p','NaN',1,'size',[md.mesh.numberofelements 1]);
@@ -45 +47 @@
                         disp(sprintf('Basal shear stress parameters: tau_b = coefficient^2 * Neff ^r * |u_b|^(s-1) * u_b * 1/f(T)\n(effective stress Neff=rho_ice*g*thickness+rho_water*g*(bed+water_layer), r=q/p and s=1/p)'));
                         fielddisplay(self,'coefficient','frictiontemp coefficient [SI]');
+                        fielddisplay(self,'f','f variable for effective pressure');
                         fielddisplay(self,'p','p exponent');
                         fielddisplay(self,'q','q exponent');
@@ -53 +56 @@
                         WriteData(fid,'enum',FrictionLawEnum,'data',5,'format','Integer');
                         WriteData(fid,'class','friction','object',self,'fieldname','coefficient','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1);
+                        WriteData(fid,'class','friction','object',self,'fieldname','f','format','Double');
                         WriteData(fid,'class','friction','object',self,'fieldname','p','format','DoubleMat','mattype',2);
                         WriteData(fid,'class','friction','object',self,'fieldname','q','format','DoubleMat','mattype',2);

issm/trunk-jpl/src/m/enum/EnumDefinitions.py

@@ -97 +97 @@
 def FrictionLawEnum(): return StringToEnum("FrictionLaw")[0]
 def FrictionGammaEnum(): return StringToEnum("FrictionGamma")[0]
+def FrictionFEnum(): return StringToEnum("FrictionF")[0]
 def FrictionWaterLayerEnum(): return StringToEnum("FrictionWaterLayer")[0]
 def FrictionEffectivePressureEnum(): return StringToEnum("FrictionEffectivePressure")[0]

issm/trunk-jpl/test/NightlyRun/test438.m

@@ -8 +8 @@
 md.friction.water_layer(:,2)=1;
 md.friction.water_layer(md.mesh.numberofvertices+1,:)=[1 2];
+md.friction.f=0.8;
 md.cluster=generic('name',oshostname(),'np',3);
 md=solve(md,TransientSolutionEnum());

issm/trunk-jpl/test/NightlyRun/test439.m

@@ -9 +9 @@
 md.friction.water_layer(:,2)=1;
 md.friction.water_layer(md.mesh.numberofvertices+1,:)=[1 2];
+md.friction.f=0.8;
 md.cluster=generic('name',oshostname(),'np',3);
 md=solve(md,TransientSolutionEnum());