[cig-commits] r6924 - in short/3D/PyLith/trunk: libsrc/feassemble modulesrc/feassemble pylith pylith/feassemble pylith/problems

brad at geodynamics.org brad at geodynamics.org
Fri May 18 11:38:04 PDT 2007


Author: brad
Date: 2007-05-18 11:37:57 -0700 (Fri, 18 May 2007)
New Revision: 6924

Added:
   short/3D/PyLith/trunk/pylith/feassemble/ImplicitElasticity.py
   short/3D/PyLith/trunk/pylith/feassemble/IntegratorImplicit.py
Modified:
   short/3D/PyLith/trunk/libsrc/feassemble/IntegratorImplicit.cc
   short/3D/PyLith/trunk/modulesrc/feassemble/feassemble.pyxe.src
   short/3D/PyLith/trunk/pylith/Makefile.am
   short/3D/PyLith/trunk/pylith/feassemble/__init__.py
   short/3D/PyLith/trunk/pylith/problems/Implicit.py
Log:
Fixed bugs in getting implicit problem to run. Need to debug.

Modified: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorImplicit.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorImplicit.cc	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorImplicit.cc	2007-05-18 18:37:57 UTC (rev 6924)
@@ -48,7 +48,7 @@
   if (_dt != -1.0)
     _dtm1 = _dt;
   else
-    _dtm1 = _dt;
+    _dtm1 = dt;
   _dt = dt;
   assert(_dt == _dtm1); // For now, don't allow variable time step
 } // timeStep

Modified: short/3D/PyLith/trunk/modulesrc/feassemble/feassemble.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/feassemble/feassemble.pyxe.src	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/modulesrc/feassemble/feassemble.pyxe.src	2007-05-18 18:37:57 UTC (rev 6924)
@@ -21,7 +21,9 @@
 
 #include "pylith/feassemble/Integrator.hh"
 #include "pylith/feassemble/IntegratorExplicit.hh"
+#include "pylith/feassemble/IntegratorImplicit.hh"
 #include "pylith/feassemble/ExplicitElasticity.hh"
+#include "pylith/feassemble/ImplicitElasticity.hh"
 
 #include "pylith/utils/petscfwd.h"
 
@@ -569,6 +571,123 @@
 
 
 # ----------------------------------------------------------------------
+cdef class IntegratorImplicit(Integrator):
+
+  def __init__(self):
+    """
+    Constructor.
+    """
+    Integrator.__init__(self)
+    return
+
+
+  def integrateResidual(self, fieldOut, fieldInT, mesh):
+    """
+    Integrate constant term (b) for quasi-static elasticity term for
+    finite elements.
+    """
+    # create shim for method 'integrateResidual'
+    #embed{ void IntegratorImplicit_integrateResidual(void* objVptr, void* fieldOutVptr, void* fieldInTVptr, void* meshVptr)
+    typedef ALE::Mesh Mesh;
+    typedef ALE::Mesh::real_section_type real_section_type;
+
+    try {
+      assert(0 != objVptr);
+      assert(0 != fieldOutVptr);
+      assert(0 != fieldInTVptr);
+      assert(0 != meshVptr);
+      ALE::Obj<Mesh>* mesh =
+        (ALE::Obj<Mesh>*) meshVptr;
+      ALE::Obj<real_section_type>* fieldOut =
+        (ALE::Obj<real_section_type>*) fieldOutVptr;
+      ALE::Obj<real_section_type>* fieldInT =
+        (ALE::Obj<real_section_type>*) fieldInTVptr;
+      ((pylith::feassemble::IntegratorImplicit*) objVptr)->integrateResidual(
+                *fieldOut, *fieldInT, *mesh);
+    } catch (const std::exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.what()));
+    } catch (const ALE::Exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.msg().c_str()));
+    } catch (...) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      "Caught unknown C++ exception.");
+    } // try/catch
+    #}embed
+    cdef void* fieldOutVptr
+    cdef void* fieldInTVptr
+    fieldOutVptr = PyCObject_AsVoidPtr(fieldOut)
+    fieldInTVptr = PyCObject_AsVoidPtr(fieldInT)
+    IntegratorImplicit_integrateResidual(self.thisptr, 
+                                         fieldOutVptr, fieldInTVptr,
+                                         ptrFromHandle(mesh))
+    return
+
+
+  def integrateJacobian(self, mat, fieldIn, mesh):
+    """
+    Compute matrix (A) associated with operator.
+    """
+    # create shim for method 'integrateJacobian'
+    #embed{ void IntegratorImplicit_integrateJacobian(void* objVptr, void* matVptr, void* fieldInVptr, void* meshVptr)
+    typedef ALE::Mesh Mesh;
+    typedef ALE::Mesh::real_section_type real_section_type;
+
+    try {
+      assert(0 != objVptr);
+      assert(0 != matVptr);
+      assert(0 != fieldInVptr);
+      assert(0 != meshVptr);
+      ALE::Obj<Mesh>* mesh =
+        (ALE::Obj<Mesh>*) meshVptr;
+      PetscMat* mat = (PetscMat*) matVptr;
+      ALE::Obj<real_section_type>* fieldIn =
+        (ALE::Obj<real_section_type>*) fieldInVptr;
+      ((pylith::feassemble::IntegratorImplicit*) objVptr)->integrateJacobian(
+                mat, *fieldIn, *mesh);
+    } catch (const std::exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.what()));
+    } catch (const ALE::Exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.msg().c_str()));
+    } catch (...) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      "Caught unknown C++ exception.");
+    } // try/catch
+    #}embed
+    cdef void* matVptr
+    cdef void* fieldInVptr
+    matVptr = PyCObject_AsVoidPtr(mat)
+    fieldInVptr = PyCObject_AsVoidPtr(fieldIn)
+    IntegratorImplicit_integrateJacobian(self.thisptr, matVptr, fieldInVptr,
+                                         ptrFromHandle(mesh))
+    return
+
+
+  property timeStep:
+    def __set__(self, dt):
+      """
+      Set timeStep.
+      """
+      # create shim for method 'timeStep'
+      #embed{ void IntegratorImplicit_timeStep_set(void* objVptr, double dt)
+      try {
+        assert(0 != objVptr);
+        ((pylith::feassemble::IntegratorImplicit*) objVptr)->timeStep(dt);
+      } catch (const std::exception& err) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        const_cast<char*>(err.what()));
+      } catch (...) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        "Caught unknown C++ exception.");
+      } // try/catch
+      #}embed
+      IntegratorImplicit_timeStep_set(self.thisptr, dt)
+
+
+# ----------------------------------------------------------------------
 cdef class ExplicitElasticity(IntegratorExplicit):
 
   def __init__(self):
