Changeset 23528

Timestamp:

12/10/18 04:07:14 (6 years ago)

Author:

bdef

Message:

CHG:minor hydro cleanup and verbose

Location:

issm/trunk-jpl/src/c

Files:

• analyses/HydrologyDCEfficientAnalysis.cpp (modified) (2 diffs)
• analyses/HydrologyDCEfficientAnalysis.h (modified) (3 diffs)
• classes/FemModel.cpp (modified) (19 diffs)
• cores/hydrology_core.cpp (modified) (4 diffs)

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

@@ -653 +653 @@
 }/*}}}*/
 
-void  HydrologyDCEfficientAnalysis::HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, Vector<IssmDouble>* recurence, int* eplzigzag_counter, Element* element){
+void  HydrologyDCEfficientAnalysis::HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, Vector<IssmDouble>* recurence, Element* element){
 
         bool        active_element;
@@ -709 +709 @@
                         }
                 }
-                if (old_active[i]==0.){
+                else if (old_active[i]==0.){
                         /*Activate if we have a residual from sediment*/
                         if(residual[i]>0.){

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

@@ -1 @@
-/*! \file HydrologyDCEfficientAnalysis.h 
+/*! \file HydrologyDCEfficientAnalysis.h
  *  \brief: header file for generic external result object
  */
@@ -22 +22 @@
                 void InitZigZagCounter(FemModel* femmodel);
                 void ResetCounter(FemModel* femmodel);
-                        
+
                 /*Finite element Analysis*/
                 void           Core(FemModel* femmodel);
@@ -40 +40 @@
                 IssmDouble EplStoring(Element* element,Gauss* gauss, Input* epl_thick_input, Input* epl_head_input, Input* base_input);
                 IssmDouble EplTransmitivity(Element* element,Gauss* gauss, Input* epl_thick_input, Input* epl_head_input, Input* base_input);
-                void HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, Vector<IssmDouble>* recurence, int* eplzigzag_counter, Element* element);
+                void HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, Vector<IssmDouble>* recurence, Element* element);
                 void HydrologyEPLGetActive(Vector<IssmDouble>* active_vec, Element* element);
                 void GetHydrologyDCInefficientHmax(IssmDouble* ph_max,Element* element, Node* innode);

issm/trunk-jpl/src/c/classes/FemModel.cpp

@@ -215 +215 @@
 
         /*Checks in debugging mode*/
-        _assert_(this->analysis_type_list); 
+        _assert_(this->analysis_type_list);
 
         /*Find analysis in list*/
@@ -929 +929 @@
                 if(profiler->Used(ESACORE)) _printf0_("   "<<setw(40)<<left<<"ESA core elapsed time:"<<setw(7)<<setprecision(6)<<profiler->TotalTime(ESACORE) << " sec\n");
                 if(profiler->Used(SLRCORE)) _printf0_("   "<<setw(40)<<left<<"SLR core elapsed time:"<<setw(7)<<setprecision(6)<<profiler->TotalTime(SLRCORE) << " sec\n");
-                if(profiler->Used(MPISERIAL)) _printf0_("   "<<setw(40)<<left<<"MPISERIAL elapsed time:"<<setw(7)<<setprecision(6)<<profiler->TotalTime(MPISERIAL) << " sec\n");                
-                
+                if(profiler->Used(MPISERIAL)) _printf0_("   "<<setw(40)<<left<<"MPISERIAL elapsed time:"<<setw(7)<<setprecision(6)<<profiler->TotalTime(MPISERIAL) << " sec\n");
+
                 /*Add to results: */
                 results->AddObject(new GenericExternalResult<IssmDouble>(results->Size()+1, ProfilingSolutionTimeEnum,  solution_time));
@@ -1188 +1188 @@
 
                 /*Update signed distance*/
-                _assert_(vertex->lid<vertices->Size()); 
+                _assert_(vertex->lid<vertices->Size());
                 distances[vertex->lid] = sqrt(dmin);
         }
@@ -2539 +2539 @@
 
         int index=AnalysisIndex(config_type);
-        _assert_(this->analysis_type_list[index]==config_type); 
+        _assert_(this->analysis_type_list[index]==config_type);
 
         /*start module: */
@@ -2719 +2719 @@
         /*Finally: interpolate all inputs and insert them into the new elements.*/
         this->InterpolateInputs(new_vertices,new_elements);
-        
+
         /*Delete old structure and set new pointers*/
         delete this->vertices;          this->vertices          = new_vertices;
