Changeset 5707

Timestamp:

09/08/10 14:42:36 (15 years ago)

Author:

Mathieu Morlighem

Message:

Final fixing of SystemMatrices for Thermal

Location:

issm/trunk/src

Files:

• c/EnumDefinitions/EnumDefinitions.h (modified) (1 diff)
• c/EnumDefinitions/EnumToString.cpp (modified) (1 diff)
• c/EnumDefinitions/StringToEnum.cpp (modified) (1 diff)
• c/objects/Elements/Penta.cpp (modified) (3 diffs)
• c/objects/Inputs/DoubleInput.cpp (modified) (1 diff)
• c/objects/Loads/Pengrid.cpp (modified) (3 diffs)
• c/solvers/solver_thermal_nonlinear.cpp (modified) (2 diffs)
• c/toolkits/petsc/patches/PetscDetermineLocalSize.cpp (modified) (1 diff)
• m/enum/AbsoluteEnum.m (modified) (1 diff)
• m/enum/AdjointEnum.m (modified) (1 diff)
• m/enum/BetaEnum.m (modified) (1 diff)
• m/enum/BoolExternalResultEnum.m (modified) (1 diff)
• m/enum/CmGradientEnum.m (modified) (1 diff)
• m/enum/CmJumpEnum.m (modified) (1 diff)
• m/enum/CmMaxEnum.m (modified) (1 diff)
• m/enum/CmMinEnum.m (modified) (1 diff)
• m/enum/ColinearEnum.m (modified) (1 diff)
• m/enum/ConnectivityEnum.m (modified) (1 diff)
• m/enum/ControlParameterEnum.m (modified) (1 diff)
• m/enum/ControlSteadyEnum.m (modified) (1 diff)
• m/enum/DakotaParameterEnum.m (modified) (1 diff)
• m/enum/DimEnum.m (modified) (1 diff)
• m/enum/DoubleElementResultEnum.m (modified) (1 diff)
• m/enum/DoubleExternalResultEnum.m (modified) (1 diff)
• m/enum/DoubleMatExternalResultEnum.m (modified) (1 diff)
• m/enum/DoubleVecExternalResultEnum.m (modified) (1 diff)
• m/enum/EnumToString.m (modified) (1 diff)
• m/enum/EpsAbsEnum.m (modified) (1 diff)
• m/enum/EpsCmEnum.m (modified) (1 diff)
• m/enum/EpsRelEnum.m (modified) (1 diff)
• m/enum/EpsResEnum.m (modified) (1 diff)
• m/enum/GradientEnum.m (modified) (1 diff)
• m/enum/GravityEnum.m (modified) (1 diff)
• m/enum/HeatCapacityEnum.m (modified) (1 diff)
• m/enum/IndexEnum.m (modified) (1 diff)
• m/enum/IndexedEnum.m (modified) (1 diff)
• m/enum/IntExternalResultEnum.m (modified) (1 diff)
• m/enum/IntersectEnum.m (modified) (1 diff)
• m/enum/IsHutterEnum.m (modified) (1 diff)
• m/enum/IsMacAyealPattynEnum.m (modified) (1 diff)
• m/enum/IsStokesEnum.m (modified) (1 diff)
• m/enum/JEnum.m (modified) (1 diff)
• m/enum/LatentHeatEnum.m (modified) (1 diff)
• m/enum/LowmemEnum.m (modified) (1 diff)
• m/enum/MaxAbsVxEnum.m (modified) (1 diff)
• m/enum/MaxAbsVyEnum.m (modified) (1 diff)
• m/enum/MaxAbsVzEnum.m (modified) (1 diff)
• m/enum/MaxIterEnum.m (modified) (1 diff)
• m/enum/MaxNonlinearIterationsEnum.m (modified) (1 diff)
• m/enum/MaxVelEnum.m (modified) (1 diff)
• m/enum/MaxVxEnum.m (modified) (1 diff)
• m/enum/MaxVyEnum.m (modified) (1 diff)
• m/enum/MaxVzEnum.m (modified) (1 diff)
• m/enum/MeltingPointEnum.m (modified) (1 diff)
• m/enum/MinMechanicalConstraintsEnum.m (modified) (1 diff)
• m/enum/MinThermalConstraintsEnum.m (modified) (1 diff)
• m/enum/MinVelEnum.m (modified) (1 diff)
• m/enum/MinVxEnum.m (modified) (1 diff)
• m/enum/MinVyEnum.m (modified) (1 diff)
• m/enum/MinVzEnum.m (modified) (1 diff)
• m/enum/MiniEnum.m (modified) (1 diff)
• m/enum/MixedLayerCapacityEnum.m (modified) (1 diff)
• m/enum/NStepsEnum.m (modified) (1 diff)
• m/enum/NdtEnum.m (modified) (1 diff)
• m/enum/NodalEnum.m (modified) (1 diff)
• m/enum/NumOutputEnum.m (modified) (1 diff)
• m/enum/NumRiftsEnum.m (modified) (1 diff)
• m/enum/NumberOfElementsEnum.m (modified) (1 diff)
• m/enum/NumberOfNodesEnum.m (modified) (1 diff)
• m/enum/NumberOfVerticesEnum.m (modified) (1 diff)
• m/enum/OldGradientEnum.m (modified) (1 diff)
• m/enum/OptScalEnum.m (modified) (1 diff)
• m/enum/OutputFilePointerEnum.m (modified) (1 diff)
• m/enum/OutputFrequencyEnum.m (modified) (1 diff)
• m/enum/P0Enum.m (modified) (1 diff)
• m/enum/P1DGEnum.m (modified) (1 diff)
• m/enum/P1Enum.m (modified) (1 diff)
• m/enum/ParameterOutputEnum.m (modified) (1 diff)
• m/enum/PatchEnum.m (modified) (1 diff)
• m/enum/PatchNodesEnum.m (modified) (1 diff)
• m/enum/PatchVerticesEnum.m (modified) (1 diff)
• m/enum/PenaltyMeltingEnum.m (modified) (1 diff)
• m/enum/PentaVertexElementResultEnum.m (modified) (1 diff)
• m/enum/PetscVecExternalResultEnum.m (modified) (1 diff)
• m/enum/QmuAnalysisEnum.m (modified) (1 diff)
• m/enum/QmuErrNameEnum.m (modified) (1 diff)
• m/enum/QmuInNameEnum.m (modified) (1 diff)
• m/enum/QmuMassFluxSegmentsEnum.m (modified) (1 diff)
• m/enum/QmuNPartEnum.m (modified) (1 diff)
• m/enum/QmuNumberOfResponsesEnum.m (modified) (1 diff)
• m/enum/QmuOutNameEnum.m (modified) (1 diff)
• m/enum/QmuPartEnum.m (modified) (1 diff)
• m/enum/QmuVxEnum.m (modified) (1 diff)
• m/enum/QmuVyEnum.m (modified) (1 diff)
• m/enum/QmuVzEnum.m (modified) (1 diff)
• m/enum/RegularEnum.m (modified) (1 diff)
• m/enum/RelativeEnum.m (modified) (1 diff)
• m/enum/ResidualEnum.m (modified) (1 diff)
• m/enum/ResponseDescriptorsEnum.m (modified) (1 diff)
• m/enum/RhoIceEnum.m (modified) (1 diff)
• m/enum/RhoWaterEnum.m (modified) (1 diff)
• m/enum/ScaledEnum.m (modified) (1 diff)
• m/enum/SeparateEnum.m (modified) (1 diff)
• m/enum/SolverStringEnum.m (modified) (1 diff)
• m/enum/SparsityEnum.m (modified) (1 diff)
• m/enum/StringExternalResultEnum.m (modified) (1 diff)
• m/enum/StringToEnum.m (modified) (1 diff)
• m/enum/TemperatureOldEnum.m (modified) (1 diff)
• m/enum/ThermalConductivityEnum.m (modified) (1 diff)
• m/enum/ThermalExchangeVelocityEnum.m (modified) (1 diff)
• m/enum/ThicknessEnum.m (modified) (1 diff)
• m/enum/ThicknessObsEnum.m (modified) (1 diff)
• m/enum/TimeEnum.m (modified) (1 diff)
• m/enum/TolXEnum.m (modified) (1 diff)
• m/enum/TriaVertexElementResultEnum.m (modified) (1 diff)
• m/enum/TypeEnum.m (modified) (1 diff)
• m/enum/VariableDescriptorsEnum.m (modified) (1 diff)
• m/enum/VelEnum.m (modified) (1 diff)
• m/enum/VelObsEnum.m (modified) (1 diff)
• m/enum/VerboseEnum.m (modified) (1 diff)
• m/enum/ViscosityOvershootEnum.m (modified) (1 diff)
• m/enum/VxAverageEnum.m (modified) (1 diff)
• m/enum/VxEnum.m (modified) (1 diff)
• m/enum/VxObsEnum.m (modified) (1 diff)
• m/enum/VxOldEnum.m (modified) (1 diff)
• m/enum/VyAverageEnum.m (modified) (1 diff)
• m/enum/VyEnum.m (modified) (1 diff)
• m/enum/VyObsEnum.m (modified) (1 diff)
• m/enum/VyOldEnum.m (modified) (1 diff)
• m/enum/VzAverageEnum.m (modified) (1 diff)
• m/enum/VzEnum.m (modified) (1 diff)
• m/enum/VzObsEnum.m (modified) (1 diff)
• m/enum/VzOldEnum.m (modified) (1 diff)
• m/enum/WaitOnLockEnum.m (modified) (1 diff)
• m/enum/WeightsEnum.m (modified) (1 diff)
• m/enum/YtsEnum.m (modified) (1 diff)
• m/solvers/solver_thermal_nonlinear.m (modified) (2 diffs)

issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

@@ -229 +229 @@
         SurfaceSlopeYEnum,
         TemperatureEnum,
+        TemporaryTemperatureEnum,
         TemperatureOldEnum,
         ThicknessEnum,

issm/trunk/src/c/EnumDefinitions/EnumToString.cpp

@@ -203 +203 @@
                 case SurfaceSlopeYEnum : return "SurfaceSlopeY";
                 case TemperatureEnum : return "Temperature";
+                case TemporaryTemperatureEnum : return "TemporaryTemperature";
                 case TemperatureOldEnum : return "TemperatureOld";
                 case ThicknessEnum : return "Thickness";

issm/trunk/src/c/EnumDefinitions/StringToEnum.cpp

@@ -201 +201 @@
         else if (strcmp(name,"SurfaceSlopeY")==0) return SurfaceSlopeYEnum;
         else if (strcmp(name,"Temperature")==0) return TemperatureEnum;
+        else if (strcmp(name,"TemporaryTemperature")==0) return TemporaryTemperatureEnum;
         else if (strcmp(name,"TemperatureOld")==0) return TemperatureOldEnum;
         else if (strcmp(name,"Thickness")==0) return ThicknessEnum;

issm/trunk/src/c/objects/Elements/Penta.cpp

@@ -5645 +5645 @@
         double       B[numdof];
         double       B_average;
