Changeset 15641

Timestamp:

07/26/13 10:51:28 (12 years ago)

Author:

Mathieu Morlighem

Message:

CHG: minor cosmetics

Location:

issm/trunk-jpl/src/c/classes/Elements

Files:

• Penta.cpp (modified) (4 diffs)
• Tria.cpp (modified) (3 diffs)

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

@@ -8083 +8083 @@
         if(!isfront) return NULL;
 
-        /*Fetch number of nodes and dof for this finite element*/
-        int         numnodes = this->NumberofNodes();
-        int         numdof   = numnodes*NDOF2;
+        /*Intermediaries*/
         IssmDouble  rho_ice,rho_water,gravity;
         IssmDouble  Jdet,surface,z_g,water_pressure,ice_pressure,air_pressure;
         IssmDouble  surface_under_water,base_under_water,pressure;
-        GaussPenta*  gauss;
-        IssmDouble* basis = xNew<IssmDouble>(numnodes);
         IssmDouble  xyz_list_front[4][3];
         IssmDouble  area_coordinates[4][3];
         IssmDouble  normal[3];
-
+        GaussPenta* gauss;
+
+        /*Fetch number of nodes and dof for this finite element*/
+        int numnodes = this->NumberofNodes();
+        int numdof   = numnodes*NDOF2;
+
+        /*Initialize Element vector and other vectors*/
         GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
-        ElementVector* pe=new ElementVector(nodes,numnodes,this->parameters,HOApproximationEnum);
+        ElementVector* pe    = new ElementVector(nodes,numnodes,this->parameters,HOApproximationEnum);
+        IssmDouble*    basis = xNew<IssmDouble>(numnodes);
+
+        /*Retrieve all inputs and parameters*/
         Input* surface_input=inputs->GetInput(SurfaceEnum); _assert_(surface_input);
         rho_water=matpar->GetRhoWater();
@@ -8105 +8110 @@
         GetQuadNormal(&normal[0],xyz_list_front);
 
-        /*Start looping on Gaussian points*/
+        /*Initialize gauss points*/
         IssmDouble zmax=xyz_list[0][2]; for(int i=1;i<6;i++) if(xyz_list[i][2]>zmax) zmax=xyz_list[i][2];
         IssmDouble zmin=xyz_list[0][2]; for(int i=1;i<6;i++) if(xyz_list[i][2]<zmin) zmin=xyz_list[i][2];
@@ -8111 +8116 @@
         else                 gauss=new GaussPenta(area_coordinates,3,3);
 
+        /* Start  looping on the number of gaussian points: */
         for(int ig=gauss->begin();ig<gauss->end();ig++){
 
@@ -8119 +8125 @@
                 GetQuadJacobianDeterminant(&Jdet,xyz_list_front,gauss);
 
-                water_pressure=rho_water*gravity*min(0.,z_g);//0 if the gaussian point is above water level
-                ice_pressure=rho_ice*gravity*(surface-z_g);
-                air_pressure=0;
-                pressure = ice_pressure + water_pressure + air_pressure;
+                water_pressure = rho_water*gravity*min(0.,z_g);//0 if the gaussian point is above water level
+                ice_pressure   = rho_ice*gravity*(surface-z_g);
+                air_pressure   = 0;
+                pressure       = ice_pressure + water_pressure + air_pressure;
 
                 for (int i=0;i<numnodes;i++){

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

@@ -3158 +3158 @@
         if(!isfront) return NULL;
 
-        /*Fetch number of nodes and dof for this finite element*/
-        int         numnodes = this->NumberofNodes();
-        int         numdof   = numnodes*NDOF2;
+        /*Intermediaries*/
         IssmDouble  rho_ice,rho_water,gravity;
         IssmDouble  Jdet,thickness,bed,water_pressure,ice_pressure,air_pressure;
         IssmDouble  surface_under_water,base_under_water,pressure;
-        GaussTria*  gauss;
-        IssmDouble* basis = xNew<IssmDouble>(numnodes);
         IssmDouble  xyz_list_front[2][3];
         IssmDouble  area_coordinates[2][3];
         IssmDouble  normal[2];
 
+        /*Fetch number of nodes and dof for this finite element*/
+        int numnodes = this->NumberofNodes();
+        int numdof   = numnodes*NDOF2;
+
+        /*Initialize Element vector and other vectors*/
         GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
-        ElementVector* pe=new ElementVector(nodes,numnodes,this->parameters,SSAApproximationEnum);
+        ElementVector* pe    = new ElementVector(nodes,numnodes,this->parameters,HOApproximationEnum);
+        IssmDouble*    basis = xNew<IssmDouble>(numnodes);
+
+        /*Retrieve all inputs and parameters*/
         Input* thickness_input=inputs->GetInput(ThicknessEnum); _assert_(thickness_input);
         Input* bed_input      =inputs->GetInput(BedEnum);       _assert_(bed_input);
-        rho_water=matpar->GetRhoWater();
-        rho_ice  =matpar->GetRhoIce();
-        gravity  =matpar->GetG();
+        rho_water = matpar->GetRhoWater();
+        rho_ice   = matpar->GetRhoIce();
+        gravity   = matpar->GetG();
         GetZeroLevelsetCoordinates(&xyz_list_front[0][0],xyz_list,MaskIcelevelsetEnum);
         GetAreaCoordinates(&area_coordinates[0][0],xyz_list_front,xyz_list,2);
@@ -3182 +3186 @@
 
         /*Start looping on Gaussian points*/
-        gauss=new GaussTria(area_coordinates,3);
-
+        GaussTria* gauss=new GaussTria(area_coordinates,3);
         for(int ig=gauss->begin();ig<gauss->end();ig++){
 
@@ -3192 +3195 @@
                 GetNodalFunctions(basis,gauss);
 
-                surface_under_water=min(0.,thickness+bed); // 0 if the top of the glacier is above water level
-                base_under_water=min(0.,bed);              // 0 if the bottom of the glacier is above water level
-                water_pressure=1.0/2.0*gravity*rho_water*(pow(surface_under_water,2) - pow(base_under_water,2));
-                ice_pressure=1.0/2.0*gravity*rho_ice*pow(thickness,2);
-                air_pressure=0;
+                surface_under_water = min(0.,thickness+bed); // 0 if the top of the glacier is above water level
+                base_under_water    = min(0.,bed);           // 0 if the bottom of the glacier is above water level
+                water_pressure = 1.0/2.0*gravity*rho_water*(pow(surface_under_water,2) - pow(base_under_water,2));
+                ice_pressure   = 1.0/2.0*gravity*rho_ice*pow(thickness,2);
+                air_pressure   = 0;
                 pressure = ice_pressure + water_pressure + air_pressure;