Changeset 24652

Timestamp:

03/11/20 11:55:27 (5 years ago)

Author:

Mathieu Morlighem

Message:

CHG: fixing PDD Sicopolis for P1DG

File:

• issm/trunk-jpl/src/c/classes/Elements/Element.cpp (modified) (3 diffs)

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

@@ -2986 +2986 @@
         const int NUM_VERTICES_MONTHS_PER_YEAR  = NUM_VERTICES * 12;
 
-        int             i;
+        int             i,vertexlids[MAXVERTICES];;
         IssmDouble* smb=xNew<IssmDouble>(NUM_VERTICES);         // surface mass balance
         IssmDouble* melt=xNew<IssmDouble>(NUM_VERTICES);                // melting comp. of surface mass balance
@@ -3003 +3003 @@
         IssmDouble inv_twelve=1./12.;                                                           //factor for monthly average
         IssmDouble time,yts,time_yr;
+
+        /*Get vertex Lids for later*/
+        this->GetVerticesLidList(&vertexlids[0]);
 
         /*Get material parameters :*/
@@ -3089 +3092 @@
                                         /*FIXME: look at other function Element::PositiveDegreeDay and propagate change! Just assert for now*/
                                         PentaInput2* temp_input = xDynamicCast<PentaInput2*>(this->GetInput2(TemperatureEnum)); _assert_(temp_input);
-                                        _assert_(temp_input->GetInputInterpolationType()==P1Enum);
-
-                                        GetInputListOnVertices(&s[0],TemperatureEnum);
-                                        yearlytemperatures[0] = s[0];
-                                        yearlytemperatures[1] = s[1];
-                                        yearlytemperatures[2] = s[2];
-                                        this->AddInput2(TemperatureEnum,&yearlytemperatures[0],P1Enum);
+                                        switch(temp_input->GetInputInterpolationType()){
+                                                case P1Enum:
+                                                        temp_input->element_values[3] = yearlytemperatures[3];
+                                                        temp_input->element_values[4] = yearlytemperatures[4];
+                                                        temp_input->element_values[5] = yearlytemperatures[5];
+                                                        temp_input->SetInput(P1Enum,NUM_VERTICES,&vertexlids[0],temp_input->element_values);
+                                                        break;
+                                                case P1DGEnum:
+                                                case P1xP2Enum:
+                                                case P1xP3Enum:
+                                                case P1xP4Enum:
+                                                        temp_input->element_values[3] = yearlytemperatures[3];
+                                                        temp_input->element_values[4] = yearlytemperatures[4];
+                                                        temp_input->element_values[5] = yearlytemperatures[5];
+                                                        temp_input->SetInput(temp_input->GetInputInterpolationType(),this->lid,this->GetNumberOfNodes(temp_input->GetInputInterpolationType()),temp_input->element_values);
+                                                        break;
+                                                default:
+                                                        _error_("Interpolation "<<EnumToStringx(temp_input->GetInputInterpolationType())<<" not supported yet");
+                                        }
+
                                         if(isenthalpy){
                                                 /*Convert that to enthalpy for the enthalpy model*/
-                                                IssmDouble enthalpy[6];
-                                                GetInputListOnVertices(&enthalpy[0],EnthalpyEnum);
-                                                ThermalToEnthalpy(&enthalpy[3],yearlytemperatures[3],0.,0.);
-                                                ThermalToEnthalpy(&enthalpy[4],yearlytemperatures[4],0.,0.);
-                                                ThermalToEnthalpy(&enthalpy[5],yearlytemperatures[5],0.,0.);
-                                                this->AddInput2(EnthalpyEnum,&enthalpy[0],P1Enum);
+                                                PentaInput2* enth_input = xDynamicCast<PentaInput2*>(this->GetInput2(EnthalpyEnum)); _assert_(enth_input);
+                                                switch(enth_input->GetInputInterpolationType()){
+                                                        case P1Enum:
+                                                                ThermalToEnthalpy(&enth_input->element_values[3],yearlytemperatures[3],0.,0.);
+                                                                ThermalToEnthalpy(&enth_input->element_values[4],yearlytemperatures[4],0.,0.);
+                                                                ThermalToEnthalpy(&enth_input->element_values[5],yearlytemperatures[5],0.,0.);
+                                                                enth_input->SetInput(P1Enum,NUM_VERTICES,&vertexlids[0],enth_input->element_values);
+                                                                break;
+                                                        case P1DGEnum:
+                                                        case P1xP2Enum:
+                                                        case P1xP3Enum:
+                                                        case P1xP4Enum:
+                                                                ThermalToEnthalpy(&enth_input->element_values[3],yearlytemperatures[3],0.,0.);
+                                                                ThermalToEnthalpy(&enth_input->element_values[4],yearlytemperatures[4],0.,0.);
+                                                                ThermalToEnthalpy(&enth_input->element_values[5],yearlytemperatures[5],0.,0.);
+                                                                enth_input->SetInput(enth_input->GetInputInterpolationType(),this->lid,this->GetNumberOfNodes(enth_input->GetInputInterpolationType()),enth_input->element_values);
+                                                                break;
+                                                        default:
+                                                                _error_("Interpolation "<<EnumToStringx(temp_input->GetInputInterpolationType())<<" not supported yet");
+                                                }
                                         }
                                 }