[cig-commits] r14988 - in short/3D/PyLith/trunk: . libsrc/feassemble modulesrc/feassemble unittests/libtests/feassemble unittests/pytests/feassemble
brad at geodynamics.org
brad at geodynamics.org
Tue May 12 14:54:55 PDT 2009
Author: brad
Date: 2009-05-12 14:54:54 -0700 (Tue, 12 May 2009)
New Revision: 14988
Modified:
short/3D/PyLith/trunk/TODO
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.icc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.hh
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc
short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.hh
short/3D/PyLith/trunk/libsrc/feassemble/QuadratureEngine.hh
short/3D/PyLith/trunk/modulesrc/feassemble/Quadrature.i
short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.cc
short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.hh
short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature0D.cc
short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadratureEngine.cc
short/3D/PyLith/trunk/unittests/pytests/feassemble/TestMeshQuadrature.py
Log:
Added computeGeometry() for single cell to eliminate need to store quadrature information. Changed interface to quadrature engine compute geometry.
Modified: short/3D/PyLith/trunk/TODO
===================================================================
--- short/3D/PyLith/trunk/TODO 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/TODO 2009-05-12 21:54:54 UTC (rev 14988)
@@ -113,6 +113,9 @@
5. Tidy up
+ Add error message if number of time steps computed by time
+ stepping scheme is 0.
+
How do we determine the orientation for a fault in a 1-D mesh? We
assume normaldir is +1.0, but cohesive cells could be created so
that the fault has a normaldir of -1.0. If the normaldir is -1.0,
Modified: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorElasticity.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -101,6 +101,7 @@
sieveMesh->getLabelStratum("material-id", materialId);
// Compute geometry for quadrature operations.
+ _quadrature->initializeGeometry();
_quadrature->computeGeometry(mesh, cells);
// Initialize material.
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -72,12 +72,77 @@
} // copy constructor
// ----------------------------------------------------------------------
+// Setup quadrature engine.
template<typename mesh_type>
void
+pylith::feassemble::Quadrature<mesh_type>::initializeGeometry(void)
+{ // initializeGeometry
+ clear();
+ assert(0 == _engine);
+
+ const int cellDim = _cellDim;
+ const int spaceDim = _spaceDim;
+
+ if (1 == spaceDim)
+ if (1 == cellDim)
+ _engine = new Quadrature1D(*this);
+ else if (0 == cellDim)
+ _engine = new Quadrature0D(*this);
+ else {
+ std::cerr << "Unknown quadrature case with cellDim '"
+ << cellDim << "' and spaceDim '" << spaceDim << "'"
+ << std::endl;
+ assert(0);
+ } // if/else
+ else if (2 == spaceDim)
+ if (2 == cellDim)
+ _engine = new Quadrature2D(*this);
+ else if (1 == cellDim)
+ _engine = new Quadrature1Din2D(*this);
+ else if (0 == cellDim)
+ _engine = new Quadrature0D(*this);
+ else {
+ std::cerr << "Unknown quadrature case with cellDim '"
+ << cellDim << "' and spaceDim '" << spaceDim << "'"
+ << std::endl;
+ assert(0);
+ } // if/else
+ else if (3 == spaceDim)
+ if (3 == cellDim)
+ _engine = new Quadrature3D(*this);
+ else if (2 == cellDim)
+ _engine = new Quadrature2Din3D(*this);
+ else if (1 == cellDim)
+ _engine = new Quadrature1Din3D(*this);
+ else if (0 == cellDim)
+ _engine = new Quadrature0D(*this);
+ else {
+ std::cerr << "Unknown quadrature case with cellDim '"
+ << cellDim << "' and spaceDim '" << spaceDim << "'"
+ << std::endl;
+ assert(0);
+ } // if/else
+ else {
+ std::cerr << "Unknown quadrature case with cellDim '"
+ << cellDim << "' and spaceDim '" << spaceDim << "'"
+ << std::endl;
+ assert(0);
+ } // if/else
+
+ assert(0 != _engine);
+ _engine->initialize();
+} // initializeGeometry
+
+// ----------------------------------------------------------------------
+// Compute geometric quantities for each cell.
+template<typename mesh_type>
+void
pylith::feassemble::Quadrature<mesh_type>::computeGeometry(
const mesh_type& mesh,
const ALE::Obj<typename mesh_type::SieveMesh::label_sequence>& cells)
-{ // precomputeGeometry
+{ // computeGeometry
+ assert(0 != _engine);
+
typedef typename mesh_type::RealSection RealSection;
typedef typename mesh_type::SieveMesh::label_sequence label_sequence;
typedef typename mesh_type::RestrictVisitor RestrictVisitor;
@@ -86,10 +151,6 @@
ALE::MemoryLogger& logger = ALE::MemoryLogger::singleton();
logger.stagePush(loggingStage);
- clear();
- _setupEngine();
- assert(0 != _engine);
-
// Allocate field and cell buffer for quadrature points
int fiberDim = _numQuadPts * _spaceDim;
_quadPtsField = new topology::Field<mesh_type>(mesh);
@@ -148,9 +209,12 @@
const int numBasis = _numBasis;
const ALE::Obj<typename mesh_type::SieveMesh>& sieveMesh = mesh.sieveMesh();
assert(!sieveMesh.isNull());
+
+ double_array coordinatesCell(numBasis*_spaceDim);
const ALE::Obj<RealSection>& coordinates =
sieveMesh->getRealSection("coordinates");
- RestrictVisitor coordsVisitor(*coordinates, numBasis*_spaceDim);
+ RestrictVisitor coordsVisitor(*coordinates,
+ coordinatesCell.size(), &coordinatesCell[0]);
const ALE::Obj<RealSection>& quadPtsSection = _quadPtsField->section();
const ALE::Obj<RealSection>& jacobianSection = _jacobianField->section();
@@ -168,9 +232,7 @@
++c_iter) {
coordsVisitor.clear();
sieveMesh->restrictClosure(*c_iter, coordsVisitor);
- const double* cellVertexCoords = coordsVisitor.getValues();
- assert(0 != cellVertexCoords);
- _engine->computeGeometry(cellVertexCoords, _spaceDim, *c_iter);
+ _engine->computeGeometry(coordinatesCell, *c_iter);
// Update fields with cell data
quadPtsSection->updatePoint(*c_iter, &quadPts[0]);
@@ -231,66 +293,5 @@
delete _basisDerivField; _basisDerivField = 0;
} // clear
-// ----------------------------------------------------------------------
-// Setup quadrature engine.
-template<typename mesh_type>
-void
-pylith::feassemble::Quadrature<mesh_type>::_setupEngine(void)
-{ // clear
- delete _engine; _engine = 0;
- const int cellDim = _cellDim;
- const int spaceDim = _spaceDim;
-
- if (1 == spaceDim)
- if (1 == cellDim)
- _engine = new Quadrature1D(*this);
- else if (0 == cellDim)
- _engine = new Quadrature0D(*this);
- else {
- std::cerr << "Unknown quadrature case with cellDim '"
- << cellDim << "' and spaceDim '" << spaceDim << "'"
- << std::endl;
- assert(0);
- } // if/else
- else if (2 == spaceDim)
- if (2 == cellDim)
- _engine = new Quadrature2D(*this);
- else if (1 == cellDim)
- _engine = new Quadrature1Din2D(*this);
- else if (0 == cellDim)
- _engine = new Quadrature0D(*this);
- else {
- std::cerr << "Unknown quadrature case with cellDim '"
- << cellDim << "' and spaceDim '" << spaceDim << "'"
- << std::endl;
- assert(0);
- } // if/else
- else if (3 == spaceDim)
- if (3 == cellDim)
- _engine = new Quadrature3D(*this);
- else if (2 == cellDim)
- _engine = new Quadrature2Din3D(*this);
- else if (1 == cellDim)
- _engine = new Quadrature1Din3D(*this);
- else if (0 == cellDim)
- _engine = new Quadrature0D(*this);
- else {
- std::cerr << "Unknown quadrature case with cellDim '"
- << cellDim << "' and spaceDim '" << spaceDim << "'"
- << std::endl;
- assert(0);
- } // if/else
- else {
- std::cerr << "Unknown quadrature case with cellDim '"
- << cellDim << "' and spaceDim '" << spaceDim << "'"
- << std::endl;
- assert(0);
- } // if/else
-
- assert(0 != _engine);
- _engine->initialize();
-} // _setupEngine
-
-
// End of file
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -57,6 +57,9 @@
*/
Quadrature(const Quadrature& q);
+ /// Setup quadrature engine.
+ void initializeGeometry(void);
+
/** Set flag for checking ill-conditioning.
*
* @param flag True to check for ill-conditioning, false otherwise.
@@ -94,7 +97,7 @@
*/
const double_array& jacobianDet(void) const;
- /** Precompute geometric quantities for each cell.
+ /** Compute geometric quantities for each cell.
*
* @param mesh Finite-element mesh
* @param cells Finite-element cells for geometry.
@@ -112,12 +115,14 @@
/// Deallocate temporary storage.
void clear(void);
-// PRIVATE METHODS //////////////////////////////////////////////////////
-private :
+ /** Precompute geometric quantities for each cell.
+ *
+ * @param coordinatesCell Coordinates of vertices in cell.
+ * @param cell Finite-element cell
+ */
+ void computeGeometry(const double_array& coordinatesCell,
+ const int cell);
- /// Setup quadrature engine.
- void _setupEngine(void);
-
// PRIVATE MEMBERS //////////////////////////////////////////////////////
private :
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.icc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature.icc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -69,6 +69,17 @@
return _engine->jacobianDet();
}
+// Precompute geometric quantities for each cell.
+template<typename mesh_type>
+inline
+void
+pylith::feassemble::Quadrature<mesh_type>::computeGeometry(const double_array& coordinatesCell,
+ const int cell) {
+ assert(0 != _engine);
+ _engine->computeGeometry(coordinatesCell, cell);
+}
+
#endif
+
// End of file
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -44,14 +44,14 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature0D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature0D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
const int numQuadPts = _quadRefCell.numQuadPts();
const int numBasis = _quadRefCell.numBasis();
+ assert(coordinatesCell.size() == numBasis*spaceDim);
const double_array& basisDerivRef = _quadRefCell.basisDerivRef();
@@ -60,10 +60,9 @@
assert(1 == numBasis);
zero();
- assert(1 == coordDim);
+ assert(numBasis*1 == coordinatesCell.size());
- for (int i=0; i < spaceDim; ++i)
- _quadPts[i] = vertCoords[i];
+ _quadPts = coordinatesCell;
_jacobian[0] = 1.0;
_jacobianDet[0] = 1.0;
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature0D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -51,12 +51,10 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell);
// PROTECTED METHODS ////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -46,13 +46,9 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature1D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature1D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
- assert(0 != vertCoords);
- assert(1 == coordDim);
-
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
const int numQuadPts = _quadRefCell.numQuadPts();
@@ -63,6 +59,8 @@
const double_array& basisDerivRef = _quadRefCell.basisDerivRef();
const CellGeometry& geometry = _quadRefCell.refGeometry();
+ assert(numBasis*spaceDim == coordinatesCell.size());
+
assert(1 == cellDim);
assert(1 == spaceDim);
zero();
@@ -75,10 +73,10 @@
// x = sum[i=0,n-1] (Ni * xi)
for (int iBasis=0; iBasis < numBasis; ++iBasis)
_quadPts[iQuadPt] +=
- basis[iQuadPt*numBasis+iBasis]*vertCoords[iBasis];
+ basis[iQuadPt*numBasis+iBasis]*coordinatesCell[iBasis];
#else
geometry.coordsRefToGlobal(&_quadPts[iQuadPt], &quadPtsRef[iQuadPt],
- vertCoords, spaceDim);
+ &coordinatesCell[0], spaceDim);
#endif
#if defined(ISOPARAMETRIC)
@@ -86,7 +84,7 @@
// J = dx/dp = sum[i=0,n-1] (dNi/dp * xi)
for (int iBasis=0; iBasis < numBasis; ++iBasis)
_jacobian[iQuadPt] +=
- basisDerivRef[iQuadPt*numBasis+iBasis] * vertCoords[iBasis];
+ basisDerivRef[iQuadPt*numBasis+iBasis] * coordinatesCell[iBasis];
// Compute determinant of Jacobian at quadrature point
// |J| = j00
@@ -97,7 +95,7 @@
// Compute Jacobian and determinant of Jacobian at quadrature point
assert(0 != _geometry);
geometry->jacobian(&_jacobian[iQuadPt], &_jacobianDet[iQuadPt],
- vertCoords, &quadPtsRef[iQuadPt], spaceDim);
+ &coordinatesCell[0], &quadPtsRef[iQuadPt], spaceDim);
_checkJacobianDet(_jacobianDet[iQuadPt], cell);
#endif
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -42,12 +42,10 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell);
// PRIVATE METHODS //////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -46,13 +46,9 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature1Din2D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature1Din2D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
- assert(0 != vertCoords);
- assert(2 == coordDim);
-
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
const int numQuadPts = _quadRefCell.numQuadPts();
@@ -63,6 +59,8 @@
const double_array& basisDerivRef = _quadRefCell.basisDerivRef();
const CellGeometry& geometry = _quadRefCell.refGeometry();
+ assert(numBasis*spaceDim == coordinatesCell.size());
+
assert(1 == cellDim);
assert(2 == spaceDim);
zero();
@@ -78,12 +76,12 @@
const double valueBasis = basis[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_quadPts[iQuadPt*spaceDim+iDim] +=
- valueBasis * vertCoords[iBasis*spaceDim+iDim];
+ valueBasis * coordinatesCell[iBasis*spaceDim+iDim];
} // for
#else
geometry.coordsRefToGlobal(&_quadPts[iQuadPt*spaceDim],
&quadPtsRef[iQuadPt*cellDim],
- vertCoords, spaceDim);
+ &coordinatesCell[0], spaceDim);
#endif
#if defined(ISOPARAMETRIC)
@@ -96,7 +94,7 @@
const double deriv = basisDerivRef[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_jacobian[iQuadPt*spaceDim+iDim] +=
- deriv * vertCoords[iBasis*spaceDim+iDim];
+ deriv * coordinatesCell[iBasis*spaceDim+iDim];
} // for
// Compute determinant of Jacobian at quadrature point
@@ -112,7 +110,7 @@
// Compute Jacobian and determinant of Jacobian at quadrature point
geometry.jacobian(&_jacobian[iQuadPt*_cellDim*spaceDim],
&_jacobianDet[iQuadPt],
- vertCoords, &quadPtsRef[iQuadPt*cellDim], spaceDim);
+ &coordinatesCell0], &quadPtsRef[iQuadPt*cellDim], spaceDim);
_checkJacobianDet(_jacobianDet[iQuadPt], cell);
#endif
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din2D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -48,12 +48,10 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell);
// PROTECTED METHODS ////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -46,12 +46,9 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature1Din3D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature1Din3D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
- assert(0 != vertCoords);
- assert(3 == coordDim);
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
@@ -63,6 +60,8 @@
const double_array& basisDerivRef = _quadRefCell.basisDerivRef();
const CellGeometry& geometry = _quadRefCell.refGeometry();
+ assert(numBasis*spaceDim == coordinatesCell.size());
+
assert(1 == cellDim);
assert(3 == spaceDim);
zero();
@@ -79,12 +78,12 @@
const double valueBasis = basis[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_quadPts[iQuadPt*spaceDim+iDim] +=
- valueBasis * vertCoords[iBasis*spaceDim+iDim];
+ valueBasis * coordinatesCell[iBasis*spaceDim+iDim];
} // for
#else
geometry.coordsRefToGlobal(&_quadPts[iQuadPt*spaceDim],
&quadPtsRef[iQuadPt*cellDim],
- vertCoords, spaceDim);
+ &coordinatesCell[0], spaceDim);
#endif
#if defined(ISOPARAMETRIC)
@@ -99,7 +98,7 @@
const double deriv = basisDerivRef[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_jacobian[iQuadPt*spaceDim+iDim] +=
- deriv * vertCoords[iBasis*spaceDim+iDim];
+ deriv * coordinatesCell[iBasis*spaceDim+iDim];
} // for
// Compute determinant of Jacobian at quadrature point
@@ -115,7 +114,7 @@
// Compute Jacobian and determinant of Jacobian at quadrature point
geometry.jacobian(&_jacobian[iQuadPt*cellDim*spaceDim],
&_jacobianDet[iQuadPt],
- vertCoords, &quadPtsRef[iQuadPt*cellDim], spaceDim);
+ &coordinatesCell[0], &quadPtsRef[iQuadPt*cellDim], spaceDim);
_checkJacobianDet(_jacobianDet[iQuadPt], cell);
#endif
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature1Din3D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -48,12 +48,10 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell);
// PROTECTED METHODS ////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -46,13 +46,9 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature2D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature2D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
- assert(0 != vertCoords);
- assert(2 == coordDim);
-
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
const int numQuadPts = _quadRefCell.numQuadPts();
@@ -63,6 +59,7 @@
const double_array& basisDerivRef = _quadRefCell.basisDerivRef();
const CellGeometry& geometry = _quadRefCell.refGeometry();
+ assert(numBasis*spaceDim == coordinatesCell.size());
assert(2 == cellDim);
assert(2 == spaceDim);
zero();
@@ -78,12 +75,12 @@
const double valueBasis = basis[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_quadPts[iQuadPt*spaceDim+iDim] +=
- valueBasis * vertCoords[iBasis*spaceDim+iDim];
+ valueBasis * coordinatesCell[iBasis*spaceDim+iDim];
} // for
#else
geometry.coordsRefToGlobal(&_quadPts[iQuadPt*spaceDim],
&quadPtsRef[iQuadPt*cellDim],
- vertCoords, spaceDim);
+ &coordinatesCell[0], spaceDim);
#endif
#if defined(ISOPARAMETRIC)
@@ -100,7 +97,7 @@
basisDerivRef[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];
+ deriv * coordinatesCell[iBasis*spaceDim+iRow];
} // for
// Compute determinant of Jacobian at quadrature point
@@ -119,7 +116,7 @@
// Compute Jacobian and determinant of Jacobian at quadrature point
geometry.jacobian(&_jacobian[iQuadPt*cellDim*spaceDim],
&_jacobianDet[iQuadPt],
- vertCoords, &quadPtsRef[iQuadPt*cellDim], spaceDim);
+ &coordinatesCell[0], &quadPtsRef[iQuadPt*cellDim], spaceDim);
_checkJacobianDet(_jacobianDet[iQuadPt], cell);
const int iJ = iQuadPt*cellDim*spaceDim;
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -49,12 +49,10 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell);
// PROTECTED METHODS ////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -48,12 +48,9 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature2Din3D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature2Din3D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
- assert(0 != vertCoords);
- assert(3 == coordDim);
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
@@ -66,6 +63,7 @@
const CellGeometry& geometry = _quadRefCell.refGeometry();
const double minJacobian = _quadRefCell.minJacobian();
+ assert(numBasis*spaceDim == coordinatesCell.size());
assert(2 == cellDim);
assert(3 == spaceDim);
zero();
@@ -82,12 +80,12 @@
const double valueBasis = basis[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_quadPts[iQuadPt*spaceDim+iDim] +=
- valueBasis * vertCoords[iBasis*spaceDim+iDim];
+ valueBasis * coordinatesCell[iBasis*spaceDim+iDim];
} // for
#else
geometry.coordsRefToGlobal(&_quadPts[iQuadPt*spaceDim],
&quadPtsRef[iQuadPt*cellDim],
- vertCoords, spaceDim);
+ &coordinatesCell[0], spaceDim);
#endif
#if defined(ISOPARAMETRIC)
@@ -107,7 +105,7 @@
basisDerivRef[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];
+ deriv * coordinatesCell[iBasis*+spaceDim+iRow];
} // for
// Compute determinant of Jacobian at quadrature point
@@ -140,7 +138,7 @@
// Compute Jacobian and determinant of Jacobian at quadrature point
geometry.jacobian(&_jacobian[iQuadPt*cellDim*spaceDim],
&_jacobianDet[iQuadPt],
- vertCoords, &quadPtsRef[iQuadPt*cellDim], spaceDim);
+ &coordinatesCell[0], &quadPtsRef[iQuadPt*cellDim], spaceDim);
_checkJacobianDet(_jacobianDet[iQuadPt], cell);
const int iJ = iQuadPt*cellDim*spaceDim;
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature2Din3D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -49,12 +49,10 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell);
// PROTECTED METHODS ////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -46,12 +46,9 @@
// ----------------------------------------------------------------------
// Compute geometric quantities for a cell at quadrature points.
void
-pylith::feassemble::Quadrature3D::computeGeometry(const double* vertCoords,
- const int coordDim,
+pylith::feassemble::Quadrature3D::computeGeometry(const double_array& coordinatesCell,
const int cell)
{ // computeGeometry
- assert(0 != vertCoords);
- assert(3 == coordDim);
const int cellDim = _quadRefCell.cellDim();
const int spaceDim = _quadRefCell.spaceDim();
@@ -63,6 +60,7 @@
const double_array& basisDerivRef = _quadRefCell.basisDerivRef();
const CellGeometry& geometry = _quadRefCell.refGeometry();
+ assert(numBasis*spaceDim == coordinatesCell.size());
assert(3 == cellDim);
assert(3 == spaceDim);
zero();
@@ -79,12 +77,12 @@
const double valueBasis = basis[iQuadPt*numBasis+iBasis];
for (int iDim=0; iDim < spaceDim; ++iDim)
_quadPts[iQuadPt*spaceDim+iDim] +=
- valueBasis * vertCoords[iBasis*spaceDim+iDim];
+ valueBasis * coordinatesCell[iBasis*spaceDim+iDim];
} // for
#else
geometry.coordsRefToGlobal(&quadPts[iQuadPt*spaceDim],
&quadPtsRef[iQuadPt*cellDim],
- vertCoords, spaceDim);
+ &coordinatesCell[0], spaceDim);
#endif
#if defined(ISOPARAMETRIC)
@@ -107,7 +105,7 @@
basisDerivRef[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];
+ deriv * coordinatesCell[iBasis*spaceDim+iRow];
} // for
// Compute determinant of Jacobian at quadrature point
@@ -137,7 +135,7 @@
// Compute Jacobian and determinant of Jacobian at quadrature point
geometry.jacobian(&_jacobian[iQuadPt*cellDim*spaceDim],
&_jacobianDet[iQuadPt],
- vertCoords, &quadPtsRef[iQuadPt*cellDim], spaceDim);
+ &coordinatesCell[0], &quadPtsRef[iQuadPt*cellDim], spaceDim);
_checkJacobianDet(_jacobianDet[iQuadPt], cell);
const int iJ = iQuadPt*cellDim*spaceDim;
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Quadrature3D.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -42,13 +42,11 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param vertCoords Coordinates of vertices of finite-element cell.
- * @param coordDim Spatial dimension of coordinate system.
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
- void computeGeometry(const double* vertCoords,
- const int coordDim,
- const int cell);
+ void computeGeometry(const double_array& coordinatesCell,
+ const int cell);
// PROTECTED METHODS ////////////////////////////////////////////////////
protected :
Modified: short/3D/PyLith/trunk/libsrc/feassemble/QuadratureEngine.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/QuadratureEngine.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/libsrc/feassemble/QuadratureEngine.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -78,13 +78,11 @@
/** Compute geometric quantities for a cell at quadrature points.
*
- * @param mesh Finite-element mesh
- * @param coordinates Section containing vertex coordinates
+ * @param coordinatesCell Coordinates of cell's vertices.
* @param cell Finite-element cell
*/
virtual
- void computeGeometry(const double* vertCoords,
- const int coordDim,
+ void computeGeometry(const double_array& coordinatesCell,
const int cell) = 0;
// PROTECTED METHODS ////////////////////////////////////////////////////
Modified: short/3D/PyLith/trunk/modulesrc/feassemble/Quadrature.i
===================================================================
--- short/3D/PyLith/trunk/modulesrc/feassemble/Quadrature.i 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/modulesrc/feassemble/Quadrature.i 2009-05-12 21:54:54 UTC (rev 14988)
@@ -48,6 +48,9 @@
* @returns True if checking for ill-conditioning, false otherwise.
*/
bool checkConditioning(void) const;
+
+ /// Setup quadrature engine.
+ void initializeGeometry(void);
/// Deallocate temporary storage.
void clear(void);
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -224,6 +224,7 @@
const ALE::Obj<SieveMesh::label_sequence>& cells = sieveMesh->heightStratum(0);
CPPUNIT_ASSERT(!cells.isNull());
+ quadrature.initializeGeometry();
quadrature.computeGeometry(mesh, cells);
size_t size = 0;
@@ -274,5 +275,110 @@
quadrature.clear();
} // testComputeGeometry
+// ----------------------------------------------------------------------
+// Test computeGeometry() will coordinates and cell.
+void
+pylith::feassemble::TestQuadrature::testComputeGeometryCell(void)
+{ // testComputeGeometryCell
+ typedef pylith::topology::Mesh::SieveMesh SieveMesh;
+ QuadratureData2DLinear data;
+ const int cellDim = data.cellDim;
+ const int numBasis = data.numBasis;
+ const int numQuadPts = data.numQuadPts;
+ const int spaceDim = data.spaceDim;
+
+ const int numCells = data.numCells;
+ double_array vertCoords(data.vertices, numBasis*spaceDim);
+ const double* quadPtsE = data.quadPts;
+ const double* jacobianE = data.jacobian;
+ const double* jacobianDetE = data.jacobianDet;
+ const double* basisDerivE = data.basisDeriv;
+
+ const double minJacobian = 1.0e-06;
+
+#if 0
+ // Create mesh with test cell
+ topology::Mesh mesh(data.cellDim);
+ const ALE::Obj<SieveMesh>& sieveMesh = mesh.sieveMesh();
+ CPPUNIT_ASSERT(!sieveMesh.isNull());
+ ALE::Obj<SieveMesh::sieve_type> sieve =
+ new SieveMesh::sieve_type(mesh.comm());
+ CPPUNIT_ASSERT(!sieve.isNull());
+
+ // Cells and vertices
+ const bool interpolate = false;
+ ALE::Obj<ALE::Mesh::sieve_type> s =
+ new ALE::Mesh::sieve_type(sieve->comm(), sieve->debug());
+
+ ALE::SieveBuilder<ALE::Mesh>::buildTopology(s, cellDim, numCells,
+ const_cast<int*>(data.cells),
+ data.numVertices,
+ interpolate, numBasis);
+ std::map<ALE::Mesh::point_type,ALE::Mesh::point_type> renumbering;
+ ALE::ISieveConverter::convertSieve(*s, *sieve, renumbering);
+ sieveMesh->setSieve(sieve);
+ sieveMesh->stratify();
+ ALE::SieveBuilder<SieveMesh>::buildCoordinates(sieveMesh, spaceDim,
+ data.vertices);
+#endif
+
+ // Setup quadrature and compute geometry
+ GeometryTri2D geometry;
+ Quadrature<topology::Mesh> quadrature;
+ quadrature.refGeometry(&geometry);
+ quadrature.minJacobian(minJacobian);
+ quadrature.initialize(data.basis, numQuadPts, numBasis,
+ data.basisDerivRef, numQuadPts, numBasis, cellDim,
+ data.quadPtsRef, numQuadPts, cellDim,
+ data.quadWts, numQuadPts,
+ spaceDim);
+
+#if 0
+ const ALE::Obj<SieveMesh::label_sequence>& cells = sieveMesh->heightStratum(0);
+ CPPUNIT_ASSERT(!cells.isNull());
+#endif
+
+ quadrature.initializeGeometry();
+ quadrature.computeGeometry(vertCoords, 0);
+
+ size_t size = 0;
+
+ // Check values from computeGeometry()
+ const double tolerance = 1.0e-06;
+
+ const double_array& quadPts = quadrature.quadPts();
+ size = numQuadPts * spaceDim;
+ CPPUNIT_ASSERT_EQUAL(size, quadPts.size());
+ for (size_t i=0; i < size; ++i)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(quadPtsE[i], quadPts[i], tolerance);
+
+ const double_array& jacobian = quadrature.jacobian();
+ size = numQuadPts * cellDim * spaceDim;
+ CPPUNIT_ASSERT_EQUAL(size, jacobian.size());
+ for (size_t i=0; i < size; ++i)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(jacobianE[i], jacobian[i], tolerance);
+
+ const double_array& jacobianDet = quadrature.jacobianDet();
+ size = numQuadPts;
+ CPPUNIT_ASSERT_EQUAL(size, jacobianDet.size());
+ for (size_t i=0; i < size; ++i)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(jacobianDetE[i], jacobianDet[i], tolerance);
+
+ const double_array& basisDeriv = quadrature.basisDeriv();
+ size = numQuadPts * numBasis * spaceDim;
+ CPPUNIT_ASSERT_EQUAL(size, basisDeriv.size());
+ for (size_t i=0; i < size; ++i)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(basisDerivE[i], basisDeriv[i], tolerance);
+
+ quadrature.clear();
+
+ CPPUNIT_ASSERT(0 == quadrature._quadPtsField);
+ CPPUNIT_ASSERT(0 == quadrature._jacobianField);
+ CPPUNIT_ASSERT(0 == quadrature._jacobianDetField);
+ CPPUNIT_ASSERT(0 == quadrature._basisDerivField);
+ CPPUNIT_ASSERT(0 == quadrature._engine);
+} // testComputeGeometryCell
+
+
// End of file
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.hh 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature.hh 2009-05-12 21:54:54 UTC (rev 14988)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testCheckConditioning );
CPPUNIT_TEST( testEngineAccessors );
CPPUNIT_TEST( testComputeGeometry );
+ CPPUNIT_TEST( testComputeGeometryCell );
CPPUNIT_TEST_SUITE_END();
@@ -59,6 +60,9 @@
/// Test computeGeometry(), retrieveGeometry(), and clear().
void testComputeGeometry(void);
+ /// Test computeGeometry() with coordinates and cell.
+ void testComputeGeometryCell(void);
+
}; // class TestQuadrature
#endif // pylith_feassemble_testquadrature_hh
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature0D.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature0D.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadrature0D.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -45,13 +45,15 @@
const int numVertices = 1;
const int numCells = 1;
- const double vertCoords[] = { 1.1 };
+ const double vertCoordsData[] = { 1.1 };
const int cells[] = { 0 };
const double quadPts[] = { 1.1 };
const double jacobian[] = { 1.0 };
const double jacobianInv[] = { 1.0 };
const double jacobianDet[] = { 1.0 };
+ double_array vertCoords(vertCoordsData, numBasis*spaceDim);
+
const double minJacobian = 1.0e-06;
QuadratureRefCell refCell;
@@ -65,7 +67,7 @@
Quadrature0D engine(refCell);
engine.initialize();
- engine.computeGeometry(vertCoords, spaceDim, 0);
+ engine.computeGeometry(vertCoords, 0);
const double tolerance = 1.0e-06;
CPPUNIT_ASSERT_DOUBLES_EQUAL(quadPts[0], engine._quadPts[0],
Modified: short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadratureEngine.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadratureEngine.cc 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/unittests/libtests/feassemble/TestQuadratureEngine.cc 2009-05-12 21:54:54 UTC (rev 14988)
@@ -158,7 +158,7 @@
const int numVertices = data.numVertices;
const int numCells = data.numCells;
- const double* vertCoords = data.vertices;
+ const double_array vertCoords(data.vertices, numBasis*spaceDim);
const int* cells = data.cells;
const double* quadPtsE = data.quadPts;
const double* jacobianE = data.jacobian;
@@ -180,9 +180,8 @@
// Check values from computeGeometry()
engine->initialize();
- engine->computeGeometry(vertCoords, spaceDim, 0);
+ engine->computeGeometry(vertCoords, 0);
-
const double tolerance = 1.0e-06;
int size = numQuadPts * spaceDim;
for (int i=0; i < size; ++i)
Modified: short/3D/PyLith/trunk/unittests/pytests/feassemble/TestMeshQuadrature.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/feassemble/TestMeshQuadrature.py 2009-05-12 19:59:54 UTC (rev 14987)
+++ short/3D/PyLith/trunk/unittests/pytests/feassemble/TestMeshQuadrature.py 2009-05-12 21:54:54 UTC (rev 14988)
@@ -128,6 +128,8 @@
quadrature.quadPtsRef()))
self.assertTrue(TestArray_checkDouble(quadWtsE.ravel(),
quadrature.quadWts()))
+
+ quadrature.initializeGeometry()
return
More information about the CIG-COMMITS
mailing list