[cig-commits] r6878 - in short/3D/PyLith/trunk: libsrc/feassemble unittests/libtests/feassemble/data

brad at geodynamics.org brad at geodynamics.org
Mon May 14 17:51:40 PDT 2007


Author: brad
Date: 2007-05-14 17:51:39 -0700 (Mon, 14 May 2007)
New Revision: 6878

Modified:
   short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc
   short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc
   short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc
   short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc
   short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureApp.py
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DLinear.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DQuadratic.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DLinear.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DQuadratic.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DLinear.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DQuadratic.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXY.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXYZ.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXZ.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearYZ.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DQuadratic.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DLinear.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc
   short/3D/PyLith/trunk/unittests/libtests/feassemble/data/feutils.py
Log:
Fixed bug in quadrature. Jacobian stuff was transposed. A few unit tests now fail and need to be fixed.

Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -67,7 +67,8 @@
     } // for
     
     // Compute Jacobian at quadrature point
-    // J = [dx/dp dy/dp]
+    // J = [dx/dp
+    //      dy/dp]
     // dx/dp = sum[i=0,n-1] (dNi/dp * xi)
     // dy/dp = sum[i=0,n-1] (dNi/dp * yi)
     for (int iBasis=0; iBasis < _numBasis; ++iBasis) {
@@ -78,7 +79,7 @@
     } // for
 
     // Compute determinant of Jacobian at quadrature point
-    // |J| = sqrt(J transpose(J))
+    // |J| = sqrt(transpose(J) J)
     double det = 0.0;
     for (int iDim=0; iDim < _spaceDim; ++iDim)
       det += _jacobian[iQuadPt*_spaceDim+iDim] * 

Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -68,7 +68,9 @@
     } // for
     
     // Compute Jacobian at quadrature point
-    // J = [dx/dp dy/dp dz/dp]
+    // J = [dx/dp
+    //      dy/dp
+    //      dz/dp]
     // dx/dp = sum[i=0,n-1] (dNi/dp * xi)
     // dy/dp = sum[i=0,n-1] (dNi/dp * yi)
     // dz/dp = sum[i=0,n-1] (dNi/dp * zi)
@@ -80,7 +82,7 @@
     } // for
 
     // Compute determinant of Jacobian at quadrature point
-    // |J| = sqrt(J transpose(J))
+    // |J| = sqrt(transpose(J) J)
     double det = 0.0;
     for (int iDim=0; iDim < _spaceDim; ++iDim)
       det += _jacobian[iQuadPt*_spaceDim+iDim] * 

Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -67,19 +67,19 @@
     } // for
     
     // Compute Jacobian at quadrature point
-    // J = [dx/dp dy/dp]
-    //     [dx/dq dy/dq]
+    // J = [dx/dp dx/dq]
+    //     [dy/dp dy/dq]
     // dx/dp = sum[i=0,n-1] (dNi/dp * xi)
+    // dx/dq = sum[i=0,n-1] (dNi/dq * xi)
     // dy/dp = sum[i=0,n-1] (dNi/dp * yi)
-    // dx/dq = sum[i=0,n-1] (dNi/dq * xi)
     // dy/dq = sum[i=0,n-1] (dNi/dq * yi)
     for (int iBasis=0; iBasis < _numBasis; ++iBasis)
-      for (int iRow=0; iRow < _cellDim; ++iRow) {
+      for (int iCol=0; iCol < _cellDim; ++iCol) {
 	const double deriv = 
-	  _basisDeriv[iQuadPt*_numBasis*_spaceDim+iBasis*_cellDim+iRow];
-	for (int iCol=0; iCol < _spaceDim; ++iCol)
-	  _jacobian[iQuadPt*_cellDim*_spaceDim+iRow*_spaceDim+iCol] +=
-	    deriv * vertCoords[iBasis*_spaceDim+iCol];
+	  _basisDeriv[iQuadPt*_numBasis*_spaceDim+iBasis*_cellDim+iCol];
+	for (int iRow=0; iRow < _spaceDim; ++iRow)
+	  _jacobian[iQuadPt*_cellDim*_spaceDim+iRow*_cellDim+iCol] +=
+	    deriv * vertCoords[iBasis*_spaceDim+iRow];
       } // for
   
     // Compute determinant of Jacobian at quadrature point

Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -71,113 +71,114 @@
     } // for
     
     // Compute Jacobian at quadrature point
-    // J = [dx/dp dy/dp dz/dp]
-    //     [dx/dq dy/dq dz/dq]
+    // J = [dx/dp dx/dq]
+    //     [dy/dp dy/dq]
+    //     [dz/dp dz/dq]
     // dx/dp = sum[i=0,n-1] (dNi/dp * xi)
+    // dx/dq = sum[i=0,n-1] (dNi/dq * xi)
     // dy/dp = sum[i=0,n-1] (dNi/dp * yi)
+    // dy/dq = sum[i=0,n-1] (dNi/dq * yi)
     // dz/dp = sum[i=0,n-1] (dNi/dp * zi)
-    // dx/dq = sum[i=0,n-1] (dNi/dq * xi)
-    // dy/dq = sum[i=0,n-1] (dNi/dq * yi)
     // dz/dq = sum[i=0,n-1] (dNi/dq * zi)
     for (int iBasis=0; iBasis < _numBasis; ++iBasis)