@@ -623,4 +742,59 @@
       ExplicitElasticity_material_set(self.thisptr, ptrFromHandle(m))
 
 
+# ----------------------------------------------------------------------
+cdef class ImplicitElasticity(IntegratorImplicit):
+
+  def __init__(self):
+    """
+    Constructor.
+    """
+    # create shim for constructor
+    #embed{ void* ImplicitElasticity_constructor()
+    void* result = 0;
+    try {
+      result = (void*)(new pylith::feassemble::ImplicitElasticity);
+      assert(0 != result);
+    } catch (const std::exception& err) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      const_cast<char*>(err.what()));
+    } catch (...) {
+      PyErr_SetString(PyExc_RuntimeError,
+                      "Caught unknown C++ exception.");
+    } // try/catch
+    return result;
+    #}embed
+
+    IntegratorImplicit.__init__(self)
+    self.thisptr = ImplicitElasticity_constructor()
+    self.handle = self._createHandle()
+    return
+
+
+  property material:
+    def __set__(self, m):
+      """
+      Set material.
+      """
+      # create shim for method 'material'
+      #embed{ void ImplicitElasticity_material_set(void* objVptr, void* mVptr)
+      try {
+        assert(0 != objVptr);
+        pylith::materials::ElasticMaterial* material =
+          (pylith::materials::ElasticMaterial*) mVptr;
+        ((pylith::feassemble::ImplicitElasticity*) objVptr)->material(material);
+      } catch (const std::exception& err) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        const_cast<char*>(err.what()));
+      } catch (...) {
+        PyErr_SetString(PyExc_RuntimeError,
+                        "Caught unknown C++ exception.");
+      } // try/catch
+      #}embed
+      if not m.name == "pylith_materials_ElasticMaterial":
+        raise TypeError, \
+              "Argument must be extension module type 'ElasticMaterial'."
+      ImplicitElasticity_material_set(self.thisptr, ptrFromHandle(m))
+
+
 # End of file 

Modified: short/3D/PyLith/trunk/pylith/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/pylith/Makefile.am	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/pylith/Makefile.am	2007-05-18 18:37:57 UTC (rev 6924)
@@ -28,10 +28,12 @@
 	faults/SlipTimeFn.py \
 	feassemble/__init__.py \
 	feassemble/ExplicitElasticity.py \
+	feassemble/ImplicitElasticity.py \
 	feassemble/FIATCell.py \
 	feassemble/FIATSimplex.py \
 	feassemble/Integrator.py \
 	feassemble/IntegratorExplicit.py \
+	feassemble/IntegratorImplicit.py \
 	feassemble/ReferenceCell.py \
 	feassemble/quadrature/Quadrature.py \
 	feassemble/quadrature/Quadrature1D.py \

