Changeset 3596

Timestamp:

04/21/10 16:11:25 (15 years ago)

Author:

Mathieu Morlighem

Message:

back to previous version for Eric

File:

• issm/trunk/src/c/objects/Tria.cpp (modified) (11 diffs)

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

@@ -110 +110 @@
 /*FUNCTION Tria::Configure {{{1*/
 void  Tria::Configure(void* ploadsin,void* pnodesin,void* pmaterialsin,void* pparametersin){
+
+        int i;
 
         DataSet* loadsin=NULL;
@@ -511 +513 @@
           /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
           GaussTria( &num_gauss, &first_gauss_area_coord, &second_gauss_area_coord, &third_gauss_area_coord, &gauss_weights, 2);
+
+#ifdef _ISSM_DEBUG_ 
+        for (i=0;i<num_gauss;i++){
+                printf("Gauss coord %i: %lf %lf %lf Weight: %lf\n",i,*(first_gauss_area_coord+i),*(second_gauss_area_coord+i),*(third_gauss_area_coord+i),*(gauss_weights+i));
+        }
+#endif
 
         /* Start  looping on the number of gaussian points: */
@@ -3798 +3806 @@
                 /* Get Jacobian determinant: */
                 GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
+#ifdef _ISSM_DEBUG_ 
+                printf("Element id %i Jacobian determinant: %g\n",GetId(),Jdet);
+#endif
 
                 /* Get nodal functions value at gaussian point:*/
@@ -3923 +3934 @@
         GetNodalFunctionsDerivatives(&dh1dh3[0][0],xyz_list, gauss_l1l2l3);
 
+        #ifdef _ISSM_DEBUG_ 
+        for (i=0;i<3;i++){
+                printf("Node %i  dh/dx=%lf dh/dy=%lf \n",i,dh1dh3[0][i],dh1dh3[1][i]);
+        }
+        #endif
+
         /*Build B: */
         for (i=0;i<numgrids;i++){
@@ -3957 +3974 @@
         /*Get dh1dh2dh3 in actual coordinate system: */
         GetNodalFunctions(&l1l2l3[0],gauss_l1l2l3);
+
+#ifdef _ISSM_DEBUG_ 
+        for (i=0;i<3;i++){
+                printf("Node %i  h=%lf \n",i,l1l2l3[i]);
+        }
+#endif
 
         /*Build B_prog: */
@@ -4535 +4558 @@
         /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
         GaussTria( &num_gauss, &first_gauss_area_coord, &second_gauss_area_coord, &third_gauss_area_coord, &gauss_weights, 4);
+#ifdef _ISSM_DEBUG_ 
+        for (i=0;i<num_gauss;i++){
+                printf("Gauss coord %i: %lf %lf %lf Weight: %lf\n",i,*(first_gauss_area_coord+i),*(second_gauss_area_coord+i),*(third_gauss_area_coord+i),*(gauss_weights+i));
+        }
+#endif
 
         /* Start  looping on the number of gaussian points: */
@@ -4564 +4592 @@
                 /* Get nodal functions value at gaussian point:*/
                 GetNodalFunctions(l1l2l3, gauss_l1l2l3);
+#ifdef _ISSM_DEBUG_
+                printf("viscositycomp %g thickness %g dvx [%g %g] dvy [%g %g]  dadjx [%g %g] dadjy[%g %g]\n",viscosity_complement,thickness,dvx[0],dvx[1],dvy[0],dvy[1],dadjx[0],dadjx[1],dadjy[0],dadjy[1]);
+#endif
 
                 /*Get nodal functions derivatives*/
@@ -4768 +4799 @@
                 GetParameterValue(&alpha_complement, &alpha_complement_list[0],gauss_l1l2l3);
                 GetParameterValue(&drag, &this->properties.k[0],gauss_l1l2l3);
+                #ifdef _ISSM_DEBUG_ 
+                        printf("Drag complement: %20.20lf Drag: %20.20lf\n",alpha_complement,drag);
+                #endif
 
                 /*recover lambda and mu: */
                 GetParameterValue(&lambda, &adjx_list[0],gauss_l1l2l3);
                 GetParameterValue(&mu, &adjy_list[0],gauss_l1l2l3);
+                #ifdef _ISSM_DEBUG_ 
+                        printf("Adjoint vector %20.20lf %20.20lf\n",lambda,mu);
+                #endif
                         
                 /*recover vx and vy: */
                 GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
+                #ifdef _ISSM_DEBUG_ 
+                        printf("Velocity vector %20.20lf %20.20lf\n",vx,vy);
+                #endif
 
                 /* Get Jacobian determinant: */
                 GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
-        
+                #ifdef _ISSM_DEBUG_ 
+                printf("Element id %i Jacobian determinant: %lf\n",GetId(),Jdet);
+                #endif
+                
                 /* Get nodal functions value at gaussian point:*/
                 GetNodalFunctions(l1l2l3, gauss_l1l2l3);
@@ -4990 +5033 @@
                 GetParameterValue(&alpha_complement, &alpha_complement_list[0],gauss_l1l2l3);
                 GetParameterValue(&drag, &this->properties.k[0],gauss_l1l2l3);
+#ifdef _ISSM_DEBUG_ 
+                printf("Drag complement: %20.20lf Drag: %20.20lf\n",alpha_complement,drag);
+#endif
 
                 /*recover lambda mu and xi: */
@@ -4995 +5041 @@
                 GetParameterValue(&mu, &adjy_list[0],gauss_l1l2l3);
                 GetParameterValue(&xi, &adjz_list[0],gauss_l1l2l3);
+#ifdef _ISSM_DEBUG_ 
+                printf("Adjoint vector %20.20lf %20.20lf\n",lambda,mu);
+#endif
 
                 /*recover vx vy and vz: */
@@ -5000 +5049 @@
                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
                 GetParameterValue(&vz, &vz_list[0],gauss_l1l2l3);
+#ifdef _ISSM_DEBUG_ 
+                printf("Velocity vector %20.20lf %20.20lf\n",vx,vy);
+
+                /*Get normal vecyor to the bed */
+                SurfaceNormal(&surface_normal[0],xyz_list);
+
+                bed_normal[0]=-surface_normal[0]; //Program is for surface, so the normal to the bed is the opposite of the result
+                bed_normal[1]=-surface_normal[1];
+                bed_normal[2]=-surface_normal[2];
+#endif
 
                 /* Get Jacobian determinant: */
                 GetJacobianDeterminant3d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
+#ifdef _ISSM_DEBUG_ 
+                printf("Element id %i Jacobian determinant: %lf\n",GetId(),Jdet);
+#endif
 
                 /* Get nodal functions value at gaussian point:*/