@@ -3383 +3383 @@
 /*}}}*/
 void FemModel::GetMesh(int** elementslist, IssmDouble** x, IssmDouble** y, int* numberofvertices, int* numberofelements){/*{{{*/
-        
-        int amrtype;    
+
+        int amrtype;
         this->parameters->FindParam(&amrtype,AmrTypeEnum);
-   
+
         switch(amrtype){
 
@@ -3401 +3401 @@
 }/*}}}*/
 void FemModel::SetMesh(int** elementslist, IssmDouble** x, IssmDouble** y, int* numberofvertices, int* numberofelements){/*{{{*/
-        
-        int amrtype;    
+
+        int amrtype;
         this->parameters->FindParam(&amrtype,AmrTypeEnum);
-   
+
         switch(amrtype){
 
@@ -3566 +3566 @@
                 }
         }
-        
+
         /*Cleanup*/
         xDelete<IssmDouble>(spc);
@@ -4644 +4644 @@
         for (int i=0;i<elements->Size();i++){
                 Element* element=xDynamicCast<Element*>(elements->GetObjectByOffset(i));
-                effanalysis->HydrologyEPLGetMask(mask,recurence,eplzigzag_counter,element);
+                effanalysis->HydrologyEPLGetMask(mask,recurence,element);
         }
         /*check for changes and increment zigzag counter, change the mask if necessary*/
@@ -4652 +4652 @@
                 Node* node=xDynamicCast<Node*>(nodes->GetObjectByOffset(i));
                 if(serial_rec[node->Sid()]==1.)eplzigzag_counter[node->Lid()] ++;
-                if(eplzigzag_counter[node->Lid()]>eplflip_lock & eplflip_lock!=0){
+                if(eplzigzag_counter[node->Lid()]>eplflip_lock && eplflip_lock!=0){
                         mask->SetValue(node->Sid(),old_active[node->Sid()],INS_VAL);
                 }
@@ -5146 +5146 @@
    Vector<IssmDouble>* vminvertexdistance = new Vector<IssmDouble>(numberofvertices);
         IssmDouble* pminvertexdistance                  = NULL;
-        IssmDouble* levelset_points                             = NULL; 
+        IssmDouble* levelset_points                             = NULL;
         IssmDouble x,y;
         IssmDouble threshold,resolution;
@@ -5174 +5174 @@
       sid= vertex->Sid();
                 minvertexdistance=INFINITY;
-      
+
                 /*Find the minimum vertex distance*/
                 for(int j=0;j<numberofpoints;j++){
@@ -5181 +5181 @@
         }
                 /*Now, insert in the vector*/
-                vminvertexdistance->SetValue(sid,minvertexdistance,INS_VAL);    
+                vminvertexdistance->SetValue(sid,minvertexdistance,INS_VAL);
    }
         /*Assemble*/
@@ -5283 +5283 @@
                                 else hmaxvertices[vid]=min(max(maxlength[i]/2.,resolution),hmaxvertices[vid]);//Try first dividing the element
                         }
-                } 
+                }
                 else {
                         /*Try unrefine the element*/
@@ -5295 +5295 @@
                                         }
                                 }
-                        } 
+                        }
                 }
         }
@@ -5306 +5306 @@
 /*}}}*/
 void FemModel::GetVerticeDistanceToZeroLevelSet(IssmDouble** pverticedistance,int levelset_type){/*{{{*/
-        
+
         //itapopo esse metodo pode ser deletado
 
@@ -5325 +5325 @@
    IssmDouble* x                                        = NULL;
    IssmDouble* y                                        = NULL;
-        int* elementslist                               = NULL; 
+        int* elementslist                               = NULL;
         int numberofpoints;
         IssmDouble distance;
@@ -5382 +5382 @@
                 this->amr->ExecuteRefinement(gl_distance,if_distance,deviatoricerror,thicknesserror,
                                                                                                 &newdatalist,&newxylist,&newelementslist);
-                if(newdatalist[0]<=0 || newdatalist[1]<=0) _error_("Error in the ReMeshNeopz.");        
+                if(newdatalist[0]<=0 || newdatalist[1]<=0) _error_("Error in the ReMeshNeopz.");
         }
 
@@ -5459 +5459 @@
         this->parameters->FindParam(&this->amr->deviatoricerror_groupthreshold,AmrDeviatoricErrorGroupThresholdEnum);
         this->parameters->FindParam(&this->amr->deviatoricerror_maximum,AmrDeviatoricErrorMaximumEnum);
-        
+
         /*Initialize NeoPZ data structure*/
         this->amr->SetMesh(&elements,&x,&y,&numberofvertices,&numberofelements);

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

@@ -50 +50 @@
                 /*intermediary: */
                 bool       isefficientlayer;
-                int        step,hydroslices;
+                int        hydrostep,hydroslices;
                 IssmDouble time,init_time,hydrotime,yts;
                 IssmDouble dt,hydrodt;
@@ -56 +56 @@
                 femmodel->parameters->FindParam(&isefficientlayer,HydrologydcIsefficientlayerEnum);
                 femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum);
-                femmodel->parameters->FindParam(&step,StepEnum);
                 femmodel->parameters->FindParam(&time,TimeEnum);
                 femmodel->parameters->FindParam(&hydroslices,HydrologyStepsPerStepEnum);
@@ -68 +67 @@
                         hydrotime=init_time;
                         hydrodt=dt/hydroslices; //computing hydro dt from dt and a divider
+                        hydrostep=0;
                         femmodel->parameters->AddObject(new DoubleParam(HydrologydtEnum,hydrodt));
                         if(hydroslices>1){
@@ -77 +77 @@
                                         femmodel->InitTransientOutputx(&stackedinput[0],4);
                                         while(hydrotime<time-(yts*DBL_EPSILON)){ //loop on hydro dts
+                                                hydrostep+=1;
                                                 hydrotime+=hydrodt;
+                                                if(VerboseSolution()) _printf0_("hydro iteration " << hydrostep << "/" << hydroslices << "  time [yr]: " << setprecision(4) << hydrotime/yts << " (time step: " << hydrodt/yts << ")\n");
+
                                                 /*save preceding timestep*/
                                                 InputDuplicatex(femmodel,SedimentHeadHydrostepEnum,SedimentHeadOldEnum);