Changeset 1717

Timestamp:

08/17/09 10:53:04 (16 years ago)

Author:

Eric.Larour

Message:

transfervel: transfer results from md.results to md.vx, md.vy
Riftfront, Icefront, TriMeshProcessRifts: new enums for water, ice and air.

Location:

issm/trunk/src

Files:

• c/objects/Icefront.cpp (modified) (7 diffs)
• c/objects/Riftfront.cpp (modified) (4 diffs)
• c/objects/Riftfront.h (modified) (3 diffs)
• m/classes/public/presolve.m (modified) (1 diff)
• m/classes/public/transfervel.m (added)
• mex/TriMeshProcessRifts/TriMeshProcessRifts.cpp (modified) (1 diff)

issm/trunk/src/c/objects/Icefront.cpp

@@ -366 +366 @@
         /*Recover material and fill parameters: */
         matpar=(Matpar*)element->GetMatPar();
-        fill=element->GetShelf();
+        if (element->GetShelf())fill=WaterEnum();
+        else fill=AirEnum();
 
         //check that the element is onbed (collapsed formulation) otherwise:pe=0
@@ -493 +494 @@
         /*Recover material and fill parameters: */
         matpar=(Matpar*)element->GetMatPar();
-        fill=element->GetShelf();
+        if (element->GetShelf())fill=WaterEnum();
+        else fill=AirEnum();
 
         /* Set pe_g to 0: */
@@ -655 +657 @@
         /*Recover material and fill parameters: */
         matpar=(Matpar*)element->GetMatPar();
-        fill=element->GetShelf();
+        if (element->GetShelf())fill=WaterEnum();
+        else fill=AirEnum();
 
         /* Set pe_g to 0: */
@@ -865 +868 @@
                 bed=b1*(1+segment_gauss_coord[ig])/2+b2*(1-segment_gauss_coord[ig])/2;
 