+        bool         converged;
 
         /*Get dof list: */
@@ -5654 +5655 @@
         }
 
-        B_average=Paterson((values[0]+values[1]+values[2]+values[3]+values[4]+values[5])/6.0);
-        for(i=0;i<numdof;i++){
-                B[i]=B_average;
-        }
-
-        /*Add thickness as inputs to the tria element: */
-        this->inputs->AddInput(new PentaVertexInput(TemperatureEnum,values));
-
-        /*Also update the rheology*/
-        this->matice->inputs->AddInput(new PentaVertexInput(RheologyBEnum,B));
+        this->inputs->GetParameterValue(&converged,ConvergedEnum);
+        if(converged){
+                this->inputs->AddInput(new PentaVertexInput(TemperatureEnum,values));
+
+                /*Update Rheology only if convreged (we must make sure that the temperature is below melting point
+                 * otherwise the rheology could be negative*/
+                B_average=Paterson((values[0]+values[1]+values[2]+values[3]+values[4]+values[5])/6.0);
+                for(i=0;i<numdof;i++) B[i]=B_average;
+                this->matice->inputs->AddInput(new PentaVertexInput(RheologyBEnum,B));
+        }
+        else{
+                this->inputs->AddInput(new PentaVertexInput(TemporaryTemperatureEnum,values));
+        }
 
         /*Free ressources:*/
@@ -5760 +5764 @@
                                 name==DragCoefficientEnum ||
                                 name==GradientEnum ||
-                                name==OldGradientEnum 
+                                name==OldGradientEnum  ||
+                                name==ConvergedEnum
                                 ) {
                 return true;

issm/trunk/src/c/objects/Inputs/DoubleInput.cpp

@@ -157 +157 @@
 /*FUNCTION DoubleInput::GetParameterValue(bool* pvalue) {{{1*/
 void DoubleInput::GetParameterValue(bool* pvalue){
-
-        ISSMERROR("cannot return a bool");
+#ifdef _SERIAL_
+        *pvalue=(bool)value;
+#else
+        ISSMERROR("Double input of enum %s cannot return a boolean",EnumToString(enum_type));
+#endif
 
 }

issm/trunk/src/c/objects/Loads/Pengrid.cpp

@@ -19 +19 @@
         
 /*Pengrid constructors and destructor*/
-/*FUNCTION Pengrid::constructor {{{1*/
+/*FUNCTION Pengrid::Pengrid(){{{1*/
 Pengrid::Pengrid(){
         this->inputs=NULL;
@@ -72 +72 @@
 }
 /*}}}*/
-/*FUNCTION Pengrid::destructor {{{1*/
+/*FUNCTION Pengrid::~Pengrid(){{{1*/
 Pengrid::~Pengrid(){
         delete inputs;
@@ -513 +513 @@
         gauss=penta->GaussFromNode(node);
         penta->inputs->GetParameterValue(&pressure,gauss,PressureEnum);
-        penta->inputs->GetParameterValue(&temperature,gauss,TemperatureEnum);
+        penta->inputs->GetParameterValue(&temperature,gauss,TemporaryTemperatureEnum);
 
         //Recover our data:

issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp

@@ -18 +18 @@
         /*intermediary: */
         Mat Kgg=NULL;
-        Mat Kgg_nopenalty=NULL;
         Mat Kff=NULL;
         Mat Kfs=NULL;
         Vec pg=NULL;
-        Vec pg_nopenalty=NULL;
         Vec pf=NULL;
 
-        int converged;
+        bool converged;
         int constraints_converged;
         int num_unstable_constraints;
@@ -45 +43 @@
 
         count=1;
-        converged=0;
+        converged=false;
 
+        if(verbose)_printf_("%s\n","starting direct shooting method");
+        InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,reset_penalties,ResetPenaltiesEnum);
+        InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,false,ConvergedEnum);
         for(;;){
 
-                if(verbose)_printf_("%s\n","starting direct shooting method");
 
-                if(count==1) reset_penalties=1; else reset_penalties=0;
-                InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,reset_penalties,ResetPenaltiesEnum);
+                SystemMatricesx(&Kgg,&pg,&melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
 
-                /*===================== DEBUGING ====================*/
-                if (count==1) SystemMatricesx(&Kgg_nopenalty,&pg_nopenalty,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,true,true,false,false);
-                MatDuplicate(Kgg_nopenalty,MAT_COPY_VALUES,&Kgg);
-                VecDuplicatePatch(&pg,pg_nopenalty);
-                PenaltySystemMatricesx(Kgg, pg,&melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag);
-                /*===================== DEBUGING ====================*/
-                //SystemMatricesx(&Kgg,&pg,&melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
-                /*===================== DEBUGING ====================*/
-
-                Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters);
-
-                /*Free ressources: */
-                MatFree(&Kgg);
-        
+                Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg);
                 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters);
 
-                //no need for pg and Kfs anymore 
                 VecFree(&pg); 
                 MatFree(&Kfs);
 
