Changeset 16592

Timestamp:

10/31/13 12:08:07 (11 years ago)

Author:

bdef

Message:

Adding a slope computation for the EPL head

Location:

issm/trunk-jpl/src

Files:

• c/classes/Elements/Tria.cpp (modified) (2 diffs)
• c/cores/hydrology_core.cpp (modified) (4 diffs)
• c/shared/Enum/EnumDefinitions.h (modified) (1 diff)
• c/shared/Enum/EnumToStringx.cpp (modified) (1 diff)
• c/shared/Enum/StringToEnumx.cpp (modified) (5 diffs)
• c/solutionsequences/solutionsequence_hydro_nonlinear.cpp (modified) (3 diffs)
• m/classes/clusters/greenplanet.m (modified) (1 diff)
• m/contrib/paraview/exportVTK.m (modified) (1 diff)
• m/enum/EnumDefinitions.py (modified) (1 diff)
• m/enum/EplHeadSlopeXEnum.m (added)
• m/enum/EplHeadSlopeYEnum.m (added)
• m/plot/plot_contour.m (modified) (1 diff)

issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

@@ -573 +573 @@
                 case BedSlopeXEnum:     input2 = inputs->GetInput(BedEnum);     _assert_(input2); break;
                 case BedSlopeYEnum:     input2 = inputs->GetInput(BedEnum);     _assert_(input2); break;
+                case EplHeadSlopeXEnum: input2 = inputs->GetInput(EplHeadEnum); _assert_(input2); break;
+                case EplHeadSlopeYEnum: input2 = inputs->GetInput(EplHeadEnum); _assert_(input2); break;
                 default: input = inputs->GetInput(input_enum);
         }
@@ -587 +589 @@
                 if(input2) input2->GetInputDerivativeValue(&slopes[0],&xyz_list[0][0],gauss);
                 switch(input_enum){
-                        case SurfaceSlopeXEnum: case BedSlopeXEnum: value = slopes[0]; break;
-                        case SurfaceSlopeYEnum: case BedSlopeYEnum: value = slopes[1]; break;
+                        case SurfaceSlopeXEnum: case BedSlopeXEnum: case EplHeadSlopeXEnum: value = slopes[0]; break;
+                        case SurfaceSlopeYEnum: case BedSlopeYEnum: case EplHeadSlopeYEnum: value = slopes[1]; break;
                         default: input->GetInputValue(&value,gauss);
                 }

issm/trunk-jpl/src/c/cores/hydrology_core.cpp

@@ -57 +57 @@
                 femmodel->parameters->SetParam(step,StepEnum);
 
+                /*Using the Shreve based Model*/
                 if (hydrology_model==HydrologyshreveEnum){
                         if(VerboseSolution()) _printf0_("   computing water column\n");
@@ -77 +78 @@
                 }
 
+                /*Using the double continuum model*/
                 else if (hydrology_model==HydrologydcEnum){
                         InputDuplicatex(femmodel,SedimentHeadEnum,SedimentHeadOldEnum);
@@ -83 +85 @@
                                 InputDuplicatex(femmodel,EplHeadEnum,EplHeadOldEnum);
                         }
-
+                        
+                        /*Proceed now to heads computations*/
                         if(VerboseSolution()) _printf0_("   computing water head\n");
                         solutionsequence_hydro_nonlinear(femmodel);
