Changeset 16976

Timestamp:

11/29/13 10:34:32 (11 years ago)

Author:

seroussi

Message:

BUG: trying to fix SSAFS coupling

File:

• issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp (modified) (4 diffs)

issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

@@ -3215 +3215 @@
 
         /*compute all stiffness matrices for this element*/
+        printf("-------------- file: StressbalanceAnalysis.cpp line: %i\n",__LINE__); 
         ElementMatrix* Ke1=CreateKMatrixFS(element);
+        printf("-------------- file: StressbalanceAnalysis.cpp line: %i\n",__LINE__); 
         int indices[3]={18,19,20};
         Ke1->StaticCondensation(3,&indices[0]);
+        printf("-------------- file: StressbalanceAnalysis.cpp line: %i\n",__LINE__); 
         int init = element->FiniteElement();
         element->SetTemporaryElementType(P1Enum); 
+        printf("-------------- file: StressbalanceAnalysis.cpp line: %i\n",__LINE__); 
         ElementMatrix* Ke2=CreateKMatrixSSA3d(element);
+        printf("-------------- file: StressbalanceAnalysis.cpp line: %i\n",__LINE__); 
         element->SetTemporaryElementType(init); 
         ElementMatrix* Ke3=CreateKMatrixCouplingSSAFS(element);
+        printf("-------------- file: StressbalanceAnalysis.cpp line: %i\n",__LINE__); 
         ElementMatrix* Ke =new ElementMatrix(Ke1,Ke2,Ke3);
 
@@ -3536 +3542 @@
         }
 
-        ElementMatrix* Ke1=element->NewElementMatrix(SSAApproximationEnum);
+        ElementMatrix* Ke1=element->NewElementMatrixCoupling(6,SSAApproximationEnum);
         ElementMatrix* Ke2=element->NewElementMatrix(FSvelocityEnum);
         ElementMatrix* Ke=new ElementMatrix(Ke1,Ke2);
@@ -3642 +3648 @@
         Node **node_list = xNew<Node*>(2*vnumnodes-1+pnumnodes);
         for(i=0;i<vnumnodes-1;i++){
-                node_list[vnumnodes+pnumnodes+i] = pentabase->GetNode(i);
-                cs_list[vnumnodes+pnumnodes+i]   = XYEnum;
+                node_list[i] = pentabase->GetNode(i);
+                cs_list[i]   = XYEnum;
         }
         for(i=0;i<vnumnodes;i++){
-                node_list[i] = element->GetNode(i);
-                cs_list[i]   = XYZEnum;
+                node_list[i+vnumnodes-1] = element->GetNode(i);
+                cs_list[i+vnumnodes-1]   = XYZEnum;
         }
         for(i=0;i<pnumnodes;i++){
-                node_list[vnumnodes+i] = element->GetNode(vnumnodes+i);
-                cs_list[vnumnodes+i]   = PressureEnum;
+                node_list[2*vnumnodes-1+i] = element->GetNode(vnumnodes+i);
+                cs_list[2*vnumnodes-1+i]   = PressureEnum;
         }
 
         /*Initialize Element matrix and return if necessary*/
-        ElementMatrix* Ke1=element->NewElementMatrix(FSvelocityEnum);
-        ElementMatrix* Ke2=pentabase->NewElementMatrixCoupling(6,SSAApproximationEnum);
+        ElementMatrix* Ke1=pentabase->NewElementMatrixCoupling(6,SSAApproximationEnum);
+        ElementMatrix* Ke2=element->NewElementMatrix(FSvelocityEnum);
         ElementMatrix* Ke =new ElementMatrix(Ke1,Ke2);
         delete Ke1; delete Ke2;
@@ -3699 +3705 @@
 
         } 
-        for(i=0;i<numdofs;i++) for(j=0;j<numdofm;j++) Ke->values[i*numdoftotal+j+numdofs]+=Ke_gg[i][j];
-        for(i=0;i<numdofm;i++) for(j=0;j<numdofs;j++) Ke->values[(i+numdofs)*numdoftotal+j]+=Ke_gg2[i][j];
+        for(i=0;i<numdofs;i++) for(j=0;j<numdofm;j++) Ke->values[(i+2*numdofm)*numdoftotal+j]+=Ke_gg[i][j];
+        for(i=0;i<numdofm;i++) for(j=0;j<numdofs;j++) Ke->values[i*numdoftotal+(j+2*numdofm)]+=Ke_gg2[i][j];
 
         /*Transform Coordinate System*/