[cig-commits] r4685 - in
short/3D/PyLith/trunk/unittests/libtests/feassemble: . data
brad at geodynamics.org
brad at geodynamics.org
Tue Oct 3 16:46:39 PDT 2006
Author: brad
Date: 2006-10-03 16:46:39 -0700 (Tue, 03 Oct 2006)
New Revision: 4685
Added:
short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc
short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.hh
Modified:
short/3D/PyLith/trunk/unittests/libtests/feassemble/Makefile.am
short/3D/PyLith/trunk/unittests/libtests/feassemble/data/Quadrature3DQuadratic.py
short/3D/PyLith/trunk/unittests/libtests/feassemble/data/generate.sh
Log:
Added unit test for quadratic element in 3-D.
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/Makefile.am 2006-10-03 22:09:38 UTC (rev 4684)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/Makefile.am 2006-10-03 23:46:39 UTC (rev 4685)
@@ -23,6 +23,7 @@
TestQuadrature1Din2D.cc \
TestQuadrature2D.cc \
TestQuadrature2Din3D.cc \
+ TestQuadrature3D.cc \
test_feassemble.cc \
data/QuadratureData.cc \
data/QuadratureData1DLinear.cc \
@@ -37,7 +38,9 @@
data/QuadratureData2Din3DLinearXY.cc \
data/QuadratureData2Din3DLinearYZ.cc \
data/QuadratureData2Din3DLinearXZ.cc \
- data/QuadratureData2Din3DQuadratic.cc
+ data/QuadratureData2Din3DQuadratic.cc \
+ data/QuadratureData3DLinear.cc \
+ data/QuadratureData3DQuadratic.cc
noinst_HEADERS = \
TestQuadrature.hh \
@@ -59,7 +62,9 @@
data/QuadratureData2Din3DLinearXY.hh \
data/QuadratureData2Din3DLinearYZ.hh \
data/QuadratureData2Din3DLinearXZ.hh \
- data/QuadratureData2Din3DQuadratic.hh
+ data/QuadratureData2Din3DQuadratic.hh \
+ data/QuadratureData3DLinear.hh \
+ data/QuadratureData3DQuadratic.hh
testfeassemble_LDFLAGS = $(PETSC_LIB)
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/Quadrature3DQuadratic.py
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/Quadrature3DQuadratic.py 2006-10-03 22:09:38 UTC (rev 4684)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/Quadrature3DQuadratic.py 2006-10-03 23:46:39 UTC (rev 4685)
@@ -21,16 +21,16 @@
# ----------------------------------------------------------------------
def N0(p):
- return (1.0-p[0]-p[1]-p[2])*(1-2.0*p[0]-2.0*p[1]-2.0*p[2])
+ return (1.0-p[0]-p[1]-p[2])*(1.0-2.0*p[0]-2.0*p[1]-2.0*p[2])
def N0p(p):
- return 3.0*p[0]+3.0*p[1]+3.0*p[2]+3.0
+ return 4.0*p[0]+4.0*p[1]+4.0*p[2]-3.0
def N0q(p):
- return 3.0*p[0]+3.0*p[1]+3.0*p[2]+3.0
+ return 4.0*p[0]+4.0*p[1]+4.0*p[2]-3.0
def N0r(p):
- return 3.0*p[0]+3.0*p[1]+3.0*p[2]+3.0
+ return 4.0*p[0]+4.0*p[1]+4.0*p[2]-3.0
def N1(p):
return p[0]*(2.0*p[0]-1.0)
@@ -105,40 +105,40 @@
return -8.0*p[2]+4.0*(1.0-p[0]-p[1])
def N7(p):
- return 4.0*p[0]*p[1]*(1.0-2.0*p[2])
+ return 4.0*p[0]*p[1]
def N7p(p):
- return 4.0*p[1]-8.0*p[1]*p[2]
+ return 4.0*p[1]
def N7q(p):
- return 4.0*p[0]-8.0*p[0]*p[2]
+ return 4.0*p[0]
def N7r(p):
- return -8.0*p[0]*p[1]
+ return 0.0
def N8(p):
- return 4.0*p[1]*p[2]*(1.0-2.0*p[0])
+ return 4.0*p[1]*p[2]
def N8p(p):
- return -8.0*p[1]*p[2]
+ return 0.0
def N8q(p):
- return 4.0*p[2]-8.0*p[0]*p[2]
+ return 4.0*p[2]
def N8r(p):
- return 4.0*p[1]-8.0*p[0]*p[1]
+ return 4.0*p[1]
def N9(p):
- return 4.0*p[0]*p[2]*(1.0-2.0*p[1])
+ return 4.0*p[0]*p[2]
def N9p(p):
- return 4.0*p[2]-8.0*p[1]*p[2]
+ return 4.0*p[2]
def N9q(p):
- return -8.0*p[0]*p[2]
+ return 0.0
def N9r(p):
- return 4.0*p[0]-8.0*p[0]*p[1]
+ return 4.0*p[0]
# ----------------------------------------------------------------------
@@ -157,18 +157,20 @@
"""
QuadratureApp.__init__(self, name)
- self.numVertices = 6
- self.spaceDim = 2
+ self.numVertices = 10
+ self.spaceDim = 3
self.numCells = 1
- self.cellDim = 2
- self.numCorners = 6
- self.numQuadPts = 3
-
- self.quadPtsRef = numpy.array( [[2.0/3.0, 1.0/6.0],
- [1.0/6.0, 2.0/3.0],
- [1.0/6.0, 1.0/6.0]],
+ self.cellDim = 3
+ self.numCorners = 10
+ self.numQuadPts = 4
+
+ # These are just approximate points used to test the quadrature routine
+ self.quadPtsRef = numpy.array( [[1.0/12.0, 1.0/12.0, 1.0/12.0],
+ [3.0/4.0, 1.0/12.0, 1.0/12.0],
+ [1.0/12.0, 3.0/4.0, 1.0/12.0],
+ [1.0/12.0, 1.0/12.0, 3.0/4.0]],
dtype=numpy.Float64)
- self.quadWts = numpy.array([1.0/6.0, 1.0/6.0, 1.0/6.0],
+ self.quadWts = numpy.array([1.0/8.0, 1.0/8.0, 1.0/8.0, 1.0/8.0],
dtype=numpy.Float64)
self.vertices = numpy.array( [[0.0, 0.0, 0.0],
[1.0, 0.0, 0.0],
@@ -181,8 +183,21 @@
[0.0, 0.5, 0.5],
[0.5, 0.0, 0.5]],
dtype=numpy.Float64)
- self.cells = numpy.array( [[0, 1, 2, 3, 4, 5]], dtype=numpy.Int32)
+ self.vertices = numpy.array( [[-0.5, -2.0, -1.0],
+ [ 2.0, -2.0, -0.5],
+ [ 1.0, 1.0, 0.0],
+ [ 0.2, 0.5, 2.0],
+ [ 0.7, -2.1, -0.8],
+ [ 0.3, -0.5, -0.5],
+ [-0.2, -0.8, 0.5],
+ [ 1.5, -0.6, -0.2],
+ [ 0.6, 0.8, 0.9],
+ [ 1.1, -0.8, 0.7]],
+ dtype=numpy.Float64)
+ self.cells = numpy.array( [[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]],
+ dtype=numpy.Int32)
+
return
@@ -202,17 +217,22 @@
iQuad = 0
for q in self.quadPtsRef:
# Basis functions at quadrature points
- basis = numpy.array([N0(q), N1(q), N2(q), N3(q), N4(q), N5(q)],
+ basis = numpy.array([N0(q), N1(q), N2(q), N3(q), N4(q),
+ N5(q), N6(q), N7(q), N8(q), N9(q)],
dtype=numpy.Float64)
self.basis[iQuad] = basis.reshape( (self.numCorners,) )
# Derivatives of basis functions at quadrature points
- deriv = numpy.array([[N0p(q), N0q(q)],
- [N1p(q), N1q(q)],
- [N2p(q), N2q(q)],
- [N3p(q), N3q(q)],
- [N4p(q), N4q(q)],
- [N5p(q), N5q(q)]],
+ deriv = numpy.array([[N0p(q), N0q(q), N0r(q)],
+ [N1p(q), N1q(q), N1r(q)],
+ [N2p(q), N2q(q), N2r(q)],
+ [N3p(q), N3q(q), N3r(q)],
+ [N4p(q), N4q(q), N4r(q)],
+ [N5p(q), N5q(q), N5r(q)],
+ [N6p(q), N6q(q), N6r(q)],
+ [N7p(q), N7q(q), N7r(q)],
+ [N8p(q), N8q(q), N8r(q)],
+ [N9p(q), N9q(q), N9r(q)]],
dtype=numpy.Float64)
self.basisDeriv[iQuad] = deriv.reshape((self.numCorners, self.cellDim))
Added: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc 2006-10-03 22:09:38 UTC (rev 4684)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.cc 2006-10-03 23:46:39 UTC (rev 4685)
@@ -0,0 +1,182 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+// DO NOT EDIT THIS FILE
+// This file was generated from python application quadrature3dquadratic.
+
+#include "QuadratureData3DQuadratic.hh"
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_numVertices = 10;
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_spaceDim = 3;
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_numCells = 1;
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_cellDim = 3;
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_numCorners = 10;
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_numQuadPts = 4;
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_vertices[] = {
+ -5.00000000e-01, -2.00000000e+00, -1.00000000e+00,
+ 2.00000000e+00, -2.00000000e+00, -5.00000000e-01,
+ 1.00000000e+00, 1.00000000e+00, 0.00000000e+00,
+ 2.00000000e-01, 5.00000000e-01, 2.00000000e+00,
+ 7.00000000e-01, -2.10000000e+00, -8.00000000e-01,
+ 3.00000000e-01, -5.00000000e-01, -5.00000000e-01,
+ -2.00000000e-01, -8.00000000e-01, 5.00000000e-01,
+ 1.50000000e+00, -6.00000000e-01, -2.00000000e-01,
+ 6.00000000e-01, 8.00000000e-01, 9.00000000e-01,
+ 1.10000000e+00, -8.00000000e-01, 7.00000000e-01,
+};
+
+const int pylith::feassemble::QuadratureData3DQuadratic::_cells[] = {
+ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9,
+};
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_quadPtsRef[] = {
+ 8.33333333e-02, 8.33333333e-02, 8.33333333e-02,
+ 7.50000000e-01, 8.33333333e-02, 8.33333333e-02,
+ 8.33333333e-02, 7.50000000e-01, 8.33333333e-02,
+ 8.33333333e-02, 8.33333333e-02, 7.50000000e-01,
+};
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_quadWts[] = {
+ 1.25000000e-01, 1.25000000e-01, 1.25000000e-01, 1.25000000e-01,
+};
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_basis[] = {
+ 3.75000000e-01, -6.94444444e-02, -6.94444444e-02,
+ -6.94444444e-02, 2.50000000e-01, 2.50000000e-01,
+ 2.50000000e-01, 2.77777778e-02, 2.77777778e-02,
+ 2.77777778e-02, -6.94444444e-02, 3.75000000e-01,
+ -6.94444444e-02, -6.94444444e-02, 2.50000000e-01,
+ 2.77777778e-02, 2.77777778e-02, 2.50000000e-01,
+ 2.77777778e-02, 2.50000000e-01, -6.94444444e-02,
+ -6.94444444e-02, 3.75000000e-01, -6.94444444e-02,
+ 2.77777778e-02, 2.50000000e-01, 2.77777778e-02,
+ 2.50000000e-01, 2.50000000e-01, 2.77777778e-02,
+ -6.94444444e-02, -6.94444444e-02, -6.94444444e-02,
+ 3.75000000e-01, 2.77777778e-02, 2.77777778e-02,
+ 2.50000000e-01, 2.77777778e-02, 2.50000000e-01,
+ 2.50000000e-01,};
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_basisDeriv[] = {
+ -2.00000000e+00, -2.00000000e+00, -2.00000000e+00,
+ -6.66666667e-01, 0.00000000e+00, 0.00000000e+00,
+ 0.00000000e+00, -6.66666667e-01, 0.00000000e+00,
+ 0.00000000e+00, 0.00000000e+00, -6.66666667e-01,
+ 2.66666667e+00, -3.33333333e-01, -3.33333333e-01,
+ -3.33333333e-01, 2.66666667e+00, -3.33333333e-01,
+ -3.33333333e-01, -3.33333333e-01, 2.66666667e+00,
+ 3.33333333e-01, 3.33333333e-01, 0.00000000e+00,
+ 0.00000000e+00, 3.33333333e-01, 3.33333333e-01,
+ 3.33333333e-01, 0.00000000e+00, 3.33333333e-01,
+ 6.66666667e-01, 6.66666667e-01, 6.66666667e-01,
+ 2.00000000e+00, 0.00000000e+00, 0.00000000e+00,
+ 0.00000000e+00, -6.66666667e-01, 0.00000000e+00,
+ 0.00000000e+00, 0.00000000e+00, -6.66666667e-01,
+ -2.66666667e+00, -3.00000000e+00, -3.00000000e+00,
+ -3.33333333e-01, 1.11022302e-16, -3.33333333e-01,
+ -3.33333333e-01, -3.33333333e-01, 1.11022302e-16,
+ 3.33333333e-01, 3.00000000e+00, 0.00000000e+00,
+ 0.00000000e+00, 3.33333333e-01, 3.33333333e-01,
+ 3.33333333e-01, 0.00000000e+00, 3.00000000e+00,
+ 6.66666667e-01, 6.66666667e-01, 6.66666667e-01,
+ -6.66666667e-01, 0.00000000e+00, 0.00000000e+00,
+ 0.00000000e+00, 2.00000000e+00, 0.00000000e+00,
+ 0.00000000e+00, 0.00000000e+00, -6.66666667e-01,
+ 1.11022302e-16, -3.33333333e-01, -3.33333333e-01,
+ -3.00000000e+00, -2.66666667e+00, -3.00000000e+00,
+ -3.33333333e-01, -3.33333333e-01, -1.11022302e-16,
+ 3.00000000e+00, 3.33333333e-01, 0.00000000e+00,
+ 0.00000000e+00, 3.33333333e-01, 3.00000000e+00,
+ 3.33333333e-01, 0.00000000e+00, 3.33333333e-01,
+ 6.66666667e-01, 6.66666667e-01, 6.66666667e-01,
+ -6.66666667e-01, 0.00000000e+00, 0.00000000e+00,
+ 0.00000000e+00, -6.66666667e-01, 0.00000000e+00,
+ 0.00000000e+00, 0.00000000e+00, 2.00000000e+00,
+ -1.11022302e-16, -3.33333333e-01, -3.33333333e-01,
+ -3.33333333e-01, -1.11022302e-16, -3.33333333e-01,
+ -3.00000000e+00, -3.00000000e+00, -2.66666667e+00,
+ 3.33333333e-01, 3.33333333e-01, 0.00000000e+00,
+ 0.00000000e+00, 3.00000000e+00, 3.33333333e-01,
+ 3.00000000e+00, 0.00000000e+00, 3.33333333e-01,
+};
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_quadPts[] = {
+ -1.20833333e-01, -1.58194444e+00, -6.40277778e-01,
+ 1.54583333e+00, -1.60416667e+00, -3.06944444e-01,
+ 9.01388889e-01, 4.40277778e-01, 2.63888889e-02,
+ 3.45833333e-01, 1.06944444e-01, 1.33750000e+00,
+};
+
+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,
+};
+
+const double pylith::feassemble::QuadratureData3DQuadratic::_jacobianDet[] = {
+ 1.77671111e+01, 1.51087407e+01, 1.53117037e+01, 2.13964444e+01,
+};
+
+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,
+};
+
+pylith::feassemble::QuadratureData3DQuadratic::QuadratureData3DQuadratic(void)
+{ // constructor
+ numVertices = _numVertices;
+ spaceDim = _spaceDim;
+ numCells = _numCells;
+ cellDim = _cellDim;
+ numCorners = _numCorners;
+ numQuadPts = _numQuadPts;
+ vertices = const_cast<double*>(_vertices);
+ cells = const_cast<int*>(_cells);
+ quadPtsRef = const_cast<double*>(_quadPtsRef);
+ quadWts = const_cast<double*>(_quadWts);
+ basis = const_cast<double*>(_basis);
+ basisDeriv = const_cast<double*>(_basisDeriv);
+ quadPts = const_cast<double*>(_quadPts);
+ jacobian = const_cast<double*>(_jacobian);
+ jacobianDet = const_cast<double*>(_jacobianDet);
+ jacobianInv = const_cast<double*>(_jacobianInv);
+} // constructor
+
+pylith::feassemble::QuadratureData3DQuadratic::~QuadratureData3DQuadratic(void)
+{}
+
+
+// End of file
Added: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.hh 2006-10-03 22:09:38 UTC (rev 4684)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/QuadratureData3DQuadratic.hh 2006-10-03 23:46:39 UTC (rev 4685)
@@ -0,0 +1,76 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+// DO NOT EDIT THIS FILE
+// This file was generated from python application quadrature3dquadratic.
+
+#if !defined(pylith_feassemble_quadraturedata3dquadratic_hh)
+#define pylith_feassemble_quadraturedata3dquadratic_hh
+
+#include "QuadratureData.hh"
+
+namespace pylith {
+ namespace feassemble {
+ class QuadratureData3DQuadratic;
+ } // pylith
+} // feassemble
+
+class pylith::feassemble::QuadratureData3DQuadratic : public QuadratureData
+{
+
+public:
+
+ /// Constructor
+ QuadratureData3DQuadratic(void);
+
+ /// Destructor
+ ~QuadratureData3DQuadratic(void);
+
+private:
+
+ static const int _numVertices;
+
+ static const int _spaceDim;
+
+ static const int _numCells;
+
+ static const int _cellDim;
+
+ static const int _numCorners;
+
+ static const int _numQuadPts;
+
+ static const double _vertices[];
+
+ static const int _cells[];
+
+ static const double _quadPtsRef[];
+
+ static const double _quadWts[];
+
+ static const double _basis[];
+
+ static const double _basisDeriv[];
+
+ static const double _quadPts[];
+
+ static const double _jacobian[];
+
+ static const double _jacobianDet[];
+
+ static const double _jacobianInv[];
+
+};
+
+#endif // pylith_feassemble_quadraturedata3dquadratic_hh
+
+// End of file
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/data/generate.sh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/data/generate.sh 2006-10-03 22:09:38 UTC (rev 4684)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/data/generate.sh 2006-10-03 23:46:39 UTC (rev 4685)
@@ -86,4 +86,9 @@
--data.object=QuadratureData3DLinear \
--data.parent=QuadratureData
+python Quadrature3DQuadratic.py \
+ --data.namespace=pylith,feassemble \
+ --data.object=QuadratureData3DQuadratic \
+ --data.parent=QuadratureData
+
# End of file
More information about the cig-commits
mailing list