Changeset 10400

Timestamp:

11/01/11 11:12:50 (13 years ago)

Author:

seroussi

Message:

bad grounding line 3d, to be improved

Location:

issm/trunk/src/c

Files:

• modules/GroundinglineMigrationx/GroundinglineMigrationx.cpp (modified) (6 diffs)
• modules/GroundinglineMigrationx/GroundinglineMigrationx.h (modified) (1 diff)
• modules/OutputResultsx/OutputResultsx.cpp (modified) (1 diff)
• objects/Elements/Penta.cpp (modified) (2 diffs)
• objects/Elements/Tria.cpp (modified) (1 diff)
• objects/Patch.h (modified) (2 diffs)
• solutions/transient_core.cpp (modified) (1 diff)

issm/trunk/src/c/modules/GroundinglineMigrationx/GroundinglineMigrationx.cpp

@@ -34 +34 @@
 
                 /*propagate ice shelf into connex areas of the ice sheet that potentially want to unground: */
-                vertices_ungrounding=PropagateFloatingiceToGrounded(elements,nodes,vertices,parameters,vertices_potentially_ungrounding);
+                vertices_ungrounding=PropagateFloatingiceToGroundedNeighbors(elements,nodes,vertices,parameters,vertices_potentially_ungrounding);
         }
 
@@ -110 +110 @@
 }
 /*}}}*/
