Changeset 10395

Timestamp:

11/01/11 07:49:56 (13 years ago)

Author:

Mathieu Morlighem

Message:

Now Stokes rotation is 1.e-10

Location:

issm/trunk/src/c/objects/Loads

Files:

• Pengrid.cpp (modified) (3 diffs)
• Pengrid.h (modified) (1 diff)

issm/trunk/src/c/objects/Loads/Pengrid.cpp

@@ -559 +559 @@
         Ke->values[2*NDOF4+2]= kmax*pow((double)10,penalty_offset);
 
+        /*Transform Coordinate System*/
+        TransformStiffnessMatrixCoord(Ke,NDOF4);
+
         /*Clean up and return*/
         return Ke;
@@ -682 +685 @@
 }
 /*}}}1*/
-/*FUNCTION Pengrid::UpdateInputs {{{1*/
-void  Pengrid::UpdateInputs(double* solution){
-        _error_("not supported yet!");
-}
-/*}}}1*/
 /*FUNCTION Pengrid::ResetConstraint {{{1*/
 void  Pengrid::ResetConstraint(void){
@@ -693 +691 @@
 }
 /*}}}1*/
+/*FUNCTION Pengrid::TransformStiffnessMatrixCoord{{{1*/
+void Pengrid::TransformStiffnessMatrixCoord(ElementMatrix* Ke,int dim){
+
+        int     i,j;
+        int     numnodes          = 1;
+        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,&node,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 Pengrid::UpdateInputs {{{1*/
+void  Pengrid::UpdateInputs(double* solution){
+        _error_("not supported yet!");
+}
+/*}}}1*/

issm/trunk/src/c/objects/Loads/Pengrid.h

@@ -90 +90 @@
                 void  UpdateInputs(double* solution);
                 void  ResetConstraint(void);
+                void  TransformStiffnessMatrixCoord(ElementMatrix* Ke,int dim);
                 /*}}}*/