Changeset 15273

Timestamp:

06/17/13 08:33:26 (12 years ago)

Author:

bdef

Message:

CHG: Minor changes in the hydrological solver

Location:

issm/trunk-jpl/src

Files:

• c/classes/Elements/Tria.cpp (modified) (6 diffs)
• c/solutionsequences/solutionsequence_hydro_nonlinear.cpp (modified) (2 diffs)
• m/classes/hydrologydc.m (modified) (3 diffs)

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

@@ -2 +2 @@
  * \brief: implementation of the Tria object
  */
-
 /*Headers:*/
 /*{{{*/
@@ -22 +21 @@
 /*Element macros*/
 #define NUMVERTICES 3
-
 /*Constructors/destructor/copy*/
 /*FUNCTION Tria::Tria(){{{*/
@@ -6506 +6504 @@
         IssmDouble wh_trans;
         IssmDouble activeEpl[numdof];
-        IssmDouble storing[numdof];
+        IssmDouble eplstoring[numdof],sedstoring[numdof];
         IssmDouble epl_head[numdof],sed_head[numdof];
 
@@ -6530 +6528 @@
 
                         this->parameters->FindParam(&leakage,HydrologydcLeakageFactorEnum);
-                        
+
                         sed_trans = matpar->GetSedimentTransmitivity();
                         sed_thick = matpar->GetSedimentThickness();
@@ -6536 +6534 @@
                         for(int i=0;i<numdof;i++){
                                 if(activeEpl[i]==0.0)continue;
-                
+
+                                eplstoring[i]=matpar->GetEplStoring();          
+                                sedstoring[i]=matpar->GetSedimentStoring();
+
                                 if(sed_head[i]>epl_head[i]){
-                                        storing[i]=matpar->GetSedimentStoring();
-                                        wh_trans=sed_trans*storing[i]*(epl_head[i]-sed_head[i])/(leakage*sed_thick);
+                                        wh_trans=sed_trans*sedstoring[i]*(epl_head[i]-sed_head[i])/(leakage*sed_thick);                         
                                 }
                                 else{
                                         this->GetHydrologyDCInefficientHmax(&h_max,nodes[i]);
-                                        storing[i]=matpar->GetEplStoring();
-                                        wh_trans=sed_trans*storing[i]*(epl_head[i]-sed_head[i])/(leakage*sed_thick);
-                                        if(sed_head[i]>h_max){
+                                        wh_trans=sed_trans*eplstoring[i]*(epl_head[i]-sed_head[i])/(leakage*sed_thick);
+                                        if(sed_head[i]>=h_max){
                                                 wh_trans=0.0;
                                         }
-                                        if((sed_head[i]+wh_trans)>h_max){
-                                                wh_trans=h_max-sed_head[i];
-                                        }
+                                        //                                      if((sed_head[i]+(wh_trans/sedstoring[i]))>h_max){
+                                        //(h_max-sed_head[i])*eplstoring[i];
+                                        //}
                                 }
                                 /*Assign output pointer*/
@@ -6623 +6622 @@
                 this->GetHydrologyDCInefficientHmax(&h_max,nodes[i]);
                 if(eplhead[i]>=h_max){
-                        //vec_mask->SetValue(nodes[i]->Sid(),1.,INS_VAL);
-                        //for(j=0;j<numdof;j++){
-
-                        //      /*Increase of the domain is on the downstream node in term of sediment head*/
-                        //      if(sedhead[j] == sedheadmin){
-                        //              vec_mask->SetValue(nodes[j]->Sid(),1.,INS_VAL);
-                        //              break;
-                        //      }
-                        //}
+                        vec_mask->SetValue(nodes[i]->Sid(),1.,INS_VAL);
+                        for(j=0;j<numdof;j++){
+                                
+                                /*Increase of the domain is on the downstream node in term of sediment head*/
+                                if(sedhead[j] == sedheadmin){
+                                        vec_mask->SetValue(nodes[j]->Sid(),1.,INS_VAL);
+                                        break;
+                                }
+                        }
                 }
         }

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

@@ -115 +115 @@
                         /*Iteration on the EPL layer*/
                         for(;;){
-                                femmodel->HydrologyTransferx();
+                                //                              femmodel->HydrologyTransferx();
                                 femmodel->SystemMatricesx(&Kff, &Kfs, &pf,&df,NULL);
                                 CreateNodalConstraintsx(&ys,femmodel->nodes,HydrologyDCEfficientAnalysisEnum);
@@ -127 +127 @@
                                 Mergesolutionfromftogx(&ug_epl,uf,ys,femmodel->nodes,femmodel->parameters); delete ys;
                                 InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,ug_epl);
-                                //femmodel->HydrologyEPLupdateDomainx();
+                                femmodel->HydrologyEPLupdateDomainx();
                                 ConstraintsStatex(&constraints_converged,&num_unstable_constraints, femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
                                 

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

@@ -27 +27 @@
                 epl_thickness            = 0;
                 epl_transmitivity        = 0;
-
   end
         methods
-                % {{{function obj = hydrologydc(varargin)
-                function obj = hydrologydc(varargin) 
-
+                function obj = hydrologydc(varargin) % {{{ 
                         switch nargin
                                 case 0
@@ -66 +63 @@
                 end 
                 % }}}
-                % {{{function md = checkconsistency(obj,md,solution,analyses) 
+                % {{{ function md = checkconsistency(obj,md,solution,analyses) 
                 function md = checkconsistency(obj,md,solution,analyses) 
                 %Early return
@@ -102 +99 @@
                 end 
                 % }}}
-                % {{{           function disp(obj)
+                % {{{ function disp(obj)
                 function disp(obj) 
                         disp(sprintf('   hydrology Dual Porous Continuum Equivalent parameters:'));