Changeset 27518 for issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp

Timestamp:

01/19/23 03:50:29 (2 years ago)

Author:

bdef

Message:

BUG: fix to the split between layers of moulin input for DoCo

File:

• issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp (modified) (5 diffs)

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

@@ -297 +297 @@
         IssmDouble dt,scalar,water_head;
         IssmDouble water_load,transfer,runoff_value;
-        IssmDouble epl_storing,epl_transmitivity;
+        IssmDouble epl_storing;  //,epl_transmitivity;
         IssmDouble Jdet,time;
         IssmDouble residual,connectivity;
+        IssmDouble active_node;
 
         IssmDouble *xyz_list             = NULL;
@@ -313 +314 @@
         ElementVector* pe    = basalelement->NewElementVector();
         IssmDouble*    basis = xNew<IssmDouble>(numnodes);
+
 
         /*Retrieve all inputs and parameters*/
@@ -353 +355 @@
                 basalelement ->JacobianDeterminant(&Jdet,xyz_list,gauss);
                 basalelement ->NodalFunctions(basis,gauss);
+                //epl_transmitivity = EplTransmitivity(basalelement,gauss,epl_thick_input);
                 epl_storing     = EplStoring(basalelement,gauss,epl_thick_input);
-                epl_transmitivity = EplTransmitivity(basalelement,gauss,epl_thick_input);
 
                 /*Loading term*/
@@ -362 +364 @@
                 scalar = Jdet*gauss->weight*(water_load+runoff_value);
                 if(dt!=0.) scalar = scalar*dt;
-                for(int i=0;i<numnodes;i++)pe->values[i]+=scalar*basis[i];
-
+                for(int i=0;i<numnodes;i++){
+                        //This is the original
+                        pe->values[i]+=scalar*basis[i];
+                        //This is the noded version
+                        /* basalelement->GetInputValue(&active_node,basalelement->nodes[i],HydrologydcMaskEplactiveNodeEnum); */
+                        /* if(!reCast<bool>(active_node)){ */
+                        /*      pe->values[i]+=scalar*basis[i]; */
+                        //}
+                        //if(basalelement->nodes[i]->Sid()==42)_printf_("EPL uni Input "<<scalar*basis[i]<<"\n");
+                }
                 /*Transient and transfer terms*/
                 if(dt!=0.){
@@ -379 +389 @@
         for(int iv=0;iv<numvertices;iv++){
                 gauss->GaussVertex(iv);
-                epl_transmitivity = EplTransmitivity(basalelement,gauss,epl_thick_input);
+                //epl_transmitivity = EplTransmitivity(basalelement,gauss,epl_thick_input);
                 connectivity = IssmDouble(basalelement->VertexConnectivity(iv));
                 residual_input->GetInputValue(&residual,gauss);