-                if (fill==1){
+                if (fill==WaterEnum()){
                         //icefront ends in water: 
                         ice_pressure=1.0/2.0*gravity*rho_ice*pow(thickness,2);
@@ -875 +878 @@
                         water_pressure=1.0/2.0*gravity*rho_water*(pow(surface_under_water,2) - pow(base_under_water,2));
                 }
-                else if (fill==0){
+                else if (fill==AirEnum()){
                         ice_pressure=1.0/2.0*gravity*rho_ice*pow(thickness,2);
                         air_pressure=0;
@@ -1129 +1132 @@
 
                         //Now deal with water pressure: 
-                        if(fill==1){ //icefront ends in water
+                        if(fill==WaterEnum()){ //icefront ends in water
                                 water_level_above_g_tria=min(0,z_g[i]);//0 if the gaussian point is above water level
                                 water_pressure_tria=rho_water*gravity*water_level_above_g_tria;
                         }
-                        else if(fill==0){
+                        else if(fill==AirEnum()){
                                 water_pressure_tria=0;
                         }
@@ -1378 +1381 @@
 
                         //Now deal with water pressure: 
-                        if(fill==1){ //icefront ends in water
+                        if(fill==WaterEnum()){ //icefront ends in water
                                 water_level_above_g_tria=min(0,z_g[i]);//0 if the gaussian point is above water level
                                 water_pressure_tria=rho_water*gravity*water_level_above_g_tria;
                         }
-                        else if(fill==0){
+                        else if(fill==AirEnum()){
                                 water_pressure_tria=0;
                         }

issm/trunk/src/c/objects/Riftfront.cpp

@@ -23 +23 @@
 }
 
-Riftfront::Riftfront(char riftfront_type[RIFTFRONTSTRING],int riftfront_id, int riftfront_node_ids[MAX_RIFTFRONT_GRIDS], int riftfront_mparid, double riftfront_h[MAX_RIFTFRONT_GRIDS],double riftfront_b[MAX_RIFTFRONT_GRIDS],double riftfront_s[MAX_RIFTFRONT_GRIDS],double riftfront_normal[2],double riftfront_length,int riftfront_fill,double riftfront_friction, double riftfront_penalty_offset, bool riftfront_penalty_lock, bool riftfront_active,int riftfront_counter,bool riftfront_prestable,bool riftfront_shelf){
+Riftfront::Riftfront(char riftfront_type[RIFTFRONTSTRING],int riftfront_id, int riftfront_node_ids[MAX_RIFTFRONT_GRIDS], int riftfront_mparid, double riftfront_h[MAX_RIFTFRONT_GRIDS],double riftfront_b[MAX_RIFTFRONT_GRIDS],double riftfront_s[MAX_RIFTFRONT_GRIDS],double riftfront_normal[2],double riftfront_length,int riftfront_fill,double riftfront_friction, double riftfront_penalty_offset, int riftfront_penalty_lock, bool riftfront_active,int riftfront_counter,bool riftfront_prestable,bool riftfront_shelf){
 
         int i;
@@ -445 +445 @@
 
                 /*Ok, now compute the pressure (in norm) that is being applied to the flanks, depending on the type of fill: */
-                if(fill==WATERFILL){
+                if(fill==WaterEnum()){
                         if(shelf){
                                 /*We are on an ice shelf, hydrostatic equilibrium is used to determine the pressure for water fill: */
@@ -455 +455 @@
                         }
                 }
-                else if(fill==AIRFILL){
+                else if(fill==AirEnum()){
                         pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2;   //icefront on an ice sheet, pressure imbalance ice vs air.
                 }
-                else if(fill==ICEFILL){ //icefront finding itself against another icefront (pressure imbalance is fully compensated, ice vs ice)
+                else if(fill==IceEnum()){ //icefront finding itself against another icefront (pressure imbalance is fully compensated, ice vs ice)
                         pressure=0;
                 }
@@ -590 +590 @@
         /*Activate or deactivate penalties: */
         if(penetration<0){
+                printf("riftfront %i is penetrating\n",this->GetId());
                 /*There is penetration, we need to active the penalty so this penetration will be NULL: */
                 activate=1;

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

@@ -14 +14 @@
 #define MAX_RIFTFRONT_GRIDS 2 //max number of grids on a rift flank, only 2 because 2d for now.
 #define RIFTFRONTSTRING 20 //max string length
-
-/*Types of rift filling: */
-#define WATERFILL 1 
-#define AIRFILL 2 
-#define ICEFILL 3
 
 class Element;
@@ -49 +44 @@
 
                 double      penalty_offset;
-                bool        penalty_lock;
+                int         penalty_lock;
 
                 /*computational: */
@@ -60 +55 @@
 
                 Riftfront();
-                Riftfront(char type[RIFTFRONTSTRING],int id, int node_ids[MAX_RIFTFRONT_GRIDS], int mparid, double h[MAX_RIFTFRONT_GRIDS],double b[MAX_RIFTFRONT_GRIDS],double s[MAX_RIFTFRONT_GRIDS],double normal[2],double length,int fill,double friction, double penalty_offset, bool penalty_lock,bool active,int counter,bool prestable,bool shelf);
+                Riftfront(char type[RIFTFRONTSTRING],int id, int node_ids[MAX_RIFTFRONT_GRIDS], int mparid, double h[MAX_RIFTFRONT_GRIDS],double b[MAX_RIFTFRONT_GRIDS],double s[MAX_RIFTFRONT_GRIDS],double normal[2],double length,int fill,double friction, double penalty_offset, int penalty_lock,bool active,int counter,bool prestable,bool shelf);
                 ~Riftfront();
 

issm/trunk/src/m/classes/public/presolve.m

@@ -6 +6 @@
 
 %first, deal with rifts.
-if isnan(md.rifts),
+if isempty(md.rifts) | isnans(md.rifts),
         md.numrifts=0;
 else

issm/trunk/src/mex/TriMeshProcessRifts/TriMeshProcessRifts.cpp

@@ -285 +285 @@
                         /*Friction and fill: set to 0 both */
                         mxSetField(pmxa_array,i,"friction",mxCreateDoubleScalar(0));
-                        mxSetField(pmxa_array,i,"fill",mxCreateDoubleScalar(1)); //default is water
+                        mxSetField(pmxa_array,i,"fill",mxCreateDoubleScalar(WaterEnum())); //default is water