-      for (int iRow=0; iRow < _cellDim; ++iRow) {
+      for (int iCol=0; iCol < _cellDim; ++iCol) {
 	const double deriv = 
-	  _basisDeriv[iQuadPt*_numBasis*_cellDim+iBasis*_cellDim+iRow];
-	for (int iCol=0; iCol < _spaceDim; ++iCol)
-	  _jacobian[iQuadPt*_cellDim*_spaceDim+iRow*_spaceDim+iCol] +=
-	    deriv * vertCoords[iBasis*+_spaceDim+iCol];
+	  _basisDeriv[iQuadPt*_numBasis*_cellDim+iBasis*_cellDim+iCol];
+	for (int iRow=0; iRow < _spaceDim; ++iRow)
+	  _jacobian[iQuadPt*_cellDim*_spaceDim+iRow*_cellDim+iCol] +=
+	    deriv * vertCoords[iBasis*+_spaceDim+iRow];
       } // for
     
     // Compute determinant of Jacobian at quadrature point
-    // |J| = sqrt(J transpose(J))
+    // |J| = sqrt(transpose(J) J)
     const int iJ = iQuadPt*_cellDim*_spaceDim;
-    const int i00 = iJ + 0*_spaceDim + 0;
-    const int i01 = iJ + 0*_spaceDim + 1;
-    const int i02 = iJ + 0*_spaceDim + 2;
-    const int i10 = iJ + 1*_spaceDim + 0;
-    const int i11 = iJ + 1*_spaceDim + 1;
-    const int i12 = iJ + 1*_spaceDim + 2;
-    // JJ = J transpose(J)
+    const int i00 = iJ + 0*_cellDim + 0;
+    const int i01 = iJ + 0*_cellDim + 1;
+    const int i10 = iJ + 1*_cellDim + 0;
+    const int i11 = iJ + 1*_cellDim + 1;
+    const int i20 = iJ + 2*_cellDim + 0;
+    const int i21 = iJ + 2*_cellDim + 1;
+    // JJ = transpose(J) J 
     const double jj00 = 
       _jacobian[i00]*_jacobian[i00] +
+      _jacobian[i10]*_jacobian[i10] +
+      _jacobian[i20]*_jacobian[i20];
+    const double jj10 =
+      _jacobian[i00]*_jacobian[i01] +
+      _jacobian[i10]*_jacobian[i11] +
+      _jacobian[i20]*_jacobian[i21];
+    const double jj01 = jj10;
+    const double jj11 = 
       _jacobian[i01]*_jacobian[i01] +
-      _jacobian[i02]*_jacobian[i02];
-    const double jj01 =
-      _jacobian[i00]*_jacobian[i10] +
-      _jacobian[i01]*_jacobian[i11] +
-      _jacobian[i02]*_jacobian[i12];
-    const double jj10 = jj01;
-    const double jj11 = 
-      _jacobian[i10]*_jacobian[i10] +
       _jacobian[i11]*_jacobian[i11] +
-      _jacobian[i12]*_jacobian[i12];
+      _jacobian[i21]*_jacobian[i21];
     const double det = sqrt(jj00*jj11 - jj01*jj10);
     _checkJacobianDet(det);
     _jacobianDet[iQuadPt] = det;
     
     // Compute inverse of Jacobian at quadrature point
-    const double d01 = 
+    const double d10 = 
       _jacobian[i00]*_jacobian[i11] - 
-      _jacobian[i10]*_jacobian[i01];
-    const double d12 = 
-      _jacobian[i01]*_jacobian[i12] - 
-      _jacobian[i11]*_jacobian[i02];
-    const double d02 = 
-      _jacobian[i00]*_jacobian[i12] - 
-      _jacobian[i10]*_jacobian[i02];
-    if (fabs(d01) > _minJacobian) {
-      // Jinv00 = 1/d01 * J11
-      // Jinv01 = 1/d01 * -J01
-      // Jinv10 = 1/d01 * -J10
-      // Jinv11 = 1/d01 * J00
-      _jacobianInv[iJ+0] =  _jacobian[i11] / d01; // Jinv00
-      _jacobianInv[iJ+1] = -_jacobian[i01] / d01; // Jinv01
-      _jacobianInv[iJ+2] = -_jacobian[i10] / d01; // Jinv10
-      _jacobianInv[iJ+3] =  _jacobian[i00] / d01; // Jinv11
-      if (fabs(d12) > _minJacobian) {
-	// Jinv20 = 1/d12 -J11
-	// Jinv21 = 1/d12 J01
-	_jacobianInv[iJ+4] = -_jacobian[i11] / d12; // Jinv20
-	_jacobianInv[iJ+5] =  _jacobian[i01] / d12; // Jinv21
+      _jacobian[i01]*_jacobian[i10];
+    const double d21 = 
+      _jacobian[i10]*_jacobian[i21] - 
+      _jacobian[i11]*_jacobian[i20];
+    const double d20 = 
+      _jacobian[i00]*_jacobian[i21] - 
+      _jacobian[i01]*_jacobian[i20];
+    if (fabs(d20) > _minJacobian) {
+      // Jinv00 = 1/d10 * J11
+      // Jinv10 = 1/d10 * -J10
+      // Jinv02 = 1/d10 * -J01
+      // Jinv11 = 1/d10 * J00
+      _jacobianInv[iJ+0] =  _jacobian[i11] / d10; // Jinv00
+      _jacobianInv[iJ+1] = -_jacobian[i10] / d10; // Jinv10
+      _jacobianInv[iJ+2] = -_jacobian[i01] / d10; // Jinv01
+      _jacobianInv[iJ+3] =  _jacobian[i00] / d10; // Jinv11
+      if (fabs(d21) > _minJacobian) {
+	// Jinv02 = 1/d21 -J11
+	// Jinv12 = 1/d21 J10
+	_jacobianInv[iJ+4] = -_jacobian[i11] / d21; // Jinv02
+	_jacobianInv[iJ+5] =  _jacobian[i10] / d21; // Jinv12
 	
-      } else if (fabs(d02) > _minJacobian) {
-	// Jinv20 = 1/d02 -J10
-	// Jinv21 = 1/d02 J00
-	_jacobianInv[iJ+4] = -_jacobian[i10] / d02; // Jinv20
-	_jacobianInv[iJ+5] =  _jacobian[i00] / d02; // Jinv21
+      } else if (fabs(d20) > _minJacobian) {
+	// Jinv02 = 1/d20 -J01
+	// Jinv12 = 1/d20 J00
+	_jacobianInv[iJ+4] = -_jacobian[i01] / d20; // Jinv02
+	_jacobianInv[iJ+5] =  _jacobian[i00] / d20; // Jinv12
       } else {
-	_jacobianInv[iJ+4] = 0.0; // Jinv20
-	_jacobianInv[iJ+5] = 0.0; // Jinv21
+	_jacobianInv[iJ+4] = 0.0; // Jinv02
+	_jacobianInv[iJ+5] = 0.0; // Jinv12
       } // if/else
-    } else if (fabs(d02) > _minJacobian) {
-      // Jinv00 = 1/d02 * J12
-      // Jinv01 = 1/d02 * -J02
-      // Jinv20 = 1/d02 * -J10
-      // Jinv21 = 1/d02 * J00
-      _jacobianInv[iJ+0] =  _jacobian[i12] / d02; // Jinv00
-      _jacobianInv[iJ+1] = -_jacobian[i02] / d02; // Jinv01
-      _jacobianInv[iJ+4] = -_jacobian[i10] / d02; // Jinv20
-      _jacobianInv[iJ+5] =  _jacobian[i00] / d02; // Jinv21
-      if (fabs(d12) > _minJacobian) {
-	// Jinv10 = 1/d12 J12
-	// Jinv11 = 1/d12 -J02
-	_jacobianInv[iJ+2] = -_jacobian[i12] / d12; // Jinv10
-	_jacobianInv[iJ+3] =  _jacobian[i02] / d12; // Jinv11
+    } else if (fabs(d20) > _minJacobian) {
+      // Jinv00 = 1/d20 * J21
+      // Jinv01 = 1/d20 * -J20
+      // Jinv20 = 1/d20 * -J01
+      // Jinv21 = 1/d20 * J00
+      _jacobianInv[iJ+0] =  _jacobian[i21] / d20; // Jinv00
+      _jacobianInv[iJ+1] = -_jacobian[i20] / d20; // Jinv10
+      _jacobianInv[iJ+4] = -_jacobian[i01] / d20; // Jinv02
+      _jacobianInv[iJ+5] =  _jacobian[i00] / d20; // Jinv12
+      if (fabs(d21) > _minJacobian) {
+	// Jinv10 = 1/d21 J21
+	// Jinv11 = 1/d21 -J20
+	_jacobianInv[iJ+2] = -_jacobian[i21] / d21; // Jinv01
+	_jacobianInv[iJ+3] =  _jacobian[i20] / d21; // Jinv11
       } else {
-	_jacobianInv[iJ+2] = 0.0; // Jinv10
+	_jacobianInv[iJ+2] = 0.0; // Jinv01
 	_jacobianInv[iJ+3] = 0.0; // Jinv11
       } // if/else
-    } else if (fabs(d12) > _minJacobian) {
+    } else if (fabs(d21) > _minJacobian) {
       _jacobianInv[iJ+0] = 0.0; // Jinv00
-      _jacobianInv[iJ+1] = 0.0; // Jinv01
-      // Jinv10 = 1/d12 J12
-      // Jinv11 = 1/d12 -J02
-      // Jinv20 = 1/d12 -J11
-      // Jinv21 = 1/d12 J01
-      _jacobianInv[iJ+2] =  _jacobian[i12] / d12; // Jinv10
-      _jacobianInv[iJ+3] = -_jacobian[i02] / d12; // Jin11
-      _jacobianInv[iJ+4] = -_jacobian[i11] / d12; // Jinv20
-      _jacobianInv[iJ+5] =  _jacobian[i01] / d12; // Jinv21
+      _jacobianInv[iJ+1] = 0.0; // Jinv10
+      // Jinv01 = 1/d21 J21
+      // Jinv11 = 1/d21 -J20
+      // Jinv02 = 1/d21 -J11
+      // Jinv12 = 1/d21 J10
+      _jacobianInv[iJ+2] =  _jacobian[i21] / d21; // Jinv01
+      _jacobianInv[iJ+3] = -_jacobian[i20] / d21; // Jin11
+      _jacobianInv[iJ+4] = -_jacobian[i11] / d21; // Jinv02
+      _jacobianInv[iJ+5] =  _jacobian[i10] / d21; // Jinv12
     } else
       throw std::runtime_error("Could not invert Jacobian.");
   } // for

Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -68,25 +68,25 @@
     } // for
     
     // Compute Jacobian at quadrature point
-    // J = [dx/dp dy/dp dz/dp]
-    //     [dx/dq dy/dq dz/dq]
-    //     [dx/dr dy/dr dz/dr]
+    // J = [dx/dp dx/dq dx/dr]
+    //     [dy/dp dy/dq dy/dr]
+    //     [dz/dp dz/dq dz/dr]
     // dx/dp = sum[i=0,n-1] (dNi/dp * xi)
+    // dx/dq = sum[i=0,n-1] (dNi/dq * xi)
+    // dx/dr = sum[i=0,n-1] (dNi/dr * xi)
     // dy/dp = sum[i=0,n-1] (dNi/dp * yi)
+    // dy/dq = sum[i=0,n-1] (dNi/dq * yi)
+    // dy/dr = sum[i=0,n-1] (dNi/dr * yi)
     // dz/dp = sum[i=0,n-1] (dNi/dp * zi)
-    // dx/dq = sum[i=0,n-1] (dNi/dq * xi)
-    // dy/dq = sum[i=0,n-1] (dNi/dq * yi)
     // dz/dq = sum[i=0,n-1] (dNi/dq * zi)
-    // dx/dr = sum[i=0,n-1] (dNi/dr * xi)
-    // dy/dr = sum[i=0,n-1] (dNi/dr * yi)
     // dz/dr = sum[i=0,n-1] (dNi/dr * zi)
     for (int iBasis=0; iBasis < _numBasis; ++iBasis)
-      for (int iRow=0; iRow < _cellDim; ++iRow) {
+      for (int iCol=0; iCol < _cellDim; ++iCol) {
 	const double deriv = 
-	  _basisDeriv[iQuadPt*_numBasis*_spaceDim+iBasis*_cellDim+iRow];
-	for (int iCol=0; iCol < _spaceDim; ++iCol)
-	  _jacobian[iQuadPt*_cellDim*_spaceDim+iRow*_spaceDim+iCol] += 
-	    deriv * vertCoords[iBasis*_spaceDim+iCol];
+	  _basisDeriv[iQuadPt*_numBasis*_spaceDim+iBasis*_cellDim+iCol];
+	for (int iRow=0; iRow < _spaceDim; ++iRow)
+	  _jacobian[iQuadPt*_cellDim*_spaceDim+iRow*_cellDim+iCol] += 
+	    deriv * vertCoords[iBasis*_spaceDim+iRow];
       } // for
 
     // Compute determinant of Jacobian at quadrature point

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureApp.py
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureApp.py	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureApp.py	2007-05-15 00:51:39 UTC (rev 6878)
@@ -163,13 +163,13 @@
                        format="%16.8e", ncols=self.cellDim)
     self.data.addArray(vtype="double", name="_jacobian",
                        values=self.jacobian,
-                       format="%16.8e", ncols=self.spaceDim)
+                       format="%16.8e", ncols=self.cellDim)
     self.data.addArray(vtype="double", name="_jacobianDet",
                        values=self.jacobianDet,
                        format="%16.8e", ncols=self.numQuadPts)
     self.data.addArray(vtype="double", name="_jacobianInv",
                        values=self.jacobianInv,
-                       format="%16.8e", ncols=self.cellDim)
+                       format="%16.8e", ncols=self.spaceDim)
       
     return
 

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DLinear.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DLinear.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DLinear.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -59,7 +59,8 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din2DLinear::_jacobian[] = {
-  4.50000000e-01,  4.00000000e-01,
+  4.50000000e-01,
+  4.00000000e-01,
 };
 
 const double pylith::feassemble::QuadratureData1Din2DLinear::_jacobianDet[] = {
@@ -67,8 +68,7 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din2DLinear::_jacobianInv[] = {
-  2.22222222e+00,
-  2.50000000e+00,
+  2.22222222e+00,  2.50000000e+00,
 };
 
 pylith::feassemble::QuadratureData1Din2DLinear::QuadratureData1Din2DLinear(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DQuadratic.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DQuadratic.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din2DQuadratic.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -70,8 +70,10 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din2DQuadratic::_jacobian[] = {
-  5.07735027e-01,  2.84529946e-01,
-  3.92264973e-01,  5.15470054e-01,
+  5.07735027e-01,
+  2.84529946e-01,
+  3.92264973e-01,
+  5.15470054e-01,
 };
 
 const double pylith::feassemble::QuadratureData1Din2DQuadratic::_jacobianDet[] = {
@@ -79,10 +81,8 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din2DQuadratic::_jacobianInv[] = {
-  1.96953125e+00,
-  3.51456855e+00,
-  2.54929721e+00,
-  1.93997691e+00,
+  1.96953125e+00,  3.51456855e+00,
+  2.54929721e+00,  1.93997691e+00,
 };
 
 pylith::feassemble::QuadratureData1Din2DQuadratic::QuadratureData1Din2DQuadratic(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DLinear.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DLinear.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DLinear.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -59,7 +59,9 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din3DLinear::_jacobian[] = {
- -7.50000000e-01,  1.75000000e+00,  2.50000000e+00,
+ -7.50000000e-01,
+  1.75000000e+00,
+  2.50000000e+00,
 };
 
 const double pylith::feassemble::QuadratureData1Din3DLinear::_jacobianDet[] = {
@@ -67,9 +69,7 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din3DLinear::_jacobianInv[] = {
- -1.33333333e+00,
-  5.71428571e-01,
-  4.00000000e-01,
+ -1.33333333e+00,  5.71428571e-01,  4.00000000e-01,
 };
 
 pylith::feassemble::QuadratureData1Din3DLinear::QuadratureData1Din3DLinear(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DQuadratic.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DQuadratic.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData1Din3DQuadratic.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -70,8 +70,12 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din3DQuadratic::_jacobian[] = {
- -6.92264973e-01,  1.80773503e+00,  2.84641016e+00,
- -8.07735027e-01,  1.69226497e+00,  2.15358984e+00,
+ -6.92264973e-01,
+  1.80773503e+00,
+  2.84641016e+00,
+ -8.07735027e-01,
+  1.69226497e+00,
+  2.15358984e+00,
 };
 
 const double pylith::feassemble::QuadratureData1Din3DQuadratic::_jacobianDet[] = {
@@ -79,12 +83,8 @@
 };
 
 const double pylith::feassemble::QuadratureData1Din3DQuadratic::_jacobianInv[] = {
- -1.44453358e+00,
-  5.53178417e-01,
-  3.51319713e-01,
- -1.23802976e+00,
-  5.90924008e-01,
-  4.64340973e-01,
+ -1.44453358e+00,  5.53178417e-01,  3.51319713e-01,
+ -1.23802976e+00,  5.90924008e-01,  4.64340973e-01,
 };
 
 pylith::feassemble::QuadratureData1Din3DQuadratic::QuadratureData1Din3DQuadratic(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DLinear.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DLinear.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DLinear.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -60,8 +60,8 @@
 };
 
 const double pylith::feassemble::QuadratureData2DLinear::_jacobian[] = {
-  1.00000000e-01,  9.00000000e-01,
- -1.20000000e+00,  2.00000000e-01,
+  1.00000000e-01, -1.20000000e+00,
+  9.00000000e-01,  2.00000000e-01,
 };
 
 const double pylith::feassemble::QuadratureData2DLinear::_jacobianDet[] = {
@@ -69,8 +69,8 @@
 };
 
 const double pylith::feassemble::QuadratureData2DLinear::_jacobianInv[] = {
-  1.81818182e-01, -8.18181818e-01,
-  1.09090909e+00,  9.09090909e-02,
+  1.81818182e-01,  1.09090909e+00,
+ -8.18181818e-01,  9.09090909e-02,
 };
 
 pylith::feassemble::QuadratureData2DLinear::QuadratureData2DLinear(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DQuadratic.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DQuadratic.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2DQuadratic.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -90,12 +90,12 @@
 };
 
 const double pylith::feassemble::QuadratureData2DQuadratic::_jacobian[] = {
-  2.33333333e+00,  1.60000000e+00,
-  8.33333333e-01,  2.03333333e+00,
-  3.33333333e+00,  1.20000000e+00,
- -1.66666667e-01,  1.63333333e+00,
-  2.33333333e+00,  8.00000000e-01,
- -1.66666667e-01,  1.63333333e+00,
+  2.33333333e+00,  8.33333333e-01,
+  1.60000000e+00,  2.03333333e+00,
+  3.33333333e+00, -1.66666667e-01,
+  1.20000000e+00,  1.63333333e+00,
+  2.33333333e+00, -1.66666667e-01,
+  8.00000000e-01,  1.63333333e+00,
 };
 
 const double pylith::feassemble::QuadratureData2DQuadratic::_jacobianDet[] = {
@@ -103,12 +103,12 @@
 };
 
 const double pylith::feassemble::QuadratureData2DQuadratic::_jacobianInv[] = {
-  5.96091205e-01, -4.69055375e-01,
- -2.44299674e-01,  6.84039088e-01,
-  2.89370079e-01, -2.12598425e-01,
-  2.95275591e-02,  5.90551181e-01,
-  4.14084507e-01, -2.02816901e-01,
-  4.22535211e-02,  5.91549296e-01,
+  5.96091205e-01, -2.44299674e-01,
+ -4.69055375e-01,  6.84039088e-01,
+  2.89370079e-01,  2.95275591e-02,
+ -2.12598425e-01,  5.90551181e-01,
+  4.14084507e-01,  4.22535211e-02,
+ -2.02816901e-01,  5.91549296e-01,
 };
 
 pylith::feassemble::QuadratureData2DQuadratic::QuadratureData2DQuadratic(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXY.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXY.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXY.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -60,8 +60,9 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXY::_jacobian[] = {
-  1.00000000e+00,  0.00000000e+00,  0.00000000e+00,
-  0.00000000e+00,  1.00000000e+00,  0.00000000e+00,
+  1.00000000e+00,  0.00000000e+00,
+  0.00000000e+00,  1.00000000e+00,
+  0.00000000e+00,  0.00000000e+00,
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXY::_jacobianDet[] = {
@@ -69,9 +70,8 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXY::_jacobianInv[] = {
-  1.00000000e+00,  0.00000000e+00,
-  0.00000000e+00,  1.00000000e+00,
-  0.00000000e+00,  0.00000000e+00,
+  1.00000000e+00,  0.00000000e+00,  0.00000000e+00,
+  1.00000000e+00,  0.00000000e+00,  0.00000000e+00,
 };
 
 pylith::feassemble::QuadratureData2Din3DLinearXY::QuadratureData2Din3DLinearXY(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXYZ.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXYZ.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXYZ.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -60,8 +60,9 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXYZ::_jacobian[] = {
-  2.50000000e+00,  2.50000000e+00,  5.00000000e-01,
- -1.50000000e+00,  4.00000000e+00,  4.50000000e+00,
+  2.50000000e+00, -1.50000000e+00,
+  2.50000000e+00,  4.00000000e+00,
+  5.00000000e-01,  4.50000000e+00,
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXYZ::_jacobianDet[] = {
@@ -69,9 +70,8 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXYZ::_jacobianInv[] = {
-  2.90909091e-01, -1.81818182e-01,
-  1.09090909e-01,  1.81818182e-01,
- -4.32432432e-01,  2.70270270e-01,
+  2.90909091e-01, -1.81818182e-01,  1.09090909e-01,
+  1.81818182e-01, -4.32432432e-01,  2.70270270e-01,
 };
 
 pylith::feassemble::QuadratureData2Din3DLinearXYZ::QuadratureData2Din3DLinearXYZ(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXZ.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXZ.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearXZ.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -60,8 +60,9 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXZ::_jacobian[] = {
- -1.00000000e+00,  0.00000000e+00,  0.00000000e+00,
-  0.00000000e+00,  0.00000000e+00,  1.00000000e+00,
+ -1.00000000e+00,  0.00000000e+00,
+  0.00000000e+00,  0.00000000e+00,
+  0.00000000e+00,  1.00000000e+00,
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXZ::_jacobianDet[] = {
@@ -69,9 +70,8 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearXZ::_jacobianInv[] = {
- -1.00000000e+00, -0.00000000e+00,
-  0.00000000e+00,  0.00000000e+00,
-  0.00000000e+00,  1.00000000e+00,
+ -1.00000000e+00,  0.00000000e+00,  0.00000000e+00,
+  0.00000000e+00, -0.00000000e+00,  1.00000000e+00,
 };
 
 pylith::feassemble::QuadratureData2Din3DLinearXZ::QuadratureData2Din3DLinearXZ(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearYZ.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearYZ.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DLinearYZ.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -60,8 +60,9 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearYZ::_jacobian[] = {
-  0.00000000e+00,  1.00000000e+00,  0.00000000e+00,
-  0.00000000e+00,  0.00000000e+00,  1.00000000e+00,
+  0.00000000e+00,  0.00000000e+00,
+  1.00000000e+00,  0.00000000e+00,
+  0.00000000e+00,  1.00000000e+00,
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearYZ::_jacobianDet[] = {
@@ -69,9 +70,8 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DLinearYZ::_jacobianInv[] = {
-  0.00000000e+00,  0.00000000e+00,
-  1.00000000e+00,  0.00000000e+00,
-  0.00000000e+00,  1.00000000e+00,
+  0.00000000e+00,  0.00000000e+00,  1.00000000e+00,
+  0.00000000e+00,  0.00000000e+00,  1.00000000e+00,
 };
 
 pylith::feassemble::QuadratureData2Din3DLinearYZ::QuadratureData2Din3DLinearYZ(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DQuadratic.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DQuadratic.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData2Din3DQuadratic.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -90,12 +90,15 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DQuadratic::_jacobian[] = {
- -1.56666667e+00,  3.56666667e+00,  5.66666667e-01,
- -2.36666667e+00,  3.66666667e-01,  2.36666667e+00,
- -1.36666667e+00,  3.36666667e+00,  3.66666667e-01,
- -2.56666667e+00,  5.66666667e-01,  2.56666667e+00,
- -1.36666667e+00,  3.36666667e+00,  3.66666667e-01,
- -2.36666667e+00,  3.66666667e-01,  2.36666667e+00,
+ -1.56666667e+00, -2.36666667e+00,
+  3.56666667e+00,  3.66666667e-01,
+  5.66666667e-01,  2.36666667e+00,
+ -1.36666667e+00, -2.56666667e+00,
+  3.36666667e+00,  5.66666667e-01,
+  3.66666667e-01,  2.56666667e+00,
+ -1.36666667e+00, -2.36666667e+00,
+  3.36666667e+00,  3.66666667e-01,
+  3.66666667e-01,  2.36666667e+00,
 };
 
 const double pylith::feassemble::QuadratureData2Din3DQuadratic::_jacobianDet[] = {
@@ -103,15 +106,12 @@
 };
 
 const double pylith::feassemble::QuadratureData2Din3DQuadratic::_jacobianInv[] = {
-  4.66101695e-02, -4.53389831e-01,
-  3.00847458e-01, -1.99152542e-01,
- -4.45344130e-02,  4.33198381e-01,
-  7.20338983e-02, -4.27966102e-01,
-  3.26271186e-01, -1.73728814e-01,
- -6.71936759e-02,  3.99209486e-01,
-  4.91071429e-02, -4.50892857e-01,
-  3.16964286e-01, -1.83035714e-01,
- -4.68085106e-02,  4.29787234e-01,
+  4.66101695e-02, -4.53389831e-01,  3.00847458e-01,
+ -1.99152542e-01, -4.45344130e-02,  4.33198381e-01,
+  7.20338983e-02, -4.27966102e-01,  3.26271186e-01,
+ -1.73728814e-01, -6.71936759e-02,  3.99209486e-01,
+  4.91071429e-02, -4.50892857e-01,  3.16964286e-01,
+ -1.83035714e-01, -4.68085106e-02,  4.29787234e-01,
 };
 
 pylith::feassemble::QuadratureData2Din3DQuadratic::QuadratureData2Din3DQuadratic(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DLinear.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DLinear.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DLinear.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -62,9 +62,9 @@
 };
 
 const double pylith::feassemble::QuadratureData3DLinear::_jacobian[] = {
-  2.50000000e+00,  5.00000000e-01,  1.00000000e-01,
-  1.50000000e+00,  9.00000000e-01,  2.00000000e-01,
-  3.00000000e-01,  1.50000000e+00,  2.50000000e+00,
+  2.50000000e+00,  1.50000000e+00,  3.00000000e-01,
+  5.00000000e-01,  9.00000000e-01,  1.50000000e+00,
+  1.00000000e-01,  2.00000000e-01,  2.50000000e+00,
 };
 
 const double pylith::feassemble::QuadratureData3DLinear::_jacobianDet[] = {
@@ -72,9 +72,9 @@
 };
 
 const double pylith::feassemble::QuadratureData3DLinear::_jacobianInv[] = {
-  6.04089219e-01, -3.40768278e-01,  3.09789343e-03,
- -1.14312268e+00,  1.92688971e+00, -1.08426270e-01,
-  6.13382900e-01, -1.11524164e+00,  4.64684015e-01,
+  6.04089219e-01, -1.14312268e+00,  6.13382900e-01,
+ -3.40768278e-01,  1.92688971e+00, -1.11524164e+00,
+  3.09789343e-03, -1.08426270e-01,  4.64684015e-01,
 };
 
 pylith::feassemble::QuadratureData3DLinear::QuadratureData3DLinear(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc	2007-05-15 00:51:39 UTC (rev 6878)
@@ -122,18 +122,18 @@
 };
 
 const double pylith::feassemble::QuadratureData3DQuadratic::_jacobian[] = {
-  2.36666667e+00, -3.00000000e-01,  3.66666667e-01,
-  1.66666667e+00,  3.03333333e+00,  1.00000000e+00,
-  5.66666667e-01,  2.40000000e+00,  2.96666667e+00,
-  2.63333333e+00,  2.33333333e-01,  6.33333333e-01,
-  1.66666667e+00,  3.03333333e+00,  1.26666667e+00,
-  8.33333333e-01,  2.66666667e+00,  2.96666667e+00,
-  2.36666667e+00, -3.00000000e-01,  6.33333333e-01,
-  1.40000000e+00,  3.03333333e+00,  1.00000000e-00,
-  5.66666667e-01,  2.66666667e+00,  2.70000000e+00,
-  2.63333333e+00, -3.33333333e-02,  3.66666667e-01,
-  1.66666667e+00,  3.30000000e+00,  7.33333333e-01,
-  8.33333333e-01,  2.66666667e+00,  2.96666667e+00,
+  2.36666667e+00,  1.66666667e+00,  5.66666667e-01,
+ -3.00000000e-01,  3.03333333e+00,  2.40000000e+00,
+  3.66666667e-01,  1.00000000e+00,  2.96666667e+00,
+  2.63333333e+00,  1.66666667e+00,  8.33333333e-01,
+  2.33333333e-01,  3.03333333e+00,  2.66666667e+00,
+  6.33333333e-01,  1.26666667e+00,  2.96666667e+00,
+  2.36666667e+00,  1.40000000e+00,  5.66666667e-01,
+ -3.00000000e-01,  3.03333333e+00,  2.66666667e+00,
+  6.33333333e-01,  1.00000000e-00,  2.70000000e+00,
+  2.63333333e+00,  1.66666667e+00,  8.33333333e-01,
+ -3.33333333e-02,  3.30000000e+00,  2.66666667e+00,
+  3.66666667e-01,  7.33333333e-01,  2.96666667e+00,
 };
 
 const double pylith::feassemble::QuadratureData3DQuadratic::_jacobianDet[] = {
@@ -141,18 +141,18 @@
 };
 
 const double pylith::feassemble::QuadratureData3DQuadratic::_jacobianInv[] = {
-  3.71410346e-01,  9.96222734e-02, -7.94851911e-02,
- -2.46397839e-01,  3.83480088e-01, -9.88092856e-02,
-  1.28389534e-01, -3.29260056e-01,  4.32196818e-01,
-  3.72043654e-01,  6.59662300e-02, -1.07590406e-01,
- -2.57393317e-01,  4.82134452e-01, -1.50906024e-01,
-  1.26858135e-01, -4.51909123e-01,  5.02946541e-01,
-  3.60726242e-01,  1.63201231e-01, -1.45059698e-01,
- -2.09861254e-01,  3.93888964e-01, -9.66580877e-02,
-  1.31562397e-01, -4.23278248e-01,  4.96279776e-01,
-  3.66156371e-01,  5.03198870e-02, -5.76939055e-02,
- -2.02525861e-01,  3.50837107e-01, -6.16924930e-02,
-  7.91928046e-02, -3.29493997e-01,  4.08738731e-01,
+  3.71410346e-01, -2.46397839e-01,  1.28389534e-01,
+  9.96222734e-02,  3.83480088e-01, -3.29260056e-01,
+ -7.94851911e-02, -9.88092856e-02,  4.32196818e-01,
+  3.72043654e-01, -2.57393317e-01,  1.26858135e-01,
+  6.59662300e-02,  4.82134452e-01, -4.51909123e-01,
+ -1.07590406e-01, -1.50906024e-01,  5.02946541e-01,
+  3.60726242e-01, -2.09861254e-01,  1.31562397e-01,
+  1.63201231e-01,  3.93888964e-01, -4.23278248e-01,
+ -1.45059698e-01, -9.66580877e-02,  4.96279776e-01,
+  3.66156371e-01, -2.02525861e-01,  7.91928046e-02,
+  5.03198870e-02,  3.50837107e-01, -3.29493997e-01,
+ -5.76939055e-02, -6.16924930e-02,  4.08738731e-01,
 };
 
 pylith::feassemble::QuadratureData3DQuadratic::QuadratureData3DQuadratic(void)

Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/feutils.py
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/feutils.py	2007-05-14 23:05:07 UTC (rev 6877)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/feutils.py	2007-05-15 00:51:39 UTC (rev 6878)
@@ -27,7 +27,7 @@
   @param vertices Coordinates of cell's vertices
   """
   jacobian = numpy.zeros( (quadrature.numQuadPts,
-                           quadrature.cellDim, quadrature.spaceDim),
+                           quadrature.spaceDim, quadrature.cellDim),
                           dtype=numpy.float64)
   jacobianInv = numpy.zeros( (quadrature.numQuadPts,
                               quadrature.spaceDim, quadrature.cellDim),
@@ -38,7 +38,7 @@
   for q in quadrature.quadPtsRef:
     # Jacobian at quadrature points
     deriv = quadrature.basisDeriv[iQuad]
-    j = numpy.dot(deriv.transpose(), vertices)
+    j = numpy.dot(vertices.transpose(), deriv)
     jacobian[iQuad] = j
 
     # Determinant of Jacobian and Jacobian inverse at quadrature points
@@ -46,56 +46,56 @@
       jacobianDet[iQuad] = numpy.linalg.det(j)
       jacobianInv[iQuad] = numpy.linalg.inv(j)
     else:
-      det = numpy.linalg.det(numpy.dot(j, j.transpose()))**0.5
+      det = numpy.linalg.det(numpy.dot(j.transpose(), j))**0.5
       jacobianDet[iQuad] = det
 
       if 1 == quadrature.cellDim:
-        jacobianInv[iQuad] = 1.0 / j.transpose()
+        jacobianInv[iQuad] = 1.0 / j
       elif 2 == quadrature.cellDim:
         minJacobian = 1.0e-06
-        jj01 = j[:,[0,1]]
-        jj12 = j[:,[1,2]]
-        jj02 = j[:,[0,2]]
-        det01 = numpy.linalg.det(jj01)
-        det12 = numpy.linalg.det(jj12)
-        det02 = numpy.linalg.det(jj02)
-        if abs(det01) > minJacobian:
-          ij01 = numpy.linalg.inv(jj01)
-          if abs(det12) > minJacobian:
-            ij12 = numpy.linalg.inv(jj12)
-            jacobianInv[iQuad] = numpy.array([ [ij01[0,0], ij01[0,1]],
-                                               [ij01[1,0], ij01[1,1]],
-                                               [ij12[1,0], ij12[1,1]] ],
+        jj10 = j[[0,1],:]
+        jj21 = j[[1,2],:]
+        jj20 = j[[0,2],:]
+        det10 = numpy.linalg.det(jj10)
+        det21 = numpy.linalg.det(jj21)
+        det20 = numpy.linalg.det(jj20)
+        if abs(det10) > minJacobian:
+          ij10 = numpy.linalg.inv(jj10)
+          if abs(det21) > minJacobian:
+            ij21 = numpy.linalg.inv(jj21)
+            jacobianInv[iQuad] = numpy.array([ [ij10[0,0], ij10[1,0]],
+                                               [ij10[0,1], ij10[1,1]],
+                                               [ij21[0,1], ij21[1,1]] ],
                                              dtype=numpy.float64)
-          elif abs(det02) > minJacobian:
-            ij02 = numpy.linalg.inv(jj02)
-            jacobianInv[iQuad] = numpy.array([ [ij01[0,0], ij01[0,1]],
-                                               [ij01[1,0], ij01[1,1]],
-                                               [ij02[1,0], ij02[1,1]] ],
+          elif abs(det20) > minJacobian:
+            ij20 = numpy.linalg.inv(jj20)
+            jacobianInv[iQuad] = numpy.array([ [ij10[0,0], ij10[1,0]],
+                                               [ij10[0,1], ij10[1,1]],
+                                               [ij20[0,1], ij20[1,1]] ],
                                              dtype=numpy.float64)
           else:
-            jacobianInv[iQuad] = numpy.array([ [ij01[0,0], ij01[0,1]],
-                                               [ij01[1,0], ij01[1,1]],
+            jacobianInv[iQuad] = numpy.array([ [ij10[0,0], ij10[0,1]],
+                                               [ij10[1,0], ij10[1,1]],
                                                [      0.0,       0.0] ],
                                              dtype=numpy.float64)
-        elif abs(det02) > minJacobian:
-          ij02 = numpy.linalg.inv(jj02)
-          if abs(det12) > minJacobian:
-            ij12 = numpy.linalg.inv(jj12)
-            jacobianInv[iQuad] = numpy.array([ [ij02[0,0], ij02[0,1]],
-                                               [ij12[0,0], ij12[0,1]],
-                                               [ij02[1,0], ij02[1,1]] ],
+        elif abs(det20) > minJacobian:
+          ij20 = numpy.linalg.inv(jj20)
+          if abs(det21) > minJacobian:
+            ij21 = numpy.linalg.inv(jj21)
+            jacobianInv[iQuad] = numpy.array([ [ij20[0,0], ij20[1,0]],
+                                               [ij21[0,0], ij21[1,0]],
+                                               [ij20[0,1], ij20[1,1]] ],
                                              dtype=numpy.float64)
           else:
-            jacobianInv[iQuad] = numpy.array([ [ij02[0,0], ij02[0,1]],
+            jacobianInv[iQuad] = numpy.array([ [ij20[0,0], ij20[1,0]],
                                                [      0.0,       0.0],
-                                               [ij02[1,0], ij02[1,1]] ],
+                                               [ij20[0,1], ij20[1,1]] ],
                                              dtype=numpy.float64)
-        elif abs(det12) > minJacobian:
-          ij12 = numpy.linalg.inv(jj12)
+        elif abs(det21) > minJacobian:
+          ij21 = numpy.linalg.inv(jj21)
           jacobianInv[iQuad] = numpy.array([ [      0.0,       0.0],
-                                             [ij12[0,0], ij12[0,1]],
-                                             [ij12[1,0], ij12[1,1]] ],
+                                             [ij21[0,0], ij21[1,0]],
+                                             [ij21[0,1], ij21[1,1]] ],
                                            dtype=numpy.float64)
         else:
           raise ValueError("Could not find inverse of Jacobian.")



More information about the cig-commits mailing list