-                /*Solve: */
                 VecFree(&tf);
                 Solverx(&tf, Kff, pf,tf_old, femmodel->parameters);
                 VecFree(&tf_old); VecDuplicatePatch(&tf_old,tf);
         
-                //no need for Kff and pf anymore
                 MatFree(&Kff);VecFree(&pf);VecFree(&tg);
 
                 Mergesolutionfromftogx(&tg, tf,femmodel->ys,femmodel->nodesets,femmodel->parameters);
-
                 InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,tg);
 
-                //Deal with penalty loads
                 PenaltyConstraintsx(&constraints_converged, &num_unstable_constraints, femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
 
                 if (!converged){
                         if(verbose)_printf_("%s%i\n","   #unstable constraints = ",num_unstable_constraints);
-                        if (num_unstable_constraints <= min_thermal_constraints)converged=1;
+                        if (num_unstable_constraints <= min_thermal_constraints)converged=true;
                 }
                 count++;
+
+                InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,converged,ConvergedEnum);
                 
-                if(converged==1)break;
+                if(converged)break;
         }
 
-        //add melting_offset to inputs:
+        InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,tg);
         InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,melting_offset,MeltingOffsetEnum);
 
         /*Free ressources: */
-        MatFree(&Kgg_nopenalty);
-        VecFree(&pg_nopenalty);
         VecFree(&tg);
         VecFree(&tf);

issm/trunk/src/c/toolkits/petsc/patches/PetscDetermineLocalSize.cpp

@@ -49 +49 @@
 
 }
-
-
- 

issm/trunk/src/m/enum/AbsoluteEnum.m

@@ -9 +9 @@
 %      macro=AbsoluteEnum()
 
-macro=243;
+macro=244;

issm/trunk/src/m/enum/AdjointEnum.m

@@ -9 +9 @@
 %      macro=AdjointEnum()
 
-macro=258;
+macro=259;

issm/trunk/src/m/enum/BetaEnum.m

@@ -9 +9 @@
 %      macro=BetaEnum()
 
-macro=253;
+macro=254;

issm/trunk/src/m/enum/BoolExternalResultEnum.m

@@ -9 +9 @@
 %      macro=BoolExternalResultEnum()
 
-macro=217;
+macro=218;

issm/trunk/src/m/enum/CmGradientEnum.m

@@ -9 +9 @@
 %      macro=CmGradientEnum()
 
-macro=254;
+macro=255;

issm/trunk/src/m/enum/CmJumpEnum.m

@@ -9 +9 @@
 %      macro=CmJumpEnum()
 
-macro=255;
+macro=256;

issm/trunk/src/m/enum/CmMaxEnum.m

@@ -9 +9 @@
 %      macro=CmMaxEnum()
 
-macro=256;
+macro=257;

issm/trunk/src/m/enum/CmMinEnum.m

@@ -9 +9 @@
 %      macro=CmMinEnum()
 
-macro=257;
+macro=258;

issm/trunk/src/m/enum/ColinearEnum.m

@@ -9 +9 @@
 %      macro=ColinearEnum()
 
-macro=251;
+macro=252;

issm/trunk/src/m/enum/ConnectivityEnum.m

@@ -9 +9 @@
 %      macro=ConnectivityEnum()
 
-macro=261;
+macro=262;

issm/trunk/src/m/enum/ControlParameterEnum.m

@@ -9 +9 @@
 %      macro=ControlParameterEnum()
 
-macro=262;
+macro=263;

issm/trunk/src/m/enum/ControlSteadyEnum.m

@@ -9 +9 @@
 %      macro=ControlSteadyEnum()
 
-macro=263;
+macro=264;

issm/trunk/src/m/enum/DakotaParameterEnum.m

@@ -9 +9 @@
 %      macro=DakotaParameterEnum()
 
-macro=264;
+macro=265;

issm/trunk/src/m/enum/DimEnum.m

@@ -9 +9 @@
 %      macro=DimEnum()
 
-macro=265;
+macro=266;

issm/trunk/src/m/enum/DoubleElementResultEnum.m

@@ -9 +9 @@
 %      macro=DoubleElementResultEnum()
 
-macro=214;
+macro=215;

issm/trunk/src/m/enum/DoubleExternalResultEnum.m

@@ -9 +9 @@
 %      macro=DoubleExternalResultEnum()
 
-macro=218;
+macro=219;

issm/trunk/src/m/enum/DoubleMatExternalResultEnum.m

@@ -9 +9 @@
 %      macro=DoubleMatExternalResultEnum()
 
-macro=220;
+macro=221;

issm/trunk/src/m/enum/DoubleVecExternalResultEnum.m

@@ -9 +9 @@
 %      macro=DoubleVecExternalResultEnum()
 
-macro=219;
+macro=220;

issm/trunk/src/m/enum/EnumToString.m

@@ -198 +198 @@
 case SurfaceSlopeYEnum(), string='SurfaceSlopeY'; return
 case TemperatureEnum(), string='Temperature'; return
+case TemporaryTemperatureEnum(), string='TemporaryTemperature'; return
 case TemperatureOldEnum(), string='TemperatureOld'; return
 case ThicknessEnum(), string='Thickness'; return

issm/trunk/src/m/enum/EpsAbsEnum.m

@@ -9 +9 @@
 %      macro=EpsAbsEnum()
 
-macro=266;
+macro=267;

issm/trunk/src/m/enum/EpsCmEnum.m

@@ -9 +9 @@
 %      macro=EpsCmEnum()
 
