[cig-commits] r14308 - in short/3D/PyLith/branches/pylith-swig: . libsrc libsrc/utils modulesrc/feassemble modulesrc/utils pylith/feassemble unittests/pytests/feassemble
brad at geodynamics.org
brad at geodynamics.org
Thu Mar 12 12:40:31 PDT 2009
Author: brad
Date: 2009-03-12 12:40:31 -0700 (Thu, 12 Mar 2009)
New Revision: 14308
Added:
short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.cc
short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.hh
short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityImplicit.i
short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/IntegratorElasticity.i
short/3D/PyLith/branches/pylith-swig/modulesrc/utils/TestArray.i
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestSubMeshQuadrature.py
Removed:
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestIntegrator.py
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.hh
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.icc
Modified:
short/3D/PyLith/branches/pylith-swig/TODO
short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am
short/3D/PyLith/branches/pylith-swig/libsrc/utils/EventLogger.hh
short/3D/PyLith/branches/pylith-swig/libsrc/utils/Makefile.am
short/3D/PyLith/branches/pylith-swig/libsrc/utils/utilsfwd.hh
short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/Makefile.am
short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/QuadratureRefCell.i
short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i
short/3D/PyLith/branches/pylith-swig/modulesrc/utils/Makefile.am
short/3D/PyLith/branches/pylith-swig/modulesrc/utils/utils.i
short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py
short/3D/PyLith/branches/pylith-swig/pylith/feassemble/Integrator.py
short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/Makefile.am
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityImplicit.py
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestMeshQuadrature.py
short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/testfeassemble.py
Log:
Started updating Python integrators and corresponding unit tests.
Modified: short/3D/PyLith/branches/pylith-swig/TODO
===================================================================
--- short/3D/PyLith/branches/pylith-swig/TODO 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/TODO 2009-03-12 19:40:31 UTC (rev 14308)
@@ -4,6 +4,10 @@
0. SWIG conversion
+ Cleanup logging. Constraints and Integrators should log at the C++
+ level using the C++ EventLogger. Add finer grain logging at C++
+ level as in ElasticityImplicit.
+
TestMesh.test_view()
TestMesh.test_checkMaterialIds()
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/Makefile.am 2009-03-12 19:40:31 UTC (rev 14308)
@@ -81,7 +81,8 @@
topology/MeshOps.cc \
topology/SubMesh.cc \
topology/SolutionFields.cc \
- utils/EventLogger.cc
+ utils/EventLogger.cc \
+ utils/TestArray.cc
# faults/BruneSlipFn.cc \
# faults/ConstRateSlipFn.cc \
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/utils/EventLogger.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/utils/EventLogger.hh 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/utils/EventLogger.hh 2009-03-12 19:40:31 UTC (rev 14308)
@@ -21,18 +21,15 @@
#if !defined(pylith_utils_eventlogger_hh)
#define pylith_utils_eventlogger_hh
+// Include directives ---------------------------------------------------
+#include "utilsfwd.hh" // forward declarations
+
#include <string> // USES std::string
#include <map> // USES std::map
#include "petsclog.h" // USES PetscLogEventBegin/End() in inline methods
-namespace pylith {
- namespace utils {
- class EventLogger;
- class TestEventLogger; // unit testing
- } // utils
-} // pylith
-
+// EventLogger ----------------------------------------------------------
class pylith::utils::EventLogger
{ // EventLogger
friend class TestEventLogger; // unit testing
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/utils/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/utils/Makefile.am 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/utils/Makefile.am 2009-03-12 19:40:31 UTC (rev 14308)
@@ -16,6 +16,7 @@
subpkginclude_HEADERS = \
EventLogger.hh \
EventLogger.icc \
+ TestArray.hh \
array.hh \
arrayfwd.hh \
constdefs.h \
Added: short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.cc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.cc (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.cc 2009-03-12 19:40:31 UTC (rev 14308)
@@ -0,0 +1,57 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+#include "TestArray.hh" // implementation of class methods
+
+#include "array.hh" // USES double_array
+
+#include <iostream> // USES std::cerr
+#include <cassert> // USES assert()
+
+// ----------------------------------------------------------------------
+// Check to make sure array of values match expected values.
+bool
+pylith::utils::TestArray::check(const double* valuesE,
+ const int nvalues,
+ const double_array& values)
+{ // check(double)
+ assert( (0 == nvalues && 0 == valuesE) ||
+ (0 < nvalues && 0 != valuesE) );
+
+ if (nvalues != values.size()) {
+ std::cerr << "Array size mismatch, expected: " << nvalues
+ << " actual: " << values.size() << std::endl;
+ return false;
+ } // if
+
+ const double tolerance = 1.0e-06;
+ bool okay = true;
+ for (int i=0; i < nvalues; ++i) {
+ okay = true;
+ if (0.0 != valuesE[i]) {
+ if (fabs(1.0 - values[i]/valuesE[i]) > tolerance)
+ okay = false;
+ } else if (fabs(values[i] - valuesE[i]) > tolerance)
+ okay = false;
+
+ if (!okay) {
+ std::cerr << "Mismatch in array at index " << i << ", expected: "
+ << valuesE[i] << ", actual: " << values[i] << std::endl;
+ return false;
+ } // if
+ } // for
+
+ return true;
+} // check(double)
+
+
+// End of file
Added: short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.hh (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/utils/TestArray.hh 2009-03-12 19:40:31 UTC (rev 14308)
@@ -0,0 +1,57 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+/**
+ * @file pylith/utils/TestArray.hh
+ *
+ * @brief C++ object for testing array values.
+ *
+ * This object is used in unit testing of SWIG interfaces where the
+ * C++ object has an accessor returning a std::valarray. The TestArray
+ * methods provide the ability to compare the array returned by the
+ * accessor against the expected values, which are supplied via a
+ * pointer and a size (number of values).
+ */
+
+#if !defined(pylith_utils_testarray_hh)
+#define pylith_utils_testarray_hh
+
+// Include directives ---------------------------------------------------
+#include "utilsfwd.hh" // forward declarations
+
+#include "arrayfwd.hh" // USES double_array
+
+// TestArray ------------------------------------------------------------
+class pylith::utils::TestArray
+{ // TestArray
+
+// PUBLIC MEMBERS ///////////////////////////////////////////////////////
+public :
+
+ /** Check to make sure array of values match expected values.
+ *
+ * @param valuesE Array of expected values.
+ * @param nvalues Array size.
+ * @param values Array of values to check.
+ */
+ static
+ bool
+ check(const double* valuesE,
+ const int nvalues,
+ const double_array& values);
+
+}; // EventLogger
+
+#endif // pylith_utils_testarray_hh
+
+
+// End of file
Modified: short/3D/PyLith/branches/pylith-swig/libsrc/utils/utilsfwd.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/libsrc/utils/utilsfwd.hh 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/libsrc/utils/utilsfwd.hh 2009-03-12 19:40:31 UTC (rev 14308)
@@ -26,6 +26,8 @@
class EventLogger;
+ class TestArray;
+
} // utils
} // pylith
Added: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityImplicit.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityImplicit.i (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityImplicit.i 2009-03-12 19:40:31 UTC (rev 14308)
@@ -0,0 +1,89 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ----------------------------------------------------------------------
+//
+
+/** @file modulesrc/feassemble/ElasticityImplicit.i
+ *
+ * @brief Python interface to C++ ElasticityImplicit object.
+ */
+
+namespace pylith {
+ namespace feassemble {
+
+ class ElasticityImplicit : public IntegratorElasticity
+ { // ElasticityImplicit
+
+ // PUBLIC MEMBERS /////////////////////////////////////////////////
+ public :
+
+ /// Constructor
+ ElasticityImplicit(void);
+
+ /// Destructor
+ ~ElasticityImplicit(void);
+
+ /** Set time step for advancing from time t to time t+dt.
+ *
+ * @param dt Time step
+ */
+ void timeStep(const double dt);
+
+ /** Get stable time step for advancing from time t to time t+dt.
+ *
+ * Default is current time step.
+ *
+ * @returns Time step
+ */
+ double stableTimeStep(void) const;
+
+ /** Set flag for setting constraints for total field solution or
+ * incremental field solution.
+ *
+ * @param flag True if using incremental solution, false otherwise.
+ */
+ void useSolnIncr(const bool flag);
+
+ /** Integrate residual part of RHS for 3-D finite elements.
+ * Includes gravity and element internal force contribution.
+ *
+ * We assume that the effects of boundary conditions are already
+ * included in the residual (tractions, concentrated nodal forces,
+ * and contributions to internal force vector due to
+ * displacement/velocity BC). This routine computes the additional
+ * external loads due to body forces plus the
+ * element internal forces for the current stress state.
+ *
+ * @param residual Field containing values for residual
+ * @param t Current time
+ * @param fields Solution fields
+ */
+ void integrateResidual(const pylith::topology::Field<pylith::topology::Mesh>& residual,
+ const double t,
+ pylith::topology::SolutionFields* const fields);
+
+ /** Integrate contributions to Jacobian matrix (A) associated with
+ * operator.
+ *
+ * @param jacobian Sparse matrix for Jacobian of system.
+ * @param t Current time
+ * @param fields Solution fields
+ */
+ void integrateJacobian(PetscMat* jacobian,
+ const double t,
+ pylith::topology::SolutionFields* const fields);
+
+ }; // ElasticityImplicit
+
+ } // feassemble
+} // pylith
+
+
+// End of file
Added: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/IntegratorElasticity.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/IntegratorElasticity.i (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/IntegratorElasticity.i 2009-03-12 19:40:31 UTC (rev 14308)
@@ -0,0 +1,82 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ----------------------------------------------------------------------
+//
+
+/** @file modulesrc/feassemble/Integrator.i
+ *
+ * @brief Python interface to C++ abstract Integrator object.
+ */
+
+%template(MeshIntegrator) pylith::feassemble::Integrator<pylith::feassemble::Quadrature<pylith::topology::Mesh> >;
+
+namespace pylith {
+ namespace feassemble {
+
+ class IntegratorElasticity :
+ public pylith::feassemble::Integrator<pylith::feassemble::Quadrature<pylith::topology::Mesh> >
+ { // IntegratorElasticity
+
+ // PUBLIC MEMBERS /////////////////////////////////////////////////
+ public :
+
+ /// Constructor
+ IntegratorElasticity(void);
+
+ /// Destructor
+ ~IntegratorElasticity(void);
+
+ /** Set material.
+ *
+ * @param m Elastic material.
+ */
+ void material(pylith::materials::ElasticMaterial* m);
+
+ /** Determine whether we need to recompute the Jacobian.
+ *
+ * @returns True if Jacobian needs to be recomputed, false otherwise.
+ */
+ bool needNewJacobian(void);
+
+ /** Set flag for setting constraints for total field solution or
+ * incremental field solution.
+ *
+ * @param flag True if using incremental solution, false otherwise.
+ */
+ void useSolnIncr(const bool flag);
+
+ /** Initialize integrator.
+ *
+ * @param mesh Finite-element mesh.
+ */
+ void initialize(const pylith::topology::Mesh& mesh);
+
+ /** Update state variables as needed.
+ *
+ * @param t Current time
+ * @param fields Solution fields
+ * @param mesh Finite-element mesh
+ */
+ void updateStateVars(const double t,
+ pylith::topology::SolutionFields* const fields);
+
+ /** Verify configuration is acceptable.
+ *
+ * @param mesh Finite-element mesh
+ */
+ void verifyConfiguration(const pylith::topology::Mesh& mesh) const;
+
+ }; // IntegratorElasticity
+
+ } // feassemble
+} // pylith
+
+
+// End of file
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/Makefile.am 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/Makefile.am 2009-03-12 19:40:31 UTC (rev 14308)
@@ -29,7 +29,9 @@
GeometryTri2D.i \
GeometryTri3D.i \
QuadratureRefCell.i \
- Quadrature.i
+ Quadrature.i \
+ Integrator.i \
+ IntegratorElasticity.i
swig_generated = \
feassemble_wrap.cxx \
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/QuadratureRefCell.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/QuadratureRefCell.i 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/QuadratureRefCell.i 2009-03-12 19:40:31 UTC (rev 14308)
@@ -99,8 +99,10 @@
const int numQuadPts4,
const int spaceDim);
%clear(const double* basis, const int numQuadPts, const int numBasis);
- %clear(const double* basisDerivRef, const int numQuadPts, const int numBasis, const int spaceDim);
- %clear(const double* quadPtsRef, const int numQuadPts, const int cellDim);
+ %clear(const double* basisDerivRef, const int numQuadPts,
+ const int numBasis, const int spaceDim);
+ %clear(const double* quadPtsRef, const int numQuadPts,
+ const int cellDim);
%clear(const double* quadWts, const int numQuadPts);
/** Set geometry associated with reference cell.
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i 2009-03-12 19:40:31 UTC (rev 14308)
@@ -33,7 +33,8 @@
#include "pylith/topology/Mesh.hh"
#include "pylith/topology/SubMesh.hh"
#include "pylith/feassemble/Quadrature.hh"
-#include "pylith/feassemble/Integrator.hh"
+#include "pylith/feassemble/ElasticityImplicit.hh"
+#include "pylith/feassemble/ElasticityExplicit.hh"
%}
@@ -75,6 +76,8 @@
%include "Quadrature.i"
%include "Integrator.i"
+%include "IntegratorElasticity.i"
+%include "ElasticityImplicit.i"
// Template instatiation
%template(MeshQuadrature) pylith::feassemble::Quadrature<pylith::topology::Mesh>;
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/utils/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/utils/Makefile.am 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/utils/Makefile.am 2009-03-12 19:40:31 UTC (rev 14308)
@@ -32,7 +32,8 @@
utils_swig_sources = \
utils.i \
- EventLogger.i
+ EventLogger.i \
+ TestArray.i
utils_swig_generated = \
utils_wrap.cxx \
@@ -73,7 +74,7 @@
endif
-INCLUDES += $(PYTHON_EGG_CPPFLAGS) -I$(PYTHON_INCDIR) $(PETSC_INCLUDE)
+INCLUDES += $(PYTHON_EGG_CPPFLAGS) -I$(NUMPY_INCDIR) -I$(PYTHON_INCDIR) $(PETSC_INCLUDE)
$(srcdir)/petsc_wrap.cxx $(srcdir)/petsc.py: $(petsc_swig_sources)
$(SWIG) -Wall -c++ -python $<
Added: short/3D/PyLith/branches/pylith-swig/modulesrc/utils/TestArray.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/utils/TestArray.i (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/utils/TestArray.i 2009-03-12 19:40:31 UTC (rev 14308)
@@ -0,0 +1,40 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+/**
+ * @file modulesrc/utils/TestArray.i
+ *
+ * @brief Python interface to C++ TestArray object.
+ */
+
+%apply(double* IN_ARRAY1, int DIM1) {
+ (const double* valuesE,
+ const int nvalues)
+ };
+%inline %{
+ /** Check to make sure array of values match expected values.
+ *
+ * @param valuesE Array of expected values.
+ * @param nvalues Array size.
+ * @param values Array of values to check.
+ */
+ bool
+ TestArray_checkDouble(const double* valuesE,
+ const int nvalues,
+ const pylith::double_array& values) {
+ pylith::utils::TestArray::check(valuesE, nvalues, values);
+ } // check(double)
+%} // inline
+%clear(const double* valuesE, const int nvalues);
+
+
+// End of file
Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/utils/utils.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/utils/utils.i 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/utils/utils.i 2009-03-12 19:40:31 UTC (rev 14308)
@@ -16,8 +16,10 @@
// Header files for module C++ code
%{
#include "pylith/utils/EventLogger.hh"
+#include "pylith/utils/TestArray.hh"
#include <petsclog.h> // USES PetscLogEventBegin/End() in inline methods
+#include "pylith/utils/arrayfwd.hh" // USES double_array
%}
%include "exception.i"
@@ -29,8 +31,20 @@
} // try/catch
} // exception
+%include "typemaps.i"
+
+// Numpy interface stuff
+%{
+#define SWIG_FILE_WITH_INIT
+%}
+%include "../include/numpy.i"
+%init %{
+import_array();
+%}
+
// Interfaces
%include "EventLogger.i"
+%include "TestArray.i"
// End of file
Modified: short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -18,11 +18,12 @@
## Factory: integrator
from IntegratorElasticity import IntegratorElasticity
+from feassemble import ElasticityImplicit as ModuleElasticityImplicit
# ElasticityImplicit class
-class ElasticityImplicit(IntegratorElasticity):
+class ElasticityImplicit(IntegratorElasticity, ModuleElasticityImplicit):
"""
- Python object for implicit time integration of dynamic elasticity
+ Python object for implicit time integration of elasticity
equation using finite-elements.
"""
@@ -33,10 +34,8 @@
Constructor.
"""
IntegratorElasticity.__init__(self, name)
+ ModuleElasticityImplicit.__init__(self)
self._loggingPrefix = "ElIm "
-
- import pylith.feassemble.feassemble as bindings
- self.cppHandle = bindings.ElasticityImplicit()
return
Modified: short/3D/PyLith/branches/pylith-swig/pylith/feassemble/Integrator.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/feassemble/Integrator.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/pylith/feassemble/Integrator.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -57,7 +57,6 @@
"""
Constructor.
"""
- self.quadrature = None
self.gravityField = None
self.mesh = None
return
@@ -67,25 +66,9 @@
"""
Do pre-initialization setup.
"""
- self._setupLogging()
-
return
- def verifyConfiguration(self):
- """
- Verify compatibility of configuration.
- """
- logEvent = "%sverify" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- self.cppHandle.verifyConfiguration(self.mesh.cppHandle)
-
- self._logger.eventEnd(logEvent)
- return
-
-
def initialize(self, totalTime, numTimeSteps, normalizer):
"""
Do initialization.
@@ -93,128 +76,14 @@
logEvent = "%sinit" % self._loggingPrefix
self._logger.eventBegin(logEvent)
- assert(None != self.cppHandle)
- self.cppHandle.normalizer = normalizer.cppHandle
+ self.normalizer(normalizer)
if None != self.gravityField:
- self.cppHandle.gravityField = self.gravityField.cppHandle
+ self.gravityField(self.gravityField)
self._logger.eventEnd(logEvent)
return
- def timeStep(self, dt):
- """
- Set time step for advancing from time t to time t+dt.
- """
- assert(None != self.cppHandle)
- self.cppHandle.timeStep = dt
- return
-
-
- def stableTimeStep(self):
- """
- Get stable time step for advancing from time t to time t+dt.
- """
- logEvent = "%stimestep" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- dt = self.cppHandle.stableTimeStep
-
- self._logger.eventEnd(logEvent)
- return dt
-
-
- def useSolnIncr(self, flag):
- """
- Set behavior for using total field solution or incremental field solution.
- """
- logEvent = "%ssolnIncr" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- self.cppHandle.useSolnIncr = flag
-
- self._logger.eventEnd(logEvent)
- return
-
-
- def needNewJacobian(self):
- """
- Returns true if we need to recompute Jacobian matrix for operator,
- false otherwise.
- """
- logEvent = "%snewJacobian" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- flag = self.cppHandle.needNewJacobian
- self._logger.eventEnd(logEvent)
- return flag
-
-
- def integrateResidual(self, residual, t, fields):
- """
- Integrate contributions to residual term at time t.
- """
- logEvent = "%sresidual" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- self.cppHandle.integrateResidual(residual, t, fields.cppHandle,
- self.mesh.cppHandle,
- self.mesh.coordsys.cppHandle)
- self._logger.eventEnd(logEvent)
- return
-
-
- def integrateJacobian(self, jacobian, t, fields):
- """
- Integrate contributions to Jacobian term at time t.
- """
- logEvent = "%sjacobian" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- self.cppHandle.integrateJacobian(jacobian, t, fields.cppHandle,
- self.mesh.cppHandle)
- self._logger.eventEnd(logEvent)
- return
-
-
- def integrateResidualAssembled(self, residual, t, fields):
- """
- Integrate contributions to residual term at time t that do not
- require assembly over cells, vertices, or processors.
- """
- logEvent = "%sresidualAs" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- self.cppHandle.integrateResidualAssembled(residual, t,
- fields.cppHandle,
- self.mesh.cppHandle,
- self.mesh.coordsys.cppHandle)
- self._logger.eventEnd(logEvent)
- return
-
-
- def integrateJacobianAssembled(self, jacobian, t, fields):
- """
- Integrate contributions to Jacobian term at time t that do not
- require assembly over cells, vertices, or processors.
- """
- logEvent = "%sjacobianAs" % self._loggingPrefix
- self._logger.eventBegin(logEvent)
-
- assert(None != self.cppHandle)
- self.cppHandle.integrateJacobianAssembled(jacobian, t,
- fields.cppHandle,
- self.mesh.cppHandle)
- self._logger.eventEnd(logEvent)
- return
-
-
def poststep(self, t, dt, totalTime, fields):
"""
Hook for doing stuff after advancing time step.
@@ -222,8 +91,7 @@
logEvent = "%spoststep" % self._loggingPrefix
self._logger.eventBegin(logEvent)
- assert(None != self.cppHandle)
- self.cppHandle.updateState(t, fields.cppHandle, self.mesh.cppHandle)
+ self.updateState(t, fields)
self._logger.eventEnd(logEvent)
return
@@ -250,15 +118,8 @@
logger.setClassName("FE Integrator")
logger.initialize()
- events = ["verify",
+ events = ["preinit",
"init",
- "timestep",
- "solnIncr",
- "newJacobian",
- "residual",
- "jacobian",
- "residualAs",
- "jacobianAs",
"poststep",
"finalize"]
for event in events:
Modified: short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -33,8 +33,6 @@
Constructor.
"""
Integrator.__init__(self)
- import journal
- self._info = journal.info(name)
self.output = None
self.availableFields = None
self.name = "Integrator Elasticity"
@@ -45,21 +43,16 @@
"""
Setup integrator.
"""
- Integrator.preinitialize(self, mesh)
-
- assert(None != self.cppHandle)
self.mesh = mesh
-
- material.preinitialize()
-
- self.quadrature = material.quadrature
- self.cppHandle.quadrature = self.quadrature.cppHandle
-
- self.material = material
- self.cppHandle.material = self.material.cppHandle
self.output = material.output
self.availableFields = material.availableFields
+ self.quadrature(material.quadrature)
+ self.material(material)
+
+ Integrator.preinitialize(self, mesh)
+ material.preinitialize()
self.output.preinitialize(self)
+
return
@@ -71,10 +64,8 @@
self._logger.eventBegin(logEvent)
Integrator.verifyConfiguration(self)
- self.material.verifyConfiguration()
self.output.verifyConfiguration(self.mesh)
-
self._logger.eventEnd(logEvent)
return
@@ -91,7 +82,6 @@
Integrator.initialize(self, totalTime, numTimeSteps, normalizer)
- self.material.initialize(self.mesh, totalTime, numTimeSteps, normalizer)
self.output.initialize(normalizer, self.quadrature)
self.output.writeInfo()
self.output.open(totalTime, numTimeSteps)
@@ -128,10 +118,9 @@
Get cell field.
"""
if None == fields:
- (field, fieldType) = self.cppHandle.cellField(name, self.mesh.cppHandle)
+ (field, fieldType) = self.cellField(name)
else:
- (field, fieldType) = self.cppHandle.cellField(name, self.mesh.cppHandle,
- fields.cppHandle)
+ (field, fieldType) = self.cellField(name, fields)
return (field, fieldType)
Modified: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/Makefile.am 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/Makefile.am 2009-03-12 19:40:31 UTC (rev 14308)
@@ -24,7 +24,7 @@
TestFIATLagrange.py \
TestFIATSimplex.py \
TestCellGeometry.py \
- TestQuadrature.py \
+ TestMeshQuadrature.py \
TestIntegrator.py \
TestElasticityExplicit.py \
TestElasticityImplicit.py
Modified: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityImplicit.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityImplicit.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityImplicit.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -39,224 +39,61 @@
"""
Test preiniitlaize().
"""
- from spatialdata.units.Nondimensional import Nondimensional
- normalizer = Nondimensional()
- normalizer.initialize()
-
- # Setup mesh
- cs = CSCart()
- cs.spaceDim = 2
- from pylith.meshio.MeshIOAscii import MeshIOAscii
- importer = MeshIOAscii()
- importer.filename = "data/tri3.mesh"
- importer.coordsys = cs
- mesh = importer.read(normalizer, debug=False, interpolate=False)
-
- # Setup material
- from pylith.feassemble.FIATSimplex import FIATSimplex
- cell = FIATSimplex()
- cell.shape = "triangle"
- cell.degree = 1
- cell.order = 1
- from pylith.feassemble.quadrature.Quadrature2D import Quadrature2D
- quadrature = Quadrature2D()
- quadrature._configure()
- quadrature.cell = cell
- minJacobian = 4.0e-02;
- quadrature.minJacobian = minJacobian
-
- from spatialdata.spatialdb.SimpleDB import SimpleDB
- from spatialdata.spatialdb.SimpleIOAscii import SimpleIOAscii
- iohandler = SimpleIOAscii()
- iohandler.filename = "data/elasticplanestrain.spatialdb"
- db = SimpleDB()
- db.label = "elastic plane strain"
- db.iohandler = iohandler
- initialStateDB = None
-
- from pylith.materials.ElasticPlaneStrain import ElasticPlaneStrain
- material = ElasticPlaneStrain()
- material.id = 0
- material.label = "elastic plane strain"
- material.db = db
- material.initialStateDB = initialStateDB
- material.quadrature = quadrature
- from pylith.meshio.OutputMatElastic import OutputMatElastic
- material.output = OutputMatElastic()
- material.output._configure()
- material.output.writer._configure()
-
- integrator = ElasticityImplicit()
- integrator.preinitialize(mesh, material)
- self.assertEqual(mesh, integrator.mesh)
- self.assertEqual(minJacobian, integrator.quadrature.minJacobian)
+ self._preinitialize()
return
-
- def test_timeStep(self):
- """
- Test timeStep().
- """
- dt = 2.3
- (mesh, integrator, fields) = self._initialize()
- integrator.timeStep(dt)
- return
-
- def test_stableTimeStep(self):
- """
- Test stableTimeStep().
- """
- (mesh, integrator, fields) = self._initialize()
-
- self.assertEqual(1.0e+30, integrator.stableTimeStep())
- return
-
-
- def test_needNewJacobian(self):
- """
- Test needNewJacobian().
- """
- (mesh, integrator, fields) = self._initialize()
- self.assertEqual(True, integrator.needNewJacobian())
- return
-
-
- def test_useSolnIncr(self):
- """
- Test useSolnIncr().
- """
- (mesh, integrator, fields) = self._initialize()
- integrator.useSolnIncr(True)
- return
-
-
- def test_integrateResidual(self):
- """
- Test integrateResidual().
-
- WARNING: This is not a rigorous test of integrateResidual() because we
- neither set the input fields or verify the results.
- """
- (mesh, integrator, fields) = self._initialize()
-
- residual = fields.getReal("residual")
- t = 3.4
- integrator.integrateResidual(residual, t, fields)
-
- # We should really add something here to check to make sure things
- # actually initialized correctly
- return
-
-
- def test_integrateJacobian(self):
- """
- Test integrateJacobian().
-
- WARNING: This is not a rigorous test of integrateJacobian() because we
- neither set the input fields or verify the results.
- """
- (mesh, integrator, fields) = self._initialize()
-
- jacobian = mesh.createMatrix(fields.getReal("residual"))
- import pylith.utils.petsc as petsc
- petsc.mat_setzero(jacobian)
- t = 7.3
- integrator.integrateJacobian(jacobian, t, fields)
- self.assertEqual(False, integrator.needNewJacobian())
-
- # We should really add something here to check to make sure things
- # actually initialized correctly
- return
-
-
- def test_poststep(self):
- """
- Test poststep().
-
- WARNING: This is not a rigorous test of poststep() because we
- neither set the input fields or verify the results.
- """
- (mesh, integrator, fields) = self._initialize()
-
- t = 0.27
-
- residual = fields.getReal("residual")
- integrator.integrateResidual(residual, t, fields)
-
- dt = 0.01
- totalTime = 10
- integrator.poststep(t, dt, totalTime, fields)
-
- # We should really add something here to check to make sure things
- # actually initialized correctly
- return
-
-
- def test_finalize(self):
- """
- Test finalize().
-
- WARNING: This is not a rigorous test of finalize() because we
- neither set the input fields or verify the results.
- """
- (mesh, integrator, fields) = self._initialize()
-
- integrator.finalize()
-
- # We should really add something here to check to make sure things
- # actually initialized correctly.
- return
-
-
# PRIVATE METHODS ////////////////////////////////////////////////////
- def _initialize(self):
+ def _preinitialize(self):
"""
- Initialize integrator.
+ Setup mesh and integrator and preinitialize integrator.
"""
- dt = 2.3
-
from spatialdata.units.Nondimensional import Nondimensional
normalizer = Nondimensional()
- normalizer.initialize()
+ normalizer._configure()
# Setup mesh
cs = CSCart()
- cs.spaceDim = 2
+ cs.inventory.spaceDim = 2
+ cs._configure()
from pylith.meshio.MeshIOAscii import MeshIOAscii
importer = MeshIOAscii()
- importer.filename = "data/tri3.mesh"
- importer.coordsys = cs
+ importer.inventory.filename = "data/tri3.mesh"
+ importer.inventory.coordsys = cs
+ importer._configure()
mesh = importer.read(normalizer, debug=False, interpolate=False)
# Setup material
from pylith.feassemble.FIATSimplex import FIATSimplex
cell = FIATSimplex()
- cell.shape = "triangle"
- cell.degree = 1
- cell.order = 1
- from pylith.feassemble.quadrature.Quadrature2D import Quadrature2D
- quadrature = Quadrature2D()
+ cell.inventory.shape = "triangle"
+ cell.inventory.degree = 1
+ cell.inventory.order = 1
+ cell._configure()
+ from pylith.feassemble.Quadrature import MeshQuadrature
+ quadrature = MeshQuadrature()
+ quadrature.inventory.cell = cell
quadrature._configure()
- quadrature.cell = cell
from spatialdata.spatialdb.SimpleDB import SimpleDB
from spatialdata.spatialdb.SimpleIOAscii import SimpleIOAscii
iohandler = SimpleIOAscii()
- iohandler.filename = "data/elasticplanestrain.spatialdb"
+ iohandler.inventory.filename = "data/elasticplanestrain.spatialdb"
+ iohandler._configure()
db = SimpleDB()
- db.label = "elastic plane strain"
- db.iohandler = iohandler
- initialStateDB = None
+ db.inventory.label = "elastic plane strain"
+ db.inventory.iohandler = iohandler
+ db._configure()
from pylith.materials.ElasticPlaneStrain import ElasticPlaneStrain
material = ElasticPlaneStrain()
- material.id = 0
- material.label = "elastic plane strain"
- material.db = db
- material.initialStateDB = initialStateDB
- material.quadrature = quadrature
+ material.inventory.label = "elastic plane strain"
+ material.inventory.id = 0
+ material.inventory.db = db
+ material.inventory.quadrature = quadrature
+ material._configure()
+
from pylith.meshio.OutputMatElastic import OutputMatElastic
material.output = OutputMatElastic()
material.output._configure()
@@ -265,24 +102,33 @@
# Setup integrator
integrator = ElasticityImplicit()
integrator.preinitialize(mesh, material)
+ return (mesh, integrator)
+
+
+ def _initialize(self, mesh, integrator):
+ """
+ Initialize integrator.
+ """
+ dt = 2.3
+
from pyre.units.time import s
integrator.initialize(totalTime=0.0*s, numTimeSteps=1, normalizer=normalizer)
integrator.timeStep(dt)
# Setup fields
- from pylith.topology.FieldsManager import FieldsManager
- fields = FieldsManager(mesh)
- fields.addReal("residual")
- fields.addReal("dispTBctpdt")
- fields.solutionField("dispTBctpdt")
- fields.setFiberDimension("residual", cs.spaceDim)
- fields.allocate("residual")
+ from pylith.topology.SolutionFields import SolutionFields
+ fields = SolutionFields(mesh)
+ fields.add("residual")
+ fields.add("disp(t), bc(t+dt)")
+ fields.solutionName("disp(t), bc(t+dt)")
+
+ residual = fields.get("residual")
+ residual.newSection(residual.VERTICES_FIELD, mesh.coordsys().spaceDim())
+ residual.allocate()
fields.copyLayout("residual")
- import pylith.topology.topology as bindings
- bindings.zeroRealSection(fields.getReal("residual"))
-
- return (mesh, integrator, fields)
+ residual.zero()
+ return
# End of file
Deleted: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestIntegrator.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestIntegrator.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestIntegrator.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -1,43 +0,0 @@
-#!/usr/bin/env python
-#
-# ======================================================================
-#
-# Brad T. Aagaard
-# U.S. Geological Survey
-#
-# {LicenseText}
-#
-# ======================================================================
-#
-
-## @file unittests/pytests/feassemble/TestIntegrator.py
-
-## @brief Unit testing of Python Integrator object.
-
-import unittest
-from pylith.feassemble.Integrator import Integrator
-
-# ----------------------------------------------------------------------
-class TestIntegrator(unittest.TestCase):
- """
- Unit testing of Python Integrator object.
- """
-
- def test_constructor(self):
- """
- Test constructor.
- """
- i = Integrator()
- return
-
-
- def test_finalize(self):
- """
- Test constructor.
- """
- i = Integrator()
- i.finalize()
- return
-
-
-# End of file
Modified: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestMeshQuadrature.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestMeshQuadrature.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestMeshQuadrature.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -10,9 +10,9 @@
# ======================================================================
#
-## @file unittests/pytests/feassemble/TestQuadrature.py
+## @file unittests/pytests/feassemble/TestMeshQuadrature.py
-## @brief Unit testing of Python Quadrature object.
+## @brief Unit testing of Python MeshQuadrature object.
import unittest
import numpy
@@ -43,7 +43,7 @@
# ----------------------------------------------------------------------
class TestMeshQuadrature(unittest.TestCase):
"""
- Unit testing of Python Quadrature object.
+ Unit testing of Python MeshQuadrature object.
"""
def test_minJacobian(self):
@@ -118,24 +118,16 @@
self.assertEqual(2, quadrature.numQuadPts())
from pylith.utils.testarray import test_double
- self.failIf(True)
- #import pylith.feassemble.testfeassemble as testmodule
+ from pylith.utils.utils import TestArray_checkDouble
- #vertices = testmodule.vertices(quadrature.cppHandle)
- #test_double(self, verticesE, numpy.array(vertices))
-
- #basis = testmodule.basis(quadrature.cppHandle)
- #test_double(self, basisE, numpy.array(basis))
-
- #basisDeriv = testmodule.basisDeriv(quadrature.cppHandle)
- #test_double(self, basisDerivE, numpy.array(basisDeriv))
-
- #quadWts = testmodule.quadWts(quadrature.cppHandle)
- #test_double(self, quadWtsE, numpy.array(quadWts))
-
- #quadPts = testmodule.quadPtsRef(quadrature.cppHandle)
- #test_double(self, quadPtsE, numpy.array(quadPts))
-
+ self.assertTrue(TestArray_checkDouble(basisE.ravel(),
+ quadrature.basis()))
+ self.assertTrue(TestArray_checkDouble(basisDerivE.ravel(),
+ quadrature.basisDerivRef()))
+ self.assertTrue(TestArray_checkDouble(quadPtsE.ravel(),
+ quadrature.quadPtsRef()))
+ self.assertTrue(TestArray_checkDouble(quadWtsE.ravel(),
+ quadrature.quadWts()))
return
Deleted: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.hh
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.hh 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.hh 2009-03-12 19:40:31 UTC (rev 14308)
@@ -1,101 +0,0 @@
-// -*- C++ -*-
-//
-// ----------------------------------------------------------------------
-//
-// Brad T. Aagaard
-// U.S. Geological Survey
-//
-// {LicenseText}
-//
-// ----------------------------------------------------------------------
-//
-
-/**
- * @file unittests/pytests/feassemble/TestQuadrature.hh
- *
- * @brief C++ TestQuadrature object
- *
- * Helper class for unit testing of Python Quadrature.
- */
-
-#if !defined(pylith_feassemble_pytestquadrature_hh)
-#define pylith_feassemble_pytestquadrature_hh
-
-#include "pylith/feassemble/Quadrature.hh"
-#include "pylith/feassemble/CellGeometry.hh"
-
-/// Namespace for spatialdata package
-namespace pylith {
- namespace feassemble {
- class TestQuadrature;
- } // feassemble
-} // pylith
-
-/// Helper class for unit testing of Python Quadrature
-class pylith::feassemble::TestQuadrature
-{ // class TestQuadrature
-
- // PUBLIC METHODS /////////////////////////////////////////////////////
-public :
-
- /** Get number of dimensions in reference cell.
- *
- * @returns Number of dimensions
- */
- static int cellDim(const Quadrature& q);
-
- /** Get number of vertices in cell.
- *
- * @returns Number of vertices
- */
- static int numBasis(const Quadrature& q);
-
- /** Get number of quadrature points.
- *
- * @returns Number of points
- */
- static int numQuadPts(const Quadrature& q);
-
- /** Get number of dimensions in coordinates of cell vertices.
- *
- * @returns Number of dimensions
- */
- static int spaceDim(const Quadrature& q);
-
- /** Get vertices in reference cell.
- *
- * @returns Array of coordinates of vertices in reference cell.
- */
- static const double* vertices(const Quadrature& q);
-
- /** Get basis functions evaluated at quadrature points.
- *
- * @returns Array of basis functions evaluated at quadrature points
- */
- static const double* basis(const Quadrature& q);
-
- /** Get derivatives of basis functions evaluated at quadrature points.
- *
- * @returns Array of derivatives of basis fns evaluated at quad pts
- */
- static const double* basisDeriv(const Quadrature& q);
-
- /** Get coordinates of quadrature points in reference cell.
- *
- * @returns Array of coordinates of quadrature points
- */
- static const double* quadPtsRef(const Quadrature& q);
-
- /** Get weights of quadrature points.
- *
- * @returns Array of weights of quadrature points
- */
- static const double* quadWts(const Quadrature& q);
-
-}; // class TestQuadrature
-
-#include "TestQuadrature.icc" // inline methods
-
-#endif // pylith_feassemble_pytestquadrature_hh
-
-// End of file
Deleted: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.icc
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.icc 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestQuadrature.icc 2009-03-12 19:40:31 UTC (rev 14308)
@@ -1,82 +0,0 @@
-// -*- C++ -*-
-//
-// ======================================================================
-//
-// Brad T. Aagaard
-// U.S. Geological Survey
-//
-// {LicenseText}
-//
-// ======================================================================
-//
-
-#if !defined(pylith_feassemble_pytestquadrature_hh)
-#error "TestQuadrature.icc must be included only from TestQuadrature.hh"
-#else
-
-// Get number of dimensions in reference cell.
-inline
-int
-pylith::feassemble::TestQuadrature::cellDim(const Quadrature& q) {
- return q._cellDim;
-}
-
-// Get number of vertices in cell.
-inline
-int
-pylith::feassemble::TestQuadrature::numBasis(const Quadrature& q) {
- return q._numBasis;
-}
-
-// Get number of quadrature points.
-inline
-int
-pylith::feassemble::TestQuadrature::numQuadPts(const Quadrature& q) {
- return q._numQuadPts;
-}
-
-// Get number of dimensions in coordinates of cell vertices.
-inline
-int
-pylith::feassemble::TestQuadrature::spaceDim(const Quadrature& q) {
- return q._spaceDim;
-}
-
-// Get coordinates of vertices in reference cell.
-inline
-const double*
-pylith::feassemble::TestQuadrature::vertices(const Quadrature& q) {
- return &q._geometry->vertices()[0];
-}
-
-// Get basis functions evaluated at quadrature points.
-inline
-const double*
-pylith::feassemble::TestQuadrature::basis(const Quadrature& q) {
- return &q._basis[0];
-}
-
-// Get derivatives of basis functions evaluated at quadrature points.
-inline
-const double*
-pylith::feassemble::TestQuadrature::basisDeriv(const Quadrature& q) {
- return &q._basisDerivRef[0];
-}
-
-// Get coordinates of quadrature points
-inline
-const double*
-pylith::feassemble::TestQuadrature::quadPtsRef(const Quadrature& q) {
- return &q._quadPtsRef[0];
-}
-
-// Get weights of quadrature points
-inline
-const double*
-pylith::feassemble::TestQuadrature::quadWts(const Quadrature& q) {
- return &q._quadWts[0];
-}
-
-#endif
-
-// End of file
Added: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestSubMeshQuadrature.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestSubMeshQuadrature.py (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestSubMeshQuadrature.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -0,0 +1,218 @@
+#!/usr/bin/env python
+#
+# ======================================================================
+#
+# Brad T. Aagaard
+# U.S. Geological Survey
+#
+# {LicenseText}
+#
+# ======================================================================
+#
+
+## @file unittests/pytests/feassemble/TestSubMeshQuadrature.py
+
+## @brief Unit testing of Python SubMeshQuadrature object.
+
+import unittest
+import numpy
+
+from pylith.feassemble.Quadrature import SubMeshQuadrature
+from pylith.feassemble.FIATSimplex import FIATSimplex
+from pylith.feassemble.FIATLagrange import FIATLagrange
+
+# ----------------------------------------------------------------------
+def N0(p):
+ return -0.5*p*(1.0-p)
+
+def N0p(p):
+ return -0.5*(1.0-p) + 0.5*p
+
+def N1(p):
+ return 0.5*p*(1.0+p)
+
+def N1p(p):
+ return +0.5*(1.0+p) + 0.5*p
+
+def N2(p):
+ return (1.0-p**2)
+
+def N2p(p):
+ return -2.0*p
+
+# ----------------------------------------------------------------------
+class TestSubMeshQuadrature(unittest.TestCase):
+ """
+ Unit testing of Python SubMeshQuadrature object.
+ """
+
+ def test_minJacobian(self):
+ """
+ Test minJacobian().
+ """
+ minJacobian = 4.0e-02;
+ q = SubMeshQuadrature()
+ q.minJacobian(minJacobian)
+ self.assertEqual(minJacobian, q.minJacobian())
+ return
+
+
+ def test_checkConditioning(self):
+ """
+ Test checkConditioning().
+ """
+ q = SubMeshQuadrature()
+
+ flag = False # default
+ self.assertEqual(flag, q.checkConditioning())
+
+ flag = True
+ q.checkConditioning(flag)
+ self.assertEqual(flag, q.checkConditioning())
+
+ flag = False
+ q.checkConditioning(flag)
+ self.assertEqual(flag, q.checkConditioning())
+
+ return
+
+
+ def test_initialize(self):
+ """
+ Test initialize().
+ """
+ cell = FIATSimplex()
+ cell.inventory.shape = "line"
+ cell.inventory.degree = 2
+ cell.inventory.order = 2
+ cell._configure()
+
+ verticesE = numpy.array([ [-1.0], [1.0], [0.0] ])
+ quadPtsE = numpy.array( [[-1.0/3**0.5],
+ [+1.0/3**0.5]],
+ dtype=numpy.float64 )
+ quadWtsE = numpy.array( [1.0, 1.0], dtype=numpy.float64 )
+
+ # Compute basis functions and derivatives at quadrature points
+ basisE = numpy.zeros( (2, 3), dtype=numpy.float64)
+ basisDerivE = numpy.zeros( (2, 3, 1), dtype=numpy.float64)
+ iQuad = 0
+ for q in quadPtsE:
+ basisE[iQuad] = numpy.array([N0(q), N1(q), N2(q)],
+ dtype=numpy.float64).reshape( (3,) )
+ deriv = numpy.array([[N0p(q)], [N1p(q)], [N2p(q)]],
+ dtype=numpy.float64)
+ basisDerivE[iQuad] = deriv.reshape((3, 1))
+ iQuad += 1
+
+ quadrature = SubMeshQuadrature()
+ quadrature.inventory.cell = cell
+ quadrature._configure()
+
+ quadrature.preinitialize(spaceDim=2)
+ quadrature.initialize()
+
+ self.assertEqual(1, quadrature.cellDim())
+ self.assertEqual(2, quadrature.spaceDim())
+ self.assertEqual(3, quadrature.numBasis())
+ self.assertEqual(2, quadrature.numQuadPts())
+
+ from pylith.utils.testarray import test_double
+ from pylith.utils.utils import TestArray_checkDouble
+
+ self.assertTrue(TestArray_checkDouble(basisE.ravel(),
+ quadrature.basis()))
+ self.assertTrue(TestArray_checkDouble(basisDerivE.ravel(),
+ quadrature.basisDerivRef()))
+ self.assertTrue(TestArray_checkDouble(quadPtsE.ravel(),
+ quadrature.quadPtsRef()))
+ self.assertTrue(TestArray_checkDouble(quadWtsE.ravel(),
+ quadrature.quadWts()))
+ return
+
+
+ def test_simplex1D(self):
+ """
+ Test setup of quadrature for simplex cells for a 1-D problem in 2-D space.
+ """
+ spaceDim = 2
+
+ cell = FIATSimplex()
+ cell.inventory.shape = "line"
+ cell._configure()
+
+ quadrature = SubMeshQuadrature()
+ quadrature.inventory.cell = cell
+ quadrature._configure()
+
+ quadrature.preinitialize(spaceDim)
+ self.assertEqual(1, quadrature.cellDim())
+ self.assertEqual(spaceDim, quadrature.spaceDim())
+ self.assertEqual(2, quadrature.numBasis())
+ return
+
+
+ def test_simplex2D(self):
+ """
+ Test setup of quadrature for simplex cells for a 2-D problem in 3-D space.
+ """
+ spaceDim = 3
+
+ cell = FIATSimplex()
+ cell.inventory.shape = "triangle"
+ cell._configure()
+
+ quadrature = SubMeshQuadrature()
+ quadrature.inventory.cell = cell
+ quadrature._configure()
+
+ quadrature.preinitialize(spaceDim)
+ self.assertEqual(2, quadrature.cellDim())
+ self.assertEqual(spaceDim, quadrature.spaceDim())
+ self.assertEqual(3, quadrature.numBasis())
+ return
+
+
+ def test_lagrange1D(self):
+ """
+ Test setup of quadrature for Lagrange cells for a 1-D problem in 2-D space.
+ """
+ spaceDim = 2
+
+ cell = FIATLagrange()
+ cell.inventory.dimension = 1
+ cell._configure()
+
+ quadrature = SubMeshQuadrature()
+ quadrature.inventory.cell = cell
+ quadrature._configure()
+
+ quadrature.preinitialize(spaceDim)
+ self.assertEqual(1, quadrature.cellDim())
+ self.assertEqual(spaceDim, quadrature.spaceDim())
+ self.assertEqual(2, quadrature.numBasis())
+ return
+
+
+ def test_lagrange2D(self):
+ """
+ Test setup of quadrature for Lagrange cells for a 2-D problem in 3-D space.
+ """
+ spaceDim = 3
+
+ cell = FIATLagrange()
+ cell.inventory.dimension = 2
+ cell._configure()
+
+ quadrature = SubMeshQuadrature()
+ quadrature.inventory.cell = cell
+ quadrature._configure()
+
+ quadrature.preinitialize(spaceDim)
+ self.assertEqual(2, quadrature.cellDim())
+ self.assertEqual(spaceDim, quadrature.spaceDim())
+ self.assertEqual(4, quadrature.numBasis())
+ return
+
+
+# End of file
Modified: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/testfeassemble.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/testfeassemble.py 2009-03-12 18:02:14 UTC (rev 14307)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/testfeassemble.py 2009-03-12 19:40:31 UTC (rev 14308)
@@ -66,15 +66,15 @@
from TestMeshQuadrature import TestMeshQuadrature
suite.addTest(unittest.makeSuite(TestMeshQuadrature))
- #from TestIntegrator import TestIntegrator
- #suite.addTest(unittest.makeSuite(TestIntegrator))
+ from TestSubMeshQuadrature import TestSubMeshQuadrature
+ suite.addTest(unittest.makeSuite(TestSubMeshQuadrature))
+ from TestElasticityImplicit import TestElasticityImplicit
+ suite.addTest(unittest.makeSuite(TestElasticityImplicit))
+
#from TestElasticityExplicit import TestElasticityExplicit
#suite.addTest(unittest.makeSuite(TestElasticityExplicit))
- #from TestElasticityImplicit import TestElasticityImplicit
- #suite.addTest(unittest.makeSuite(TestElasticityImplicit))
-
return suite
More information about the CIG-COMMITS
mailing list