@@ -89 +92 @@
                                 if(VerboseSolution()) _printf0_("   saving results \n");
                                 if(isefficientlayer){
-                                        int outputs[4] = {SedimentHeadEnum,SedimentHeadResidualEnum,EplHeadEnum,HydrologydcMaskEplactiveEnum};
-                                        femmodel->RequestedOutputsx(&femmodel->results,&outputs[0],4);
+                                        int outputs[6] = {SedimentHeadEnum,SedimentHeadResidualEnum,EplHeadEnum,HydrologydcMaskEplactiveEnum,EplHeadSlopeXEnum,EplHeadSlopeYEnum};
+                                        femmodel->RequestedOutputsx(&femmodel->results,&outputs[0],6);
                                 }
                                 else{

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

@@ -100 +100 @@
         EplHeadEnum,
         EplHeadOldEnum,
-  HydrologydcRelTolEnum,
+  EplHeadSlopeXEnum,
+        EplHeadSlopeYEnum,
+        HydrologydcRelTolEnum,
         HydrologydcSpcsedimentHeadEnum,
         HydrologydcSedimentCompressibilityEnum,

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

@@ -108 +108 @@
                 case EplHeadEnum : return "EplHead";
                 case EplHeadOldEnum : return "EplHeadOld";
+                case EplHeadSlopeXEnum : return "EplHeadSlopeX";
+                case EplHeadSlopeYEnum : return "EplHeadSlopeY";
                 case HydrologydcRelTolEnum : return "HydrologydcRelTol";
                 case HydrologydcSpcsedimentHeadEnum : return "HydrologydcSpcsedimentHead";

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

@@ -108 +108 @@
               else if (strcmp(name,"EplHead")==0) return EplHeadEnum;
               else if (strcmp(name,"EplHeadOld")==0) return EplHeadOldEnum;
+              else if (strcmp(name,"EplHeadSlopeX")==0) return EplHeadSlopeXEnum;
+              else if (strcmp(name,"EplHeadSlopeY")==0) return EplHeadSlopeYEnum;
               else if (strcmp(name,"HydrologydcRelTol")==0) return HydrologydcRelTolEnum;
               else if (strcmp(name,"HydrologydcSpcsedimentHead")==0) return HydrologydcSpcsedimentHeadEnum;
@@ -135 +137 @@
               else if (strcmp(name,"InversionControlParameters")==0) return InversionControlParametersEnum;
               else if (strcmp(name,"InversionCostFunctionThreshold")==0) return InversionCostFunctionThresholdEnum;
-              else if (strcmp(name,"InversionCostFunctionsCoefficients")==0) return InversionCostFunctionsCoefficientsEnum;
-              else if (strcmp(name,"InversionCostFunctions")==0) return InversionCostFunctionsEnum;
          else stage=2;
    }
    if(stage==2){
-              if (strcmp(name,"InversionGradientScaling")==0) return InversionGradientScalingEnum;
+              if (strcmp(name,"InversionCostFunctionsCoefficients")==0) return InversionCostFunctionsCoefficientsEnum;
+              else if (strcmp(name,"InversionCostFunctions")==0) return InversionCostFunctionsEnum;
+              else if (strcmp(name,"InversionGradientScaling")==0) return InversionGradientScalingEnum;
               else if (strcmp(name,"InversionIscontrol")==0) return InversionIscontrolEnum;
               else if (strcmp(name,"InversionTao")==0) return InversionTaoEnum;
@@ -258 +260 @@
               else if (strcmp(name,"Surface")==0) return SurfaceEnum;
               else if (strcmp(name,"SurfaceforcingsPrecipitation")==0) return SurfaceforcingsPrecipitationEnum;
-              else if (strcmp(name,"SurfaceforcingsMassBalance")==0) return SurfaceforcingsMassBalanceEnum;
-              else if (strcmp(name,"SurfaceforcingsIspdd")==0) return SurfaceforcingsIspddEnum;
          else stage=3;
    }
    if(stage==3){
-              if (strcmp(name,"SurfaceforcingsDesfac")==0) return SurfaceforcingsDesfacEnum;
+              if (strcmp(name,"SurfaceforcingsMassBalance")==0) return SurfaceforcingsMassBalanceEnum;
+              else if (strcmp(name,"SurfaceforcingsIspdd")==0) return SurfaceforcingsIspddEnum;
+              else if (strcmp(name,"SurfaceforcingsDesfac")==0) return SurfaceforcingsDesfacEnum;
               else if (strcmp(name,"SurfaceforcingsS0p")==0) return SurfaceforcingsS0pEnum;
               else if (strcmp(name,"SurfaceforcingsIssmbgradients")==0) return SurfaceforcingsIssmbgradientsEnum;
@@ -381 +383 @@
               else if (strcmp(name,"FileParam")==0) return FileParamEnum;
               else if (strcmp(name,"Input")==0) return InputEnum;
-              else if (strcmp(name,"IntInput")==0) return IntInputEnum;
-              else if (strcmp(name,"InputToExtrude")==0) return InputToExtrudeEnum;
          else stage=4;
    }
    if(stage==4){
-              if (strcmp(name,"InputToL2Project")==0) return InputToL2ProjectEnum;
+              if (strcmp(name,"IntInput")==0) return IntInputEnum;
+              else if (strcmp(name,"InputToExtrude")==0) return InputToExtrudeEnum;
+              else if (strcmp(name,"InputToL2Project")==0) return InputToL2ProjectEnum;
               else if (strcmp(name,"IntParam")==0) return IntParamEnum;
               else if (strcmp(name,"IntVecParam")==0) return IntVecParamEnum;
@@ -504 +506 @@
               else if (strcmp(name,"HydrologyWaterVy")==0) return HydrologyWaterVyEnum;
               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 stage=5;
    }
    if(stage==5){
-              if (strcmp(name,"StressTensorxz")==0) return StressTensorxzEnum;
+              if (strcmp(name,"StressTensorxx")==0) return StressTensorxxEnum;
+              else 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;

issm/trunk-jpl/src/c/solutionsequences/solutionsequence_hydro_nonlinear.cpp

@@ -7 +7 @@
 #include "../shared/shared.h"
 #include "../modules/modules.h"
+/*FIXME, dirty hack to get the solutionsequence linear needed to compute the slopes*/
+#include "../solutionsequences/solutionsequences.h"
 
 void solutionsequence_hydro_nonlinear(FemModel* femmodel){
@@ -92 +94 @@
                         InputUpdateFromSolutionx(femmodel,ug_sed);
                         ConstraintsStatex(&constraints_converged,&num_unstable_constraints,femmodel);
-
                         if (!sedconverged){
                                 if(VerboseConvergence()) _printf0_("   # Sediment unstable constraints = " << num_unstable_constraints << "\n");
@@ -113 +114 @@
                 /*Second layer*/
                 if(isefficientlayer){
+
+                        /*Start by retrieving the EPL head slopes*/
+                        if(VerboseSolution()) _printf0_("computing EPL Head slope...\n");
+                        femmodel->SetCurrentConfiguration(L2ProjectionBaseAnalysisEnum);
+                        femmodel->parameters->SetParam(EplHeadSlopeXEnum,InputToL2ProjectEnum);
+                        solutionsequence_linear(femmodel);
+                        femmodel->parameters->SetParam(EplHeadSlopeYEnum,InputToL2ProjectEnum);
+                        solutionsequence_linear(femmodel);
+                        
                         femmodel->SetCurrentConfiguration(HydrologyDCEfficientAnalysisEnum);
                         InputUpdateFromConstantx(femmodel,true,ResetPenaltiesEnum);

issm/trunk-jpl/src/m/classes/clusters/greenplanet.m

@@ -147 +147 @@
 
                          disp('launching solution sequence on remote cluster');
-                         launchcommand=['cd ' cluster.executionpath ' && rm -rf ./' dirname ' && mkdir ' dirname ...
-                                 ' && cd ' dirname ' && mv ../' dirname '.tar.gz ./ && tar -zxf ' dirname '.tar.gz  && qsub ' modelname '.queue '];
+                         % launchcommand=['cd ' cluster.executionpath ' && rm -rf ./' dirname ' && mkdir ' dirname ...
+                         %       ' && cd ' dirname ' && mv ../' dirname '.tar.gz ./ && tar -zxf ' dirname '.tar.gz  && hostname && qsub ' modelname '.queue '];
                          issmssh(cluster.name,cluster.login,cluster.port,launchcommand);
                  end %}}}

issm/trunk-jpl/src/m/contrib/paraview/exportVTK.m

@@ -118 +118 @@
                 fieldnames=fields(res_struct);
                 num_of_fields=length(fieldnames);
-                for k=1:num_of_sols
+                for k=1:num_of_fields
                         if ((length(res_struct.(fieldnames{k})))==num_of_points);
                                 %paraview does not like NaN, replacing

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

@@ -100 +100 @@
 def EplHeadEnum(): return StringToEnum("EplHead")[0]
 def EplHeadOldEnum(): return StringToEnum("EplHeadOld")[0]
+def EplHeadSlopeXEnum(): return StringToEnum("EplHeadSlopeX")[0]
+def EplHeadSlopeYEnum(): return StringToEnum("EplHeadSlopeY")[0]
 def HydrologydcRelTolEnum(): return StringToEnum("HydrologydcRelTol")[0]
 def HydrologydcSpcsedimentHeadEnum(): return StringToEnum("HydrologydcSpcsedimentHead")[0]

issm/trunk-jpl/src/m/plot/plot_contour.m

@@ -230 +230 @@
                         h=patch('Xdata',[xc;NaN],'Ydata',[yc;NaN],'facecolor','none','edgecolor',color,'linewidth',linewidth);
                         c = horzcat([level, xc'; length(xc), yc']);
-                        clabel(c,h,'FontSize',10,'labelspacing',20000,'color',color);
+                        %clabel(c,h,'FontSize',10,'labelspacing',20000,'color',color);
                         hold on
                 end