Changeset 10375

Timestamp:

10/31/11 11:48:12 (13 years ago)

Author:

Mathieu Morlighem

Message:

extended local CS to MacAyeal3d

Location:

issm/trunk/src/c/objects/Elements

Files:

• Penta.cpp (modified) (4 diffs)
• Penta.h (modified) (2 diffs)
• Tria.h (modified) (1 diff)

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

@@ -4806 +4806 @@
 }
 /*}}}*/
-
 #endif
 
@@ -7076 +7075 @@
         int     i;
         double  rho_ice,g;
+        double  values0[numdof];
         double  values[numdof];
         double  vx[NUMVERTICES];
@@ -7091 +7091 @@
 
         /*Use the dof list to index into the solution vector: */
-        for(i=0;i<numdof;i++) values[i]=solution[doflist[i]];
+        for(i=0;i<numdof;i++) values0[i]=solution[doflist[i]];
+
+        /*Transform solution in Cartesian Space*/
+        TransformSolutionCoord(&values[0],&values0[0],2,true);
 
         /*Ok, we have vx and vy in values, fill in vx and vy arrays and extrude */
@@ -7709 +7712 @@
 }
 /*}}}*/
+/*FUNCTION Penta::TransformStiffnessMatrixCoord{{{1*/
+void Penta::TransformStiffnessMatrixCoord(ElementMatrix* Ke,int dim){
+
+        int     i,j;
+        int     numnodes          = NUMVERTICES;
+        double *transform         = NULL;
+        double *values            = NULL;
+
+        /*Copy current stiffness matrix*/
+        values=(double*)xmalloc(Ke->nrows*Ke->ncols*sizeof(double));
+        for(i=0;i<Ke->nrows;i++) for(j=0;j<Ke->ncols;j++) values[i*Ke->ncols+j]=Ke->values[i*Ke->ncols+j];
+
+        /*Get Coordinate Systems transform matrix*/
+        CoordinateSystemTransform(&transform,nodes,numnodes,dim);
+
+        /*Transform matrix: T*Ke*T^t */
+        TripleMultiply(transform,numnodes*dim,numnodes*dim,1,
+                                values,Ke->nrows,Ke->ncols,0,
+                                transform,numnodes*dim,numnodes*dim,0,
+                                &Ke->values[0],0);
+
+        /*Free Matrix*/
+        xfree((void**)&transform);
+        xfree((void**)&values);
+}
+/*}}}*/
+/*FUNCTION Penta::TransformLoadVectorCoord{{{1*/
+void Penta::TransformLoadVectorCoord(ElementVector* pe,int dim){
+
+        int     i,j;
+        int     numnodes          = NUMVERTICES;
+        double *transform         = NULL;
+        double *values            = NULL;
+
+        /*Copy current load vector*/
+        values=(double*)xmalloc(pe->nrows*sizeof(double));
+        for(i=0;i<pe->nrows;i++) values[i]=pe->values[i];
+
+        /*Get Coordinate Systems transform matrix*/
+        CoordinateSystemTransform(&transform,nodes,numnodes,dim);
+
+        /*Transform matrix: T^t*pe */
+        MatrixMultiply(transform,numnodes*dim,numnodes*dim,1,
+                                values,pe->nrows,1,0,
+                                &pe->values[0],0);
+
+        /*Free Matrix*/
+        xfree((void**)&transform);
+        xfree((void**)&values);
+}
+/*}}}*/
+/*FUNCTION Penta::TransformSolutionCoord{{{1*/
+void Penta::TransformSolutionCoord(double* solution,double* solution0,int dim,bool is2d){
+
+        int     i,j;
+        int     numnodes;
+        double *transform         = NULL;
+
+        /*Get Coordinate Systems transform matrix*/
+        if(is2d) numnodes=NUMVERTICES;
+        else     numnodes=NUMVERTICES2D;
+        CoordinateSystemTransform(&transform,nodes,numnodes,dim);
+
+        /*Transform matrix: T*U */
+        MatrixMultiply(transform,numnodes*dim,numnodes*dim,0,
+                                solution0,numnodes*dim,1,0,
+                                &solution[0],0);
+
+        /*Free Matrix*/
+        xfree((void**)&transform);
+}
+/*}}}*/
 #endif
 

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

@@ -235 +235 @@
                 void           InputUpdateFromSolutionDiagnosticVert( double* solutiong);
                 void           InputUpdateFromSolutionDiagnosticStokes( double* solutiong);
-                void           GetSolutionFromInputsDiagnosticHoriz(Vec solutiong);
-                void           GetSolutionFromInputsDiagnosticHutter(Vec solutiong);
-                void           GetSolutionFromInputsDiagnosticStokes(Vec solutiong);
-                void           GetSolutionFromInputsDiagnosticVert(Vec solutiong);
+                void             GetSolutionFromInputsDiagnosticHoriz(Vec solutiong);
+                void             GetSolutionFromInputsDiagnosticHutter(Vec solutiong);
+                void             GetSolutionFromInputsDiagnosticStokes(Vec solutiong);
+                void             GetSolutionFromInputsDiagnosticVert(Vec solutiong);
                 ElementVector* CreatePVectorCouplingMacAyealStokes(void);
                 ElementVector* CreatePVectorCouplingMacAyealStokesViscous(void);
@@ -258 +258 @@
                 ElementVector* CreatePVectorDiagnosticVertVolume(void);
                 ElementVector* CreatePVectorDiagnosticVertBase(void);
+                void           TransformStiffnessMatrixCoord(ElementMatrix* Ke,int dim);
+                void           TransformLoadVectorCoord(ElementVector* pe,int dim);
+                void           TransformSolutionCoord(double* solution,double* solution0,int dim,bool is2d=false);
                 #endif
 

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

@@ -192 +192 @@
                 void             SetClone(int* minranks);
                 void             SurfaceNormal(double* surface_normal, double xyz_list[3][3]);
-
                 
                 #ifdef _HAVE_DIAGNOSTIC_