-macro=267;
+macro=268;

issm/trunk/src/m/enum/EpsRelEnum.m

@@ -9 +9 @@
 %      macro=EpsRelEnum()
 
-macro=268;
+macro=269;

issm/trunk/src/m/enum/EpsResEnum.m

@@ -9 +9 @@
 %      macro=EpsResEnum()
 
-macro=269;
+macro=270;

issm/trunk/src/m/enum/GradientEnum.m

@@ -9 +9 @@
 %      macro=GradientEnum()
 
-macro=259;
+macro=260;

issm/trunk/src/m/enum/GravityEnum.m

@@ -9 +9 @@
 %      macro=GravityEnum()
 
-macro=246;
+macro=247;

issm/trunk/src/m/enum/HeatCapacityEnum.m

@@ -9 +9 @@
 %      macro=HeatCapacityEnum()
 
-macro=270;
+macro=271;

issm/trunk/src/m/enum/IndexEnum.m

@@ -9 +9 @@
 %      macro=IndexEnum()
 
-macro=302;
+macro=303;

issm/trunk/src/m/enum/IndexedEnum.m

@@ -9 +9 @@
 %      macro=IndexedEnum()
 
-macro=300;
+macro=301;

issm/trunk/src/m/enum/IntExternalResultEnum.m

@@ -9 +9 @@
 %      macro=IntExternalResultEnum()
 
-macro=221;
+macro=222;

issm/trunk/src/m/enum/IntersectEnum.m

@@ -9 +9 @@
 %      macro=IntersectEnum()
 
-macro=250;
+macro=251;

issm/trunk/src/m/enum/IsHutterEnum.m

@@ -9 +9 @@
 %      macro=IsHutterEnum()
 
-macro=271;
+macro=272;

issm/trunk/src/m/enum/IsMacAyealPattynEnum.m

@@ -9 +9 @@
 %      macro=IsMacAyealPattynEnum()
 
-macro=272;
+macro=273;

issm/trunk/src/m/enum/IsStokesEnum.m

@@ -9 +9 @@
 %      macro=IsStokesEnum()
 
-macro=273;
+macro=274;

issm/trunk/src/m/enum/JEnum.m

@@ -9 +9 @@
 %      macro=JEnum()
 
-macro=224;
+macro=225;

issm/trunk/src/m/enum/LatentHeatEnum.m

@@ -9 +9 @@
 %      macro=LatentHeatEnum()
 
-macro=274;
+macro=275;

issm/trunk/src/m/enum/LowmemEnum.m

@@ -9 +9 @@
 %      macro=LowmemEnum()
 
-macro=275;
+macro=276;

issm/trunk/src/m/enum/MaxAbsVxEnum.m

@@ -9 +9 @@
 %      macro=MaxAbsVxEnum()
 
-macro=234;
+macro=235;

issm/trunk/src/m/enum/MaxAbsVyEnum.m

@@ -9 +9 @@
 %      macro=MaxAbsVyEnum()
 
-macro=237;
+macro=238;

issm/trunk/src/m/enum/MaxAbsVzEnum.m

@@ -9 +9 @@
 %      macro=MaxAbsVzEnum()
 
-macro=240;
+macro=241;

issm/trunk/src/m/enum/MaxIterEnum.m

@@ -9 +9 @@
 %      macro=MaxIterEnum()
 
-macro=276;
+macro=277;

issm/trunk/src/m/enum/MaxNonlinearIterationsEnum.m

@@ -9 +9 @@
 %      macro=MaxNonlinearIterationsEnum()
 
-macro=277;
+macro=278;

issm/trunk/src/m/enum/MaxVelEnum.m

@@ -9 +9 @@
 %      macro=MaxVelEnum()
 
-macro=231;
+macro=232;

issm/trunk/src/m/enum/MaxVxEnum.m

@@ -9 +9 @@
 %      macro=MaxVxEnum()
 
-macro=233;
+macro=234;

issm/trunk/src/m/enum/MaxVyEnum.m

@@ -9 +9 @@
 %      macro=MaxVyEnum()
 
-macro=236;
+macro=237;

issm/trunk/src/m/enum/MaxVzEnum.m

@@ -9 +9 @@
 %      macro=MaxVzEnum()
 
-macro=239;
+macro=240;

issm/trunk/src/m/enum/MeltingPointEnum.m

@@ -9 +9 @@
 %      macro=MeltingPointEnum()
 
-macro=278;
+macro=279;

issm/trunk/src/m/enum/MinMechanicalConstraintsEnum.m

@@ -9 +9 @@
 %      macro=MinMechanicalConstraintsEnum()
 
-macro=279;
+macro=280;

issm/trunk/src/m/enum/MinThermalConstraintsEnum.m

@@ -9 +9 @@
 %      macro=MinThermalConstraintsEnum()
 
-macro=280;
+macro=281;

issm/trunk/src/m/enum/MinVelEnum.m

@@ -9 +9 @@
 %      macro=MinVelEnum()
 
-macro=230;
+macro=231;

issm/trunk/src/m/enum/MinVxEnum.m

@@ -9 +9 @@
 %      macro=MinVxEnum()
 
-macro=232;
+macro=233;

issm/trunk/src/m/enum/MinVyEnum.m

@@ -9 +9 @@
 %      macro=MinVyEnum()
 
-macro=235;
+macro=236;

issm/trunk/src/m/enum/MinVzEnum.m

@@ -9 +9 @@
 %      macro=MinVzEnum()
 