Added: short/3D/PyLith/trunk/pylith/feassemble/ImplicitElasticity.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/ImplicitElasticity.py	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/pylith/feassemble/ImplicitElasticity.py	2007-05-18 18:37:57 UTC (rev 6924)
@@ -0,0 +1,62 @@
+#!/usr/bin/env python
+#
+# ----------------------------------------------------------------------
+#
+#                           Brad T. Aagaard
+#                        U.S. Geological Survey
+#
+# <LicenseText>
+#
+# ----------------------------------------------------------------------
+#
+
+## @file pylith/feassemble/ImplicitElasticity.py
+##
+## @brief Python object for implicit time integration of dynamic
+## elasticity equation using finite-elements.
+##
+## Factory: integrator
+
+from IntegratorImplicit import IntegratorImplicit
+
+# ImplicitElasticity class
+class ImplicitElasticity(IntegratorImplicit):
+  """
+  Python object for implicit time integration of dynamic elasticity
+  equation using finite-elements.
+  """
+
+  # PUBLIC METHODS /////////////////////////////////////////////////////
+
+  def __init__(self, name="implicitelasticity"):
+    """
+    Constructor.
+    """
+    IntegratorImplicit.__init__(self, name)
+
+    import pylith.feassemble.feassemble as bindings
+    self.cppHandle = bindings.ImplicitElasticity()
+    return
+
+
+  def initMaterial(self, mesh, material):
+    """
+    Initialize material properties.
+    """
+    self._info.log("Initializing integrator material '%s'." % material.label)
+    material.initialize(mesh)
+    self.material = material
+    self.cppHandle.material = self.material.cppHandle
+    return
+  
+  
+# FACTORIES ////////////////////////////////////////////////////////////
+
+def integrator():
+  """
+  Factory associated with ImplicitElasticity.
+  """
+  return ImplicitElasticity()
+
+
+# End of file 

Added: short/3D/PyLith/trunk/pylith/feassemble/IntegratorImplicit.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/IntegratorImplicit.py	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/pylith/feassemble/IntegratorImplicit.py	2007-05-18 18:37:57 UTC (rev 6924)
@@ -0,0 +1,72 @@
+#!/usr/bin/env python
+#
+# ----------------------------------------------------------------------
+#
+#                           Brad T. Aagaard
+#                        U.S. Geological Survey
+#
+# <LicenseText>
+#
+# ----------------------------------------------------------------------
+#
+
+## @file pylith/feassemble/IntegratorImplicit.py
+##
+## @brief Python object for implicit time integration of actions with
+## finite-elements.
+##
+## Factory: integrator
+
+from Integrator import Integrator
+
+# IntegratorInertia class
+class IntegratorImplicit(Integrator):
+  """
+  Python object for implicit integration of operator actions with
+  finite-elements.
+
+  Factory: integrator.
+  """
+
+  # PUBLIC METHODS /////////////////////////////////////////////////////
+
+  def __init__(self, name="integratorimplicit"):
+    """
+    Constructor.
+    """
+    Integrator.__init__(self, name)
+    return
+
+
+  def timeStep(self, t):
+    """
+    Set time step for advancing from time t to time t+dt.
+    """
+    self.cppHandle.timeStep = t.value
+    return
+
+
+  def stableTimeStep(self):
+    """
+    Get stable time step for advancing from time t to time t+dt.
+    """
+    return self.cppHandle.getStableTimeStep()
+
+
+  def integrateResidual(self, fieldOut, fieldInT):
+    """
+    Integrate residual term for quasi-static terms for finite-elements.
+    """
+    self.cppHandle.integrateResidual(fieldOut, fieldInT, self.mesh.cppHandle)
+    return
+
+
+  def integrateJacobian(self, jacobian, fieldInT):
+    """
+    Integrate Jacobian term for quasi-static terms for finite-elements.
+    """
+    self.cppHandle.integrateJacobian(jacobian, fieldInT, self.mesh.cppHandle)
+    return
+
+
+# End of file 

Modified: short/3D/PyLith/trunk/pylith/feassemble/__init__.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/__init__.py	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/pylith/feassemble/__init__.py	2007-05-18 18:37:57 UTC (rev 6924)
@@ -18,8 +18,10 @@
            'FIATCell',
            'FIATLagrange',
            'FIATSimplex',
+           'ImplicitElasticity',
            'Integrator',
            'IntegratorExplicit',
+           'IntegratorImplicit',
            'ReferenceCell',
            'quadrature']
 

Modified: short/3D/PyLith/trunk/pylith/problems/Implicit.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/Implicit.py	2007-05-18 17:43:21 UTC (rev 6923)
+++ short/3D/PyLith/trunk/pylith/problems/Implicit.py	2007-05-18 18:37:57 UTC (rev 6924)
@@ -101,6 +101,8 @@
     mesh.allocateRealSection(self.dispTpdt)
     mesh.allocateRealSection(self.residual)
 
+    self.jacobian = mesh.createMatrix(self.residual)
+
     self._info.log("Integrating Jacobian of operator.")
     for integrator in self.integrators:
       integrator.timeStep(dt)



More information about the cig-commits mailing list