[cig-commits] r14342 - in short/3D/PyLith/branches/pylith-swig: modulesrc/feassemble pylith/feassemble unittests/pytests/feassemble

brad at geodynamics.org brad at geodynamics.org
Mon Mar 16 12:41:55 PDT 2009


Author: brad
Date: 2009-03-16 12:41:55 -0700 (Mon, 16 Mar 2009)
New Revision: 14342

Added:
   short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityExplicit.i
Modified:
   short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/Makefile.am
   short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i
   short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityExplicit.py
   short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py
   short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py
   short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityExplicit.py
   short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityImplicit.py
   short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/testfeassemble.py
Log:
Worked on Python tests of elasticity integrators.

Added: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityExplicit.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityExplicit.i	                        (rev 0)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/ElasticityExplicit.i	2009-03-16 19:41:55 UTC (rev 14342)
@@ -0,0 +1,73 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+//                           Brad T. Aagaard
+//                        U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ----------------------------------------------------------------------
+//
+
+/** @file modulesrc/feassemble/ElasticityExplicit.i
+ *
+ * @brief Python interface to C++ ElasticityExplicit object.
+ */
+
+namespace pylith {
+  namespace feassemble {
+
+    class ElasticityExplicit : public IntegratorElasticity
+    { // ElasticityExplicit
+
+      // PUBLIC MEMBERS /////////////////////////////////////////////////
+    public :
+      
+      /// Constructor
+      ElasticityExplicit(void);
+      
+      /// Destructor
+      ~ElasticityExplicit(void);
+      
+      /** Set time step for advancing from time t to time t+dt.
+       *
+       * @param dt Time step
+       */
+      void timeStep(const double dt);
+      
+      /** 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 contributions to residual term (r) for operator.
+       *
+       * @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(pylith::topology::Jacobian* jacobian,
+			     const double t,
+			     pylith::topology::SolutionFields* const fields);
+
+    }; // ElasticityExplicit
+
+  } // 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-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/Makefile.am	2009-03-16 19:41:55 UTC (rev 14342)
@@ -31,7 +31,9 @@
 	QuadratureRefCell.i \
 	Quadrature.i \
 	Integrator.i \
-	IntegratorElasticity.i
+	IntegratorElasticity.i \
+	ElasticityImplicit.i \
+	ElasticityExplicit.i
 
 swig_generated = \
 	feassemble_wrap.cxx \

Modified: short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i
===================================================================
--- short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i	2009-03-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/modulesrc/feassemble/feassemble.i	2009-03-16 19:41:55 UTC (rev 14342)
@@ -78,6 +78,7 @@
 %include "Integrator.i"
 %include "IntegratorElasticity.i"
 %include "ElasticityImplicit.i"
+%include "ElasticityExplicit.i"
 
 // Template instatiation
 %template(MeshQuadrature) pylith::feassemble::Quadrature<pylith::topology::Mesh>;

Modified: short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityExplicit.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityExplicit.py	2009-03-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityExplicit.py	2009-03-16 19:41:55 UTC (rev 14342)
@@ -18,9 +18,10 @@
 ## Factory: integrator
 
 from IntegratorElasticity import IntegratorElasticity
+from feassemble import ElasticityExplicit as ModuleElasticityExplicit
 
 # ElasticityExplicit class
-class ElasticityExplicit(IntegratorElasticity):
+class ElasticityExplicit(IntegratorElasticity, ModuleElasticityExplicit):
   """
   Python object for explicit time integration of dynamic elasticity
   equation using finite-elements.
@@ -33,10 +34,22 @@
     Constructor.
     """
     IntegratorElasticity.__init__(self, name)
+    ModuleElasticityExplicit.__init__(self)
     self._loggingPrefix = "ElEx "
+    return
 
-    import pylith.feassemble.feassemble as bindings
-    self.cppHandle = bindings.ElasticityExplicit()
+
+  def initialize(self, totalTime, numTimeSteps, normalizer):
+    """
+    Do initialization.
+    """
+    logEvent = "%sinit" % self._loggingPrefix
+    self._logger.eventBegin(logEvent)
+
+    IntegratorElasticity.initialize(self, totalTime, numTimeSteps, normalizer)
+    ModuleElasticityExplicit.initialize(self, self.mesh)
+    
+    self._logger.eventEnd(logEvent)
     return
 
 

Modified: short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py	2009-03-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/pylith/feassemble/ElasticityImplicit.py	2009-03-16 19:41:55 UTC (rev 14342)
@@ -46,11 +46,8 @@
     logEvent = "%sinit" % self._loggingPrefix
     self._logger.eventBegin(logEvent)
 
-    print "AAA"
     IntegratorElasticity.initialize(self, totalTime, numTimeSteps, normalizer)
-    print "BBB"
     ModuleElasticityImplicit.initialize(self, self.mesh)
-    print "CCC"
     
     self._logger.eventEnd(logEvent)
     return

Modified: short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py	2009-03-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/pylith/feassemble/IntegratorElasticity.py	2009-03-16 19:41:55 UTC (rev 14342)
@@ -86,9 +86,7 @@
     self._info.log("Initializing integrator for material '%s'." % \
                    self.materialObj.label)
 
-    print "DDD"
     Integrator.initialize(self, totalTime, numTimeSteps, normalizer)
-    print "EEE"
 
     #self.output.initialize(normalizer, self.materialObj.quadrature)
     #self.output.writeInfo()

Modified: short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityExplicit.py
===================================================================
--- short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityExplicit.py	2009-03-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityExplicit.py	2009-03-16 19:41:55 UTC (rev 14342)
@@ -38,100 +38,99 @@
   def test_preinitialize(self):
     """
     Test preiniitlaize().
+
+    WARNING: This is not a rigorous test of preinitialize() because we
+    neither set the input fields or verify the results.
     """
-    from spatialdata.units.Nondimensional import Nondimensional
-    normalizer = Nondimensional()
-    normalizer.initialize()
+    (mesh, integrator) = self._preinitialize()
 
-    # 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)
+    # No test of result.
+    return
 
-    # 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.quadrature = quadrature
-    material.initialStateDB = initialStateDB
-    from pylith.meshio.OutputMatElastic import OutputMatElastic
-    material.output = OutputMatElastic()
-    material.output._configure()
-    material.output.writer._configure()
+  def test_verifyConfiguration(self):
+    """
+    Test verifyConfiguration().
 
-    integrator = ElasticityExplicit()
-    integrator.preinitialize(mesh, material)
-    self.assertEqual(mesh, integrator.mesh)
-    self.assertEqual(minJacobian, integrator.quadrature.minJacobian)
+    WARNING: This is not a rigorous test of verifyConfiguration()
+    because we neither set the input fields or verify the results.
+    """
+    (mesh, integrator) = self._preinitialize()
+    integrator.verifyConfiguration()
+
+    # No test of result.
     return
-    
 
+
+  def test_initialize(self):
+    """
+    Test initialize().
+
+    WARNING: This is not a rigorous test of initialize() because we
+    neither set the input fields or verify the results.
+    """
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
+
+    # No test of result.
+    return
+
+
   def test_timeStep(self):
     """
     Test timeStep().
+
+    WARNING: This is not a rigorous test of timeStep() because we
+    neither set the input fields or verify the results.
     """
     dt = 2.3
-    (mesh, integrator, fields) = self._initialize()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
     integrator.timeStep(dt)
+
+    # No test of result.
     return
 
-  
+
   def test_stableTimeStep(self):
     """
     Test stableTimeStep().
     """
-    (mesh, integrator, fields) = self._initialize()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
 
     self.assertEqual(1.0e+30, integrator.stableTimeStep())
     return
 
-  
+
   def test_needNewJacobian(self):
     """
     Test needNewJacobian().
     """
-    (mesh, integrator, fields) = self._initialize()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
+
     self.assertEqual(True, integrator.needNewJacobian())
     return
 
-  
+
   def test_useSolnIncr(self):
     """
     Test useSolnIncr().
+
+    WARNING: This is not a rigorous test of useSolnIncr() because we
+    neither set the input fields or verify the results.
     """
-    (mesh, integrator, fields) = self._initialize()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
+
     try:
       integrator.useSolnIncr(True)
       self.failIf(True)
     except:
       self.failIf(False)
+
+    # No test of result.
     return
 
 
@@ -139,42 +138,42 @@
     """
     Test integrateResidual().
 
-    WARNING: This is not a rigorous test of integrateResidual() because we
-    neither set the input fields or verify the results.
+    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()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
 
-    residual = fields.getReal("residual")
-
-    t = 0.45
+    residual = fields.get("residual")
+    t = 3.4
     integrator.integrateResidual(residual, t, fields)
 
-    # We should really add something here to check to make sure things
-    # actually initialized correctly    
+    # No test of result.
     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.
+    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()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
 
-    jacobian = mesh.createMatrix(fields.getReal("residual"))
-    import pylith.utils.petsc as petsc
-    petsc.mat_setzero(jacobian)
-    t = 0.145
+    from pylith.topology.Jacobian import Jacobian
+    jacobian = Jacobian(fields)
+    jacobian.zero()
+    t = 7.3
+    self.assertEqual(True, integrator.needNewJacobian())
     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    
+    
+    # No test of result.
     return
 
-  
+
   def test_poststep(self):
     """
     Test poststep().
@@ -182,86 +181,69 @@
     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()
+    (mesh, integrator) = self._preinitialize()
+    fields = self._initialize(mesh, integrator)
 
-    t = 3.45
-
-    residual = fields.getReal("residual")
-    integrator.integrateResidual(residual, t, fields)
-
-    dt = 0.02
-    totalTime = 5.0
+    t = 7.3
+    dt = 0.1
+    totalTime = 23.0
     integrator.poststep(t, dt, totalTime, fields)
 
-    # We should really add something here to check to make sure things
-    # actually initialized correctly    
+    # No test of result
     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.dbProperties = db
+    material.inventory.quadrature = quadrature
+    material._configure()
+    
     from pylith.meshio.OutputMatElastic import OutputMatElastic
     material.output = OutputMatElastic()
     material.output._configure()
@@ -270,27 +252,40 @@
     # Setup integrator
     integrator = ElasticityExplicit()
     integrator.preinitialize(mesh, material)
+    return (mesh, integrator)
+
+
+  def _initialize(self, mesh, integrator):
+    """
+    Initialize integrator.
+    """
+    dt = 2.3
+    
+    from spatialdata.units.Nondimensional import Nondimensional
+    normalizer = Nondimensional()
+    normalizer._configure()
+
     from pyre.units.time import s
-    integrator.initialize(totalTime=0.0*s, numTimeSteps=1, normalizer=normalizer)
+    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("solution")
-    fields.addReal("dispT")
-    fields.addReal("dispTmdt")
-    fields.createHistory(["solution", "dispT", "dispTmdt"])
-    fields.solutionField("solution")
-    fields.setFiberDimension("residual", cs.spaceDim)
-    fields.allocate("residual")
+    from pylith.topology.SolutionFields import SolutionFields
+    fields = SolutionFields(mesh)
+    fields.add("residual")
+    fields.add("disp(t+dt)")
+    fields.add("disp(t)")
+    fields.add("disp(t-dt)")
+    fields.solutionName("disp(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 fields
 
 
 # End of file 

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-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/TestElasticityImplicit.py	2009-03-16 19:41:55 UTC (rev 14342)
@@ -162,6 +162,7 @@
     jacobian = Jacobian(fields)
     jacobian.zero()
     t = 7.3
+    self.assertEqual(True, integrator.needNewJacobian())
     integrator.integrateJacobian(jacobian, t, fields)
     self.assertEqual(False, integrator.needNewJacobian())
     
@@ -261,11 +262,9 @@
     normalizer._configure()
 
     from pyre.units.time import s
-    print "AA"
-    integrator.initialize(totalTime=0.0*s, numTimeSteps=1, normalizer=normalizer)
-    print "BB"
+    integrator.initialize(totalTime=0.0*s, numTimeSteps=1,
+                          normalizer=normalizer)
     integrator.timeStep(dt)
-    print "CC"
 
     # Setup fields
     from pylith.topology.SolutionFields import SolutionFields

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-16 19:00:56 UTC (rev 14341)
+++ short/3D/PyLith/branches/pylith-swig/unittests/pytests/feassemble/testfeassemble.py	2009-03-16 19:41:55 UTC (rev 14342)
@@ -72,8 +72,8 @@
     from TestElasticityImplicit import TestElasticityImplicit
     suite.addTest(unittest.makeSuite(TestElasticityImplicit))
 
-    #from TestElasticityExplicit import TestElasticityExplicit
-    #suite.addTest(unittest.makeSuite(TestElasticityExplicit))
+    from TestElasticityExplicit import TestElasticityExplicit
+    suite.addTest(unittest.makeSuite(TestElasticityExplicit))
 
     return suite
 



More information about the CIG-COMMITS mailing list