-macro=238;
+macro=239;

issm/trunk/src/m/enum/MiniEnum.m

@@ -9 +9 @@
 %      macro=MiniEnum()
 
-macro=213;
+macro=214;

issm/trunk/src/m/enum/MixedLayerCapacityEnum.m

@@ -9 +9 @@
 %      macro=MixedLayerCapacityEnum()
 
-macro=248;
+macro=249;

issm/trunk/src/m/enum/NStepsEnum.m

@@ -9 +9 @@
 %      macro=NStepsEnum()
 
-macro=281;
+macro=282;

issm/trunk/src/m/enum/NdtEnum.m

@@ -9 +9 @@
 %      macro=NdtEnum()
 
-macro=282;
+macro=283;

issm/trunk/src/m/enum/NodalEnum.m

@@ -9 +9 @@
 %      macro=NodalEnum()
 
-macro=304;
+macro=305;

issm/trunk/src/m/enum/NumOutputEnum.m

@@ -9 +9 @@
 %      macro=NumOutputEnum()
 
-macro=283;
+macro=284;

issm/trunk/src/m/enum/NumRiftsEnum.m

@@ -9 +9 @@
 %      macro=NumRiftsEnum()
 
-macro=284;
+macro=285;

issm/trunk/src/m/enum/NumberOfElementsEnum.m

@@ -9 +9 @@
 %      macro=NumberOfElementsEnum()
 
-macro=285;
+macro=286;

issm/trunk/src/m/enum/NumberOfNodesEnum.m

@@ -9 +9 @@
 %      macro=NumberOfNodesEnum()
 
-macro=286;
+macro=287;

issm/trunk/src/m/enum/NumberOfVerticesEnum.m

@@ -9 +9 @@
 %      macro=NumberOfVerticesEnum()
 
-macro=287;
+macro=288;

issm/trunk/src/m/enum/OldGradientEnum.m

@@ -9 +9 @@
 %      macro=OldGradientEnum()
 
-macro=260;
+macro=261;

issm/trunk/src/m/enum/OptScalEnum.m

@@ -9 +9 @@
 %      macro=OptScalEnum()
 
-macro=288;
+macro=289;

issm/trunk/src/m/enum/OutputFilePointerEnum.m

@@ -9 +9 @@
 %      macro=OutputFilePointerEnum()
 
-macro=289;
+macro=290;

issm/trunk/src/m/enum/OutputFrequencyEnum.m

@@ -9 +9 @@
 %      macro=OutputFrequencyEnum()
 
-macro=229;
+macro=230;

issm/trunk/src/m/enum/P0Enum.m

@@ -9 +9 @@
 %      macro=P0Enum()
 
-macro=210;
+macro=211;

issm/trunk/src/m/enum/P1DGEnum.m

@@ -9 +9 @@
 %      macro=P1DGEnum()
 
-macro=212;
+macro=213;

issm/trunk/src/m/enum/P1Enum.m

@@ -9 +9 @@
 %      macro=P1Enum()
 
-macro=211;
+macro=212;

issm/trunk/src/m/enum/ParameterOutputEnum.m

@@ -9 +9 @@
 %      macro=ParameterOutputEnum()
 
-macro=290;
+macro=291;

issm/trunk/src/m/enum/PatchEnum.m

@@ -9 +9 @@
 %      macro=PatchEnum()
 
-macro=225;
+macro=226;

issm/trunk/src/m/enum/PatchNodesEnum.m

@@ -9 +9 @@
 %      macro=PatchNodesEnum()
 
-macro=227;
+macro=228;

issm/trunk/src/m/enum/PatchVerticesEnum.m

@@ -9 +9 @@
 %      macro=PatchVerticesEnum()
 
-macro=226;
+macro=227;

issm/trunk/src/m/enum/PenaltyMeltingEnum.m

@@ -9 +9 @@
 %      macro=PenaltyMeltingEnum()
 
-macro=291;
+macro=292;

issm/trunk/src/m/enum/PentaVertexElementResultEnum.m

@@ -9 +9 @@
 %      macro=PentaVertexElementResultEnum()
 
-macro=216;
+macro=217;

issm/trunk/src/m/enum/PetscVecExternalResultEnum.m

@@ -9 +9 @@
 %      macro=PetscVecExternalResultEnum()
 
-macro=222;
+macro=223;

issm/trunk/src/m/enum/QmuAnalysisEnum.m

@@ -9 +9 @@
 %      macro=QmuAnalysisEnum()
 
-macro=292;
+macro=293;

issm/trunk/src/m/enum/QmuErrNameEnum.m

@@ -9 +9 @@
 %      macro=QmuErrNameEnum()
 
-macro=293;
+macro=294;

issm/trunk/src/m/enum/QmuInNameEnum.m

@@ -9 +9 @@
 %      macro=QmuInNameEnum()
 
-macro=294;
+macro=295;

issm/trunk/src/m/enum/QmuMassFluxSegmentsEnum.m

@@ -9 +9 @@
 %      macro=QmuMassFluxSegmentsEnum()
 
-macro=295;
+macro=296;

issm/trunk/src/m/enum/QmuNPartEnum.m

@@ -9 +9 @@
 %      macro=QmuNPartEnum()
 
-macro=296;
+macro=297;

issm/trunk/src/m/enum/QmuNumberOfResponsesEnum.m

@@ -9 +9 @@
 %      macro=QmuNumberOfResponsesEnum()
 
