Changeset 8414

Timestamp:

05/24/11 12:07:28 (14 years ago)

Author:

seroussi

Message:

moved PVector basevert into Penta

Location:

issm/trunk/src/c/objects/Elements

Files:

• Penta.cpp (modified) (1 diff)
• Tria.cpp (modified) (1 diff)
• Tria.h (modified) (1 diff)

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

@@ -2827 +2827 @@
 ElementVector* Penta::CreatePVectorDiagnosticVertBase(void){
 
+
+        /*Constants*/
+        const int    numdof=NDOF1*NUMVERTICES;
+
+        /*Intermediaries */
+        int        i,j,ig;
+        int        approximation;
+        double     xyz_list[NUMVERTICES][3];
+        double     xyz_list_tria[NUMVERTICES2D][3];
+        double     Jdet2d;
+        double     vx,vy,vz,dbdx,dbdy,basalmeltingvalue;
+        double     slope[2];
+        double     l1l6[NUMVERTICES];
+        GaussPenta* gauss=NULL;
+
         if (!IsOnBed()) return NULL;
 
-        /*Call Tria function*/
-        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
-        ElementVector* pe=tria->CreatePVectorDiagnosticBaseVert();
-        delete tria->matice; delete tria;
+        /*Initialize Element vector*/
+        ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters);
+
+        /*Retrieve all inputs and parameters*/
+        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
+        for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) xyz_list_tria[i][j]=xyz_list[i][j];
+        inputs->GetParameterValue(&approximation,ApproximationEnum);
+        Input* bed_input=inputs->GetInput(BedEnum);             _assert_(bed_input);
+        Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input);
+        Input* vx_input=inputs->GetInput(VxEnum);               _assert_(vx_input);
+        Input* vy_input=inputs->GetInput(VyEnum);               _assert_(vy_input);
+        Input* vzstokes_input=NULL;
+        if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
+                vzstokes_input=inputs->GetInput(VzStokesEnum);       _assert_(vzstokes_input);
+        }
+
+        /* Start  looping on the number of gaussian points: */
+        gauss=new GaussPenta(0,1,2,2);
+        for(ig=gauss->begin();ig<gauss->end();ig++){
+
+                gauss->GaussPoint(ig);
+
+                basal_melting_input->GetParameterValue(&basalmeltingvalue, gauss);
+                bed_input->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],gauss);
+                vx_input->GetParameterValue(&vx, gauss);
+                vy_input->GetParameterValue(&vy, gauss);
+                if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
+                        vzstokes_input->GetParameterValue(&vz, gauss);
+                }
+                else vz=0;
+
+                dbdx=slope[0];
+                dbdy=slope[1];
+
+                GetTriaJacobianDeterminand(&Jdet2d, &xyz_list_tria[0][0],gauss);
+                GetNodalFunctionsP1(&l1l6[0], gauss);
+
+                for(i=0;i<numdof;i++) pe->values[i]+=-Jdet2d*gauss->weight*(vx*dbdx+vy*dbdy-vz-basalmeltingvalue)*l1l6[i];
+        }
 
         /*Clean up and return*/
+        delete gauss;
         return pe;
 }

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

@@ -1717 +1717 @@
 
                 for(i=0;i<numdof;i++) pe->values[i]+=Jdettria*gauss->weight*(surface_mass_balance_g-basal_melting_g)*L[i];
-        }
-
-        /*Clean up and return*/
-        delete gauss;
-        return pe;
-}
-/*}}}*/
-/*FUNCTION Tria::CreatePVectorDiagnosticBaseVert {{{1*/
-ElementVector* Tria::CreatePVectorDiagnosticBaseVert(void){
-
-        /*Constants*/
-        const int    numdof=NDOF1*NUMVERTICES;
-
-        /*Intermediaries */
-        int        i,j,ig;
-        int        approximation;
-        double     xyz_list[NUMVERTICES][3];
-        double     Jdet;
-        double     vx,vy,vz,dbdx,dbdy,basalmeltingvalue;
-        double     slope[2];
-        double     L[NUMVERTICES];
-        GaussTria* gauss=NULL;
-
-        /*Initialize Element vector*/
-        ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters);
-
-        /*Retrieve all inputs and parameters*/
-        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
-        inputs->GetParameterValue(&approximation,ApproximationEnum);
-        Input* bed_input=inputs->GetInput(BedEnum);             _assert_(bed_input);
-        Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input);
-        Input* vx_input=inputs->GetInput(VxEnum);               _assert_(vx_input);
-        Input* vy_input=inputs->GetInput(VyEnum);               _assert_(vy_input);
-        Input* vzstokes_input=NULL;
-        if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
-                vzstokes_input=inputs->GetInput(VzStokesEnum);       _assert_(vzstokes_input);
-        }
-
-        /* Start  looping on the number of gaussian points: */
-        gauss=new GaussTria(2);
-        for(ig=gauss->begin();ig<gauss->end();ig++){
-
-                gauss->GaussPoint(ig);
-
-                basal_melting_input->GetParameterValue(&basalmeltingvalue, gauss);
-                bed_input->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],gauss);
-                vx_input->GetParameterValue(&vx, gauss);
-                vy_input->GetParameterValue(&vy, gauss);
-                if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
-                        vzstokes_input->GetParameterValue(&vz, gauss);
-                }
-                else vz=0;
-
-                dbdx=slope[0];
-                dbdy=slope[1];
-
-                GetJacobianDeterminant3d(&Jdet, &xyz_list[0][0],gauss);
-                GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);
-
-                for(i=0;i<numdof;i++) pe->values[i]+=-Jdet*gauss->weight*(vx*dbdx+vy*dbdy-vz-basalmeltingvalue)*L[i];
         }
 

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

@@ -164 +164 @@
                 ElementVector* CreatePVectorBalancethickness_CG(void);
                 ElementVector* CreatePVectorBalancevelocities(void);
-                ElementVector* CreatePVectorDiagnosticBaseVert(void);
                 ElementVector* CreatePVectorDiagnosticMacAyeal(void);
                 ElementVector* CreatePVectorAdjointHoriz(void);