[cig-commits] r15887 - in short/3D/PyLith/branches/pylith-friction: . libsrc/feassemble pylith/problems unittests/libtests/feassemble
brad at geodynamics.org
brad at geodynamics.org
Tue Oct 27 16:00:37 PDT 2009
Author: brad
Date: 2009-10-27 16:00:36 -0700 (Tue, 27 Oct 2009)
New Revision: 15887
Modified:
short/3D/PyLith/branches/pylith-friction/TODO
short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Integrator.cc
short/3D/PyLith/branches/pylith-friction/pylith/problems/ExplicitLumped.py
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.cc
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DLinear.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DQuadratic.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DLinear.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DQuadratic.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DLinear.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DQuadratic.hh
short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestIntegrator.cc
Log:
Added unit test for lumped Jacobian.
Modified: short/3D/PyLith/branches/pylith-friction/TODO
===================================================================
--- short/3D/PyLith/branches/pylith-friction/TODO 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/TODO 2009-10-27 23:00:36 UTC (rev 15887)
@@ -11,7 +11,8 @@
Field split
uniform global refinement for tets with faults
customized SNES line search
-
+ restrictOperator() for sparse matrix
+ updateAllVisitor for Mesh
----------------------------------------------------------------------
FRICTION
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Integrator.cc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Integrator.cc 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/feassemble/Integrator.cc 2009-10-27 23:00:36 UTC (rev 15887)
@@ -169,8 +169,9 @@
for (int iBasis=0; iBasis < numBasis; ++iBasis)
for (int iDim=0; iDim < spaceDim; ++iDim) {
value = 0.0;
- for (int jBasis=0, indexI=iBasis*spaceDim; jBasis < numBasis; ++jBasis)
- value += _cellMatrix[indexI+jBasis*spaceDim+iDim];
+ const int indexR = (iBasis*spaceDim+iDim) * numBasis * spaceDim;
+ for (int jBasis=0; jBasis < numBasis; ++jBasis)
+ value += _cellMatrix[indexR+jBasis*spaceDim+iDim];
_cellVector[iBasis*spaceDim+iDim] = value;
} // for
Modified: short/3D/PyLith/branches/pylith-friction/pylith/problems/ExplicitLumped.py
===================================================================
--- short/3D/PyLith/branches/pylith-friction/pylith/problems/ExplicitLumped.py 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/pylith/problems/ExplicitLumped.py 2009-10-27 23:00:36 UTC (rev 15887)
@@ -121,6 +121,7 @@
jacobian.label("jacobian")
solution.vectorFieldType(solution.VECTOR)
jacobian.newSection(jacobian.VERTICES_FIELD, dimension)
+ jacobian.allocate()
self.jacobian = jacobian
self._debug.log(resourceUsageString())
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.cc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.cc 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.cc 2009-10-27 23:00:36 UTC (rev 15887)
@@ -248,6 +248,68 @@
} // testIntegrateJacobian
// ----------------------------------------------------------------------
+// Test integrateJacobian().
+void
+pylith::feassemble::TestElasticityExplicit::testIntegrateJacobianLumped(void)
+{ // testIntegrateJacobian
+ CPPUNIT_ASSERT(0 != _data);
+
+ topology::Mesh mesh;
+ ElasticityExplicit integrator;
+ topology::SolutionFields fields(mesh);
+ _initialize(&mesh, &integrator, &fields);
+ integrator._needNewJacobian = true;
+
+ topology::Field<topology::Mesh> jacobian(mesh);
+ jacobian.label("Jacobian");
+ jacobian.vectorFieldType(topology::FieldBase::VECTOR);
+ jacobian.newSection(topology::FieldBase::VERTICES_FIELD, _data->spaceDim);
+ jacobian.allocate();
+
+ const double t = 1.0;
+ integrator.integrateJacobian(jacobian, t, &fields);
+ CPPUNIT_ASSERT_EQUAL(false, integrator.needNewJacobian());
+ jacobian.complete();
+
+ const double* valsMatrixE = _data->valsJacobian;
+ const int sizeE = _data->numVertices * _data->spaceDim;
+ double_array valsE(sizeE);
+ const int spaceDim = _data->spaceDim;
+ const int numBasis = _data->numVertices;
+ for (int iBasis=0; iBasis < numBasis; ++iBasis)
+ for (int iDim=0; iDim < spaceDim; ++iDim) {
+ const int indexRow = (iBasis*spaceDim+iDim)*numBasis*spaceDim;
+ double value = 0.0;
+ for (int jBasis=0; jBasis < numBasis; ++jBasis)
+ value += valsMatrixE[indexRow + jBasis*spaceDim+iDim];
+ valsE[iBasis*spaceDim+iDim] = value;
+ } // for
+
+ // TEMPORARY
+ jacobian.view("JACOBIAN");
+ std::cout << "\n\nJACOBIAN FULL" << std::endl;
+ const int n = numBasis*spaceDim;
+ for (int r=0; r < n; ++r) {
+ for (int c=0; c < n; ++c)
+ std::cout << " " << valsMatrixE[r*n+c];
+ std::cout << "\n";
+ } // for
+
+ const ALE::Obj<RealSection>& jacobianSection = jacobian.section();
+ CPPUNIT_ASSERT(!jacobianSection.isNull());
+ const double* vals = jacobianSection->restrictSpace();
+ const int size = jacobianSection->sizeWithBC();
+ CPPUNIT_ASSERT_EQUAL(sizeE, size);
+
+ const double tolerance = 1.0e-06;
+ for (int i=0; i < size; ++i)
+ if (fabs(valsE[i]) > 1.0)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, vals[i]/valsE[i], tolerance);
+ else
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(valsE[i], vals[i], tolerance);
+} // testIntegrateJacobianLumped
+
+// ----------------------------------------------------------------------
// Test updateStateVars().
void
pylith::feassemble::TestElasticityExplicit::testUpdateStateVars(void)
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -89,6 +89,9 @@
/// Test integrateJacobian().
void testIntegrateJacobian(void);
+ /// Test integrateJacobianLumped().
+ void testIntegrateJacobianLumped(void);
+
/// Test updateStateVars().
void testUpdateStateVars(void);
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DLinear.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DLinear.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DLinear.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testInitialize );
CPPUNIT_TEST( testIntegrateResidual );
CPPUNIT_TEST( testIntegrateJacobian );
+ CPPUNIT_TEST( testIntegrateJacobianLumped );
CPPUNIT_TEST( testUpdateStateVars );
CPPUNIT_TEST( testStableTimeStep );
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DQuadratic.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DQuadratic.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit1DQuadratic.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testInitialize );
CPPUNIT_TEST( testIntegrateResidual );
CPPUNIT_TEST( testIntegrateJacobian );
+ CPPUNIT_TEST( testIntegrateJacobianLumped );
CPPUNIT_TEST( testUpdateStateVars );
CPPUNIT_TEST( testStableTimeStep );
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DLinear.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DLinear.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DLinear.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testInitialize );
CPPUNIT_TEST( testIntegrateResidual );
CPPUNIT_TEST( testIntegrateJacobian );
+ CPPUNIT_TEST( testIntegrateJacobianLumped );
CPPUNIT_TEST( testUpdateStateVars );
CPPUNIT_TEST( testStableTimeStep );
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DQuadratic.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DQuadratic.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit2DQuadratic.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testInitialize );
CPPUNIT_TEST( testIntegrateResidual );
CPPUNIT_TEST( testIntegrateJacobian );
+ CPPUNIT_TEST( testIntegrateJacobianLumped );
CPPUNIT_TEST( testUpdateStateVars );
CPPUNIT_TEST( testStableTimeStep );
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DLinear.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DLinear.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DLinear.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testInitialize );
CPPUNIT_TEST( testIntegrateResidual );
CPPUNIT_TEST( testIntegrateJacobian );
+ CPPUNIT_TEST( testIntegrateJacobianLumped );
CPPUNIT_TEST( testUpdateStateVars );
CPPUNIT_TEST( testStableTimeStep );
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DQuadratic.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DQuadratic.hh 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestElasticityExplicit3DQuadratic.hh 2009-10-27 23:00:36 UTC (rev 15887)
@@ -41,6 +41,7 @@
CPPUNIT_TEST( testInitialize );
CPPUNIT_TEST( testIntegrateResidual );
CPPUNIT_TEST( testIntegrateJacobian );
+ CPPUNIT_TEST( testIntegrateJacobianLumped );
CPPUNIT_TEST( testUpdateStateVars );
CPPUNIT_TEST( testStableTimeStep );
Modified: short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestIntegrator.cc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestIntegrator.cc 2009-10-27 22:15:06 UTC (rev 15886)
+++ short/3D/PyLith/branches/pylith-friction/unittests/libtests/feassemble/TestIntegrator.cc 2009-10-27 23:00:36 UTC (rev 15887)
@@ -208,9 +208,10 @@
const int numBasis = quadrature.numBasis();
const int spaceDim = quadrature.spaceDim();
- for (int iBasis=0, index=0; iBasis < numBasis; ++iBasis)
- for (int iDim=0; iDim < spaceDim; ++iDim, ++index) {
+ for (int iBasis=0; iBasis < numBasis; ++iBasis)
+ for (int iDim=0; iDim < spaceDim; ++iDim) {
double value = 0;
+ const int index = (iBasis*spaceDim+iDim)*numBasis*spaceDim;
for (int jBasis=0; jBasis < numBasis; ++jBasis)
value += 1.23 + 1.2*(index+jBasis*spaceDim+iDim);
CPPUNIT_ASSERT_EQUAL(value, integrator._cellVector[iBasis*spaceDim+iDim]);
More information about the CIG-COMMITS
mailing list