-macro=299;
+macro=300;

issm/trunk/src/m/enum/QmuOutNameEnum.m

@@ -9 +9 @@
 %      macro=QmuOutNameEnum()
 
-macro=297;
+macro=298;

issm/trunk/src/m/enum/QmuPartEnum.m

@@ -9 +9 @@
 %      macro=QmuPartEnum()
 
-macro=298;
+macro=299;

issm/trunk/src/m/enum/QmuVxEnum.m

@@ -9 +9 @@
 %      macro=QmuVxEnum()
 
-macro=198;
+macro=199;

issm/trunk/src/m/enum/QmuVyEnum.m

@@ -9 +9 @@
 %      macro=QmuVyEnum()
 
-macro=203;
+macro=204;

issm/trunk/src/m/enum/QmuVzEnum.m

@@ -9 +9 @@
 %      macro=QmuVzEnum()
 
-macro=208;
+macro=209;

issm/trunk/src/m/enum/RegularEnum.m

@@ -9 +9 @@
 %      macro=RegularEnum()
 
-macro=301;
+macro=302;

issm/trunk/src/m/enum/RelativeEnum.m

@@ -9 +9 @@
 %      macro=RelativeEnum()
 
-macro=241;
+macro=242;

issm/trunk/src/m/enum/ResidualEnum.m

@@ -9 +9 @@
 %      macro=ResidualEnum()
 
-macro=242;
+macro=243;

issm/trunk/src/m/enum/ResponseDescriptorsEnum.m

@@ -9 +9 @@
 %      macro=ResponseDescriptorsEnum()
 
-macro=305;
+macro=306;

issm/trunk/src/m/enum/RhoIceEnum.m

@@ -9 +9 @@
 %      macro=RhoIceEnum()
 
-macro=244;
+macro=245;

issm/trunk/src/m/enum/RhoWaterEnum.m

@@ -9 +9 @@
 %      macro=RhoWaterEnum()
 
-macro=245;
+macro=246;

issm/trunk/src/m/enum/ScaledEnum.m

@@ -9 +9 @@
 %      macro=ScaledEnum()
 
-macro=303;
+macro=304;

issm/trunk/src/m/enum/SeparateEnum.m

@@ -9 +9 @@
 %      macro=SeparateEnum()
 
-macro=252;
+macro=253;

issm/trunk/src/m/enum/SolverStringEnum.m

@@ -9 +9 @@
 %      macro=SolverStringEnum()
 
-macro=306;
+macro=307;

issm/trunk/src/m/enum/SparsityEnum.m

@@ -9 +9 @@
 %      macro=SparsityEnum()
 
-macro=307;
+macro=308;

issm/trunk/src/m/enum/StringExternalResultEnum.m

@@ -9 +9 @@
 %      macro=StringExternalResultEnum()
 
-macro=223;
+macro=224;

issm/trunk/src/m/enum/StringToEnum.m

@@ -196 +196 @@
                         elseif (strcmpi(name,'SurfaceSlopeY')), enum=SurfaceSlopeYEnum(); return
                         elseif (strcmpi(name,'Temperature')), enum=TemperatureEnum(); return
+                        elseif (strcmpi(name,'TemporaryTemperature')), enum=TemporaryTemperatureEnum(); return
                         elseif (strcmpi(name,'TemperatureOld')), enum=TemperatureOldEnum(); return
                         elseif (strcmpi(name,'Thickness')), enum=ThicknessEnum(); return

issm/trunk/src/m/enum/TemperatureOldEnum.m

@@ -9 +9 @@
 %      macro=TemperatureOldEnum()
 
-macro=187;
+macro=188;

issm/trunk/src/m/enum/ThermalConductivityEnum.m

@@ -9 +9 @@
 %      macro=ThermalConductivityEnum()
 
-macro=247;
+macro=248;

issm/trunk/src/m/enum/ThermalExchangeVelocityEnum.m

@@ -9 +9 @@
 %      macro=ThermalExchangeVelocityEnum()
 
-macro=249;
+macro=250;

issm/trunk/src/m/enum/ThicknessEnum.m

@@ -9 +9 @@
 %      macro=ThicknessEnum()
 
-macro=188;
+macro=189;

issm/trunk/src/m/enum/ThicknessObsEnum.m

@@ -9 +9 @@
 %      macro=ThicknessObsEnum()
 
-macro=189;
+macro=190;

issm/trunk/src/m/enum/TimeEnum.m

@@ -9 +9 @@
 %      macro=TimeEnum()
 
-macro=228;
+macro=229;

issm/trunk/src/m/enum/TolXEnum.m

@@ -9 +9 @@
 %      macro=TolXEnum()
 
-macro=308;
+macro=309;

issm/trunk/src/m/enum/TriaVertexElementResultEnum.m

@@ -9 +9 @@
 %      macro=TriaVertexElementResultEnum()
 
-macro=215;
+macro=216;

issm/trunk/src/m/enum/TypeEnum.m

@@ -9 +9 @@
 %      macro=TypeEnum()
 
-macro=190;
+macro=191;

issm/trunk/src/m/enum/VariableDescriptorsEnum.m

@@ -9 +9 @@
 %      macro=VariableDescriptorsEnum()
 
-macro=309;
+macro=310;

issm/trunk/src/m/enum/VelEnum.m

@@ -9 +9 @@
 %      macro=VelEnum()
 
-macro=191;
+macro=192;