-/*FUNCTION PropagateFloatingiceToGrounded {{{1*/
-double*    PropagateFloatingiceToGrounded(Elements* elements,Nodes* nodes,Vertices* vertices,Parameters* parameters,double* vertices_potentially_ungrounding){ 
+/*FUNCTION PropagateFloatingiceToGroundedNeighbors {{{1*/
+double*    PropagateFloatingiceToGroundedNeighbors(Elements* elements,Nodes* nodes,Vertices* vertices,Parameters* parameters,double* vertices_potentially_ungrounding){ 
 
         int      i,analysis_type;
@@ -117 +117 @@
         int      nflipped,local_nflipped;
         double*  nodes_on_floatingice                  = NULL;
-        double*  elements_neighbouring_floatingce      = NULL;
-        Vec      vec_elements_neighbouring_floatingice = NULL;
+        double*  elements_neighboring_floatingce      = NULL;
+        Vec      vec_elements_neighboring_floatingice = NULL;
         Vec      vec_nodes_on_floatingice              = NULL;
         Node*    node                                  = NULL;
@@ -136 +136 @@
                 
                 /*Vector of size number of elements*/
-                vec_elements_neighbouring_floatingice=NewVec(elements->NumberOfElements(),true);
+                vec_elements_neighboring_floatingice=NewVec(elements->NumberOfElements(),true);
 
                 /*Figure out if any of the nodes of the element will be floating -> elements neighbouting the floating ice*/
                 for(i=0;i<elements->Size();i++){
                         element=(Element*)elements->GetObjectByOffset(i);
-                        VecSetValue(vec_elements_neighbouring_floatingice,element->Sid(),element->IsNodeOnShelfFromFlags(nodes_on_floatingice)?1.0:0.0,INSERT_VALUES);
+                        VecSetValue(vec_elements_neighboring_floatingice,element->Sid(),element->IsNodeOnShelfFromFlags(nodes_on_floatingice)?1.0:0.0,INSERT_VALUES);
                 }
 
                 /*Assemble vector and serialize: */
-                VecAssemblyBegin(vec_elements_neighbouring_floatingice);
-                VecAssemblyEnd(vec_elements_neighbouring_floatingice);
-                VecToMPISerial(&elements_neighbouring_floatingce,vec_elements_neighbouring_floatingice);
+                VecAssemblyBegin(vec_elements_neighboring_floatingice);
+                VecAssemblyEnd(vec_elements_neighboring_floatingice);
+                VecToMPISerial(&elements_neighboring_floatingce,vec_elements_neighboring_floatingice);
 
-                /*Go through elements_neighbouring_floatingce, and update vector of the nodes that will start floating*/
+                /*Go through elements_neighboring_floatingce, and update vector of the nodes that will start floating*/
                 local_nflipped=0;
                 for(i=0;i<elements->Size();i++){
                         element=(Element*)elements->GetObjectByOffset(i);
-                        if(elements_neighbouring_floatingce[element->Sid()]){
+                        if(elements_neighboring_floatingce[element->Sid()]){
                                 local_nflipped+=element->UpdatePotentialSheetUngrounding(vertices_potentially_ungrounding,vec_nodes_on_floatingice,nodes_on_floatingice);
                         }
@@ -164 +164 @@
 
                 /*Avoid leaks: */
-                xfree((void**)&elements_neighbouring_floatingce);
+                xfree((void**)&elements_neighboring_floatingce);
                 xfree((void**)&nodes_on_floatingice);
 
                 /*Assemble and serialize:*/
-                VecFree(&vec_elements_neighbouring_floatingice);
+                VecFree(&vec_elements_neighboring_floatingice);
                 VecToMPISerial(&nodes_on_floatingice,vec_nodes_on_floatingice); 
         }
@@ -174 +174 @@
         /*Free ressources:*/
         VecFree(&vec_nodes_on_floatingice);
-        xfree((void**)&elements_neighbouring_floatingce);
+        xfree((void**)&elements_neighboring_floatingce);
 
         return nodes_on_floatingice;

issm/trunk/src/c/modules/GroundinglineMigrationx/GroundinglineMigrationx.h

@@ -16 +16 @@
 Vec        CreateNodesOnFloatingIce(Nodes* nodes,int configuration_type);
 double*    PotentialSheetUngrounding(Elements* elements,Vertices* vertices,Parameters* parameters);
-double*    PropagateFloatingiceToGrounded(Elements* elements,Nodes* nodes,Vertices* vertices,Parameters* parameters,double* vertices_potentially_ungrounding);
+double*    PropagateFloatingiceToGroundedNeighbors(Elements* elements,Nodes* nodes,Vertices* vertices,Parameters* parameters,double* vertices_potentially_ungrounding);
 #endif  /* _GROUNDINGLINEMIGRATIONX_H */

issm/trunk/src/c/modules/OutputResultsx/OutputResultsx.cpp

@@ -26 +26 @@
         char       *outputfilename          = NULL;
         char        cpu_outputfilename[100];        //easier to convert an integer with sprintf
-        bool        io_gather,results_on_vertices;
+        bool        io_gather;
         int         solutiontype;
         Results    *results                 = NULL;
-        bool        dakota_analysis            = false;
+        bool        dakota_analysis         = false;
         
         #ifdef _SERIAL_

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

@@ -2038 +2038 @@
 /*FUNCTION Penta::MigrateGroundingLine{{{1*/
 void  Penta::MigrateGroundingLine(double* old_floating_ice,double* sheet_ungrounding){
-        _error_("not supported yet!");
+
+        int     i,migration_style,unground;
+        bool    elementonshelf = false;
+        double  bed_hydro,yts,gl_melting_rate;
+        double  rho_water,rho_ice,density;
+        double  melting[NUMVERTICES];
+        double  h[NUMVERTICES],s[NUMVERTICES],b[NUMVERTICES],ba[NUMVERTICES];
+
+        /*Recover info at the vertices: */
+        parameters->FindParam(&migration_style,GroundinglineMigrationEnum);
+        parameters->FindParam(&yts,ConstantsYtsEnum);
+        GetInputListOnVertices(&h[0],ThicknessEnum);
+        GetInputListOnVertices(&s[0],SurfaceEnum);
+        GetInputListOnVertices(&b[0],BedEnum);
+        GetInputListOnVertices(&ba[0],BathymetryEnum);
+        rho_water=matpar->GetRhoWater();
+        rho_ice=matpar->GetRhoIce();
+        density=rho_ice/rho_water;
+        
+        /*go through vertices, and update inputs, considering them to be PentaVertex type: */
+        for(i=0;i<NUMVERTICES;i++){
+                /*Ice shelf: if bed below bathymetry, impose it at the bathymetry and update surface, elso do nothing */
+                if(old_floating_ice[nodes[i]->Sid()]){
+                        if(b[i]<=ba[i]){ 
+                                b[i]=ba[i];
+                                s[i]=b[i]+h[i];
+                                nodes[i]->inputs->AddInput(new BoolInput(MaskVertexonfloatingiceEnum,false));
+                                nodes[i]->inputs->AddInput(new BoolInput(MaskVertexongroundediceEnum,true));
+                        }
+                }
+                /*Ice sheet: if hydrostatic bed above bathymetry, ice sheet starts to unground, elso do nothing */
+                /*Change only if AgressiveMigration or if the ice sheet is in contact with the ocean*/
+                else{
+                        bed_hydro=-density*h[i];
+                        if (bed_hydro>ba[i]){
+                                /*Unground only if the element is connected to the ice shelf*/
+                                if(migration_style==AgressiveMigrationEnum){
+                                        s[i]=(1-density)*h[i];
+                                        b[i]=-density*h[i];
+                                        nodes[i]->inputs->AddInput(new BoolInput(MaskVertexonfloatingiceEnum,true));
+                                        nodes[i]->inputs->AddInput(new BoolInput(MaskVertexongroundediceEnum,false));
+                                }
+                                else if(migration_style==SoftMigrationEnum && sheet_ungrounding[nodes[i]->Sid()]){
+                                        s[i]=(1-density)*h[i];
+                                        b[i]=-density*h[i];
+                                        nodes[i]->inputs->AddInput(new BoolInput(MaskVertexonfloatingiceEnum,true));
+                                        nodes[i]->inputs->AddInput(new BoolInput(MaskVertexongroundediceEnum,false));
+                                }
+                        }
+                }
+        }
+
+        /*If at least one vertex is now floating, the element is now floating*/
+        for(i=0;i<NUMVERTICES;i++){
+                if(nodes[i]->IsFloating()){
+                        elementonshelf=true;
+                        break;
+                }
+        }
+        
+   /*Add basal melting rate if element just ungrounded*/
+        if(!this->IsFloating() && elementonshelf==true){
+                for(i=0;i<NUMVERTICES;i++)melting[i]=gl_melting_rate/yts;
+                this->inputs->AddInput(new PentaVertexInput(BasalforcingsMeltingRateEnum,&melting[0]));
+        } 
+
+        /*Update inputs*/
+   this->inputs->AddInput(new BoolInput(MaskElementonfloatingiceEnum,elementonshelf));
+
+        /*Update inputs*/    
+        this->inputs->AddInput(new PentaVertexInput(SurfaceEnum,&s[0]));
+        this->inputs->AddInput(new PentaVertexInput(BedEnum,&b[0]));
 }
 /*}}}*/
@@ -2100 +2171 @@
 
         int     i;
-        int     numrows  = 0;
-        int     numnodes = 0;
-        int     temp_numnodes=0;
+        int     numrows       = 0;
+        int     numnodes      = 0;
+        int     temp_numnodes = 0;
 
         /*Go through all the results objects, and update the counters: */

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

@@ -2035 +2035 @@
 
         int     i;
-        int     numrows     = 0;
-        int     numnodes    = 0;
-        int     temp_numnodes    = 0;
+        int     numrows       = 0;
+        int     numnodes      = 0;
+        int     temp_numnodes = 0;
 
         /*Go through all the results objects, and update the counters: */

issm/trunk/src/c/objects/Patch.h

@@ -25 +25 @@
 
 /*Headers:*/
-/*{{{1*/
-/*}}}*/
-
 class Patch{
 
@@ -47 +44 @@
 
 };
-
 #endif //ifndef _PATCH_H_

issm/trunk/src/c/solutions/transient_core.cpp

@@ -97 +97 @@
 
                 if(isgroundingline){
-                        if(dim==3) _error_("Grounding line migration not implemented in 3d");
+                //      if(dim==3) _error_("Grounding line migration not implemented in 3d");
                         _printf_(VerboseSolution(),"   computing new grounding line position\n");
                         #ifdef _HAVE_GROUNDINGLINE_