issm/trunk/src/m/enum/VelObsEnum.m

@@ -9 +9 @@
 %      macro=VelObsEnum()
 
-macro=192;
+macro=193;

issm/trunk/src/m/enum/VerboseEnum.m

@@ -9 +9 @@
 %      macro=VerboseEnum()
 
-macro=310;
+macro=311;

issm/trunk/src/m/enum/ViscosityOvershootEnum.m

@@ -9 +9 @@
 %      macro=ViscosityOvershootEnum()
 
-macro=193;
+macro=194;

issm/trunk/src/m/enum/VxAverageEnum.m

@@ -9 +9 @@
 %      macro=VxAverageEnum()
 
-macro=194;
+macro=195;

issm/trunk/src/m/enum/VxEnum.m

@@ -9 +9 @@
 %      macro=VxEnum()
 
-macro=195;
+macro=196;

issm/trunk/src/m/enum/VxObsEnum.m

@@ -9 +9 @@
 %      macro=VxObsEnum()
 
-macro=196;
+macro=197;

issm/trunk/src/m/enum/VxOldEnum.m

@@ -9 +9 @@
 %      macro=VxOldEnum()
 
-macro=197;
+macro=198;

issm/trunk/src/m/enum/VyAverageEnum.m

@@ -9 +9 @@
 %      macro=VyAverageEnum()
 
-macro=199;
+macro=200;

issm/trunk/src/m/enum/VyEnum.m

@@ -9 +9 @@
 %      macro=VyEnum()
 
-macro=200;
+macro=201;

issm/trunk/src/m/enum/VyObsEnum.m

@@ -9 +9 @@
 %      macro=VyObsEnum()
 
-macro=201;
+macro=202;

issm/trunk/src/m/enum/VyOldEnum.m

@@ -9 +9 @@
 %      macro=VyOldEnum()
 
-macro=202;
+macro=203;

issm/trunk/src/m/enum/VzAverageEnum.m

@@ -9 +9 @@
 %      macro=VzAverageEnum()
 
-macro=204;
+macro=205;

issm/trunk/src/m/enum/VzEnum.m

@@ -9 +9 @@
 %      macro=VzEnum()
 
-macro=205;
+macro=206;

issm/trunk/src/m/enum/VzObsEnum.m

@@ -9 +9 @@
 %      macro=VzObsEnum()
 
-macro=206;
+macro=207;

issm/trunk/src/m/enum/VzOldEnum.m

@@ -9 +9 @@
 %      macro=VzOldEnum()
 
-macro=207;
+macro=208;

issm/trunk/src/m/enum/WaitOnLockEnum.m

@@ -9 +9 @@
 %      macro=WaitOnLockEnum()
 
-macro=311;
+macro=312;

issm/trunk/src/m/enum/WeightsEnum.m

@@ -9 +9 @@
 %      macro=WeightsEnum()
 
-macro=209;
+macro=210;

issm/trunk/src/m/enum/YtsEnum.m

@@ -9 +9 @@
 %      macro=YtsEnum()
 
-macro=315;
+macro=316;

issm/trunk/src/m/solvers/solver_thermal_nonlinear.m

@@ -13 +13 @@
 
         displaystring(femmodel.parameters.Verbose,'\n%s',['   starting direct shooting method']);
+
+        %Reset penalties and initialize convergence as false
+        [femmodel.elements loads]=InputUpdateFromConstant(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,0,ConvergedEnum);
+        [femmodel.elements femmodel.loads]=InputUpdateFromConstant(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,1,ResetPenaltiesEnum);
+
         while(~converged),
 
-                %Update inputs in datasets
-                if count==1 reset_penalties=0; else reset_penalties=1; end
-                [femmodel.elements femmodel.loads]=InputUpdateFromConstant(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,reset_penalties,ResetPenaltiesEnum);
-
-                % ================= DEBUG FOR NOW ====================
-                if count==1 
-                        [K_gg_nopenalty,p_g_nopenalty,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,true,true,false,false);
-                end
-                [K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-                % ================= DEBUG FOR NOW ====================
-                %[K_gg,p_g,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-                % ================= DEBUG FOR NOW ====================
+                [K_gg,p_g,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
 
                 [K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
@@ -33 +27 @@
                 displaystring(femmodel.parameters.Verbose,'%s%g','   condition number of stiffness matrix: ',condest(K_ff));
                 t_f=Solver(K_ff,p_f,[],femmodel.parameters);
+                t_g= Mergesolutionfromftog( t_f, femmodel.ys, femmodel.nodesets,femmodel.parameters); 
 
-                t_g= Mergesolutionfromftog( t_f, femmodel.ys, femmodel.nodesets,femmodel.parameters); 
                 [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,t_g);
-
                 [femmodel.loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads, femmodel.materials,femmodel.parameters);
         
                 if ~converged,
                         displaystring(femmodel.parameters.Verbose,'%s%i','   #unstable constraints ',num_unstable_constraints);
-                        
                         if num_unstable_constraints<=femmodel.parameters.MinThermalConstraints,
                                 converged=1;
                         end
                 end
-
+                [femmodel.elements loads]=InputUpdateFromConstant(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,double(converged),ConvergedEnum);
                 count=count+1;
         end
 
-        %add melting_offset  into inputs
+        [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,t_g);
         [femmodel.elements femmodel.loads]=InputUpdateFromConstant(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,melting_offset,MeltingOffsetEnum);
 end
-