[cig-commits] r15792 - in short/3D/PyLith/branches/pylith-friction: libsrc/faults modulesrc/faults playpen/friction pylith pylith/faults
brad at geodynamics.org
brad at geodynamics.org
Fri Oct 9 15:17:01 PDT 2009
Author: brad
Date: 2009-10-09 15:17:00 -0700 (Fri, 09 Oct 2009)
New Revision: 15792
Added:
short/3D/PyLith/branches/pylith-friction/modulesrc/faults/FaultCohesiveDynL.i
short/3D/PyLith/branches/pylith-friction/pylith/faults/FaultCohesiveDynL.py
Modified:
short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.cc
short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.hh
short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.icc
short/3D/PyLith/branches/pylith-friction/modulesrc/faults/Makefile.am
short/3D/PyLith/branches/pylith-friction/modulesrc/faults/faults.i
short/3D/PyLith/branches/pylith-friction/playpen/friction/twoquad4-axial.cfg
short/3D/PyLith/branches/pylith-friction/pylith/Makefile.am
Log:
Added module and Python files for FaultCohesiveDynL.
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.cc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.cc 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.cc 2009-10-09 22:17:00 UTC (rev 15792)
@@ -628,9 +628,6 @@
assert(0 != _quadrature);
assert(0 != _fields);
- // Cohesive cells with normal vertices i and j, and constraint
- // vertex k.
-
// Get cell information and setup storage for cell data
const int spaceDim = _quadrature->spaceDim();
const int numBasis = _quadrature->numBasis();
@@ -643,6 +640,7 @@
double_array dispTIncrCell(numCorners*spaceDim);
double_array dispTpdtCell(numCorners*spaceDim);
double_array tractionVertex(spaceDim);
+ double_array tractionStickVertex(spaceDim);
double_array slipVertex(spaceDim);
// Total fault area associated with each vertex (assembled over all cells).
@@ -689,6 +687,7 @@
for (SieveMesh::label_sequence::iterator c_iter=cellsCohesiveBegin;
c_iter != cellsCohesiveEnd;
++c_iter) {
+ // Get fault cell label associated with cohesive cell
const SieveMesh::point_type c_fault = _cohesiveToFault[*c_iter];
// Get area at fault cell's vertices.
@@ -717,11 +716,18 @@
assert(areaCell[iConstraint] > 0);
const double areaVertex = areaCell[iConstraint];
- // :KLUDGE: (TEMPORARY) Compute traction by dividing force by area
+ // :KLUDGE: (TEMPORARY) Solution at Lagrange constraint vertices
+ // is the Lagrange multiplier value, which is currently the
+ // force. Compute traction by dividing force by area
for (int i=0; i < spaceDim; ++i)
tractionVertex[i] = dispTpdtCell[indexK+i] / areaVertex;
- // Compute the current slip from current displacements.
+ // Tractions required for sticking (if not exceeded slip
+ // increment is zero).
+ tractionStickVertex = tractionVertex;
+
+ // Compute the current slip at the constraint vertex from
+ // current displacements at the normal cohesive vertices.
for (int i=0; i < spaceDim; ++i)
slipVertex[i] = dispTpdtCell[indexJ+i] - dispTpdtCell[indexI+i];
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.hh
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.hh 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.hh 2009-10-09 22:17:00 UTC (rev 15792)
@@ -16,8 +16,8 @@
* (dynamic) slip implemented with cohesive elements.
*/
-#if !defined(pylith_faults_faultcohesivedyn_hh)
-#define pylith_faults_faultcohesivedyn_hh
+#if !defined(pylith_faults_faultcohesivedynl_hh)
+#define pylith_faults_faultcohesivedynl_hh
// Include directives ---------------------------------------------------
#include "FaultCohesive.hh" // ISA FaultCohesive
@@ -287,7 +287,7 @@
#include "FaultCohesiveDynL.icc" // inline methods
-#endif // pylith_faults_faultcohesivedyn_hh
+#endif // pylith_faults_faultcohesivedynl_hh
// End of file
Modified: short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.icc
===================================================================
--- short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.icc 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/libsrc/faults/FaultCohesiveDynL.icc 2009-10-09 22:17:00 UTC (rev 15792)
@@ -10,7 +10,7 @@
// ----------------------------------------------------------------------
//
-#if !defined(pylith_faults_faultcohesivedyn_hh)
+#if !defined(pylith_faults_faultcohesivedynl_hh)
#error "FaultCohesiveDynL.icc can only be included from FaultCohesiveDynL.hh"
#endif
Added: short/3D/PyLith/branches/pylith-friction/modulesrc/faults/FaultCohesiveDynL.i
===================================================================
--- short/3D/PyLith/branches/pylith-friction/modulesrc/faults/FaultCohesiveDynL.i (rev 0)
+++ short/3D/PyLith/branches/pylith-friction/modulesrc/faults/FaultCohesiveDynL.i 2009-10-09 22:17:00 UTC (rev 15792)
@@ -0,0 +1,165 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ----------------------------------------------------------------------
+//
+
+/** @file modulesrc/faults/FaultCohesiveDynL.i
+ *
+ * @brief Python interface to C++ FaultCohesiveDynL object.
+ */
+
+namespace pylith {
+ namespace faults {
+
+ class FaultCohesiveDynL : public FaultCohesive
+ { // class FaultCohesiveDynL
+
+ // PUBLIC METHODS /////////////////////////////////////////////////
+ public :
+
+ /// Default constructor.
+ FaultCohesiveDynL(void);
+
+ /// Destructor.
+ virtual
+ ~FaultCohesiveDynL(void);
+
+ /// Deallocate PETSc and local data structures.
+ virtual
+ void deallocate(void);
+
+ /** Sets the spatial database for the inital tractions.
+ *
+ * @param db spatial database for initial tractions
+ */
+ void dbInitial(spatialdata::spatialdb::SpatialDB* db);
+
+ /** Initialize fault. Determine orientation and setup boundary
+ * condition parameters.
+ *
+ * @param mesh Finite-element mesh.
+ * @param upDir Direction perpendicular to along-strike direction that is
+ * not collinear with fault normal (usually "up" direction but could
+ * be up-dip direction; only applies to fault surfaces in a 3-D domain).
+ * @param normalDir General preferred direction for fault normal
+ * (used to pick which of two possible normal directions for
+ * interface; only applies to fault surfaces in a 3-D domain).
+ */
+ void initialize(const pylith::topology::Mesh& mesh,
+ const double upDir[3],
+ const double normalDir[3]);
+
+ /** Split solution field for separate preconditioning.
+ *
+ * @param field Solution field.
+ */
+ void splitField(pylith::topology::Field<pylith::topology::Mesh>* field);
+
+ /** Integrate contributions to residual term (r) for operator that
+ * require assembly across processors.
+ *
+ * @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 residual term (r) for operator that
+ * do not require assembly across cells, vertices, or processors.
+ *
+ * @param residual Field containing values for residual
+ * @param t Current time
+ * @param fields Solution fields
+ */
+ void integrateResidualAssembled(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 that do not require assembly across cells, vertices, or
+ * processors.
+ *
+ * @param mat Sparse matrix
+ * @param t Current time
+ * @param fields Solution fields
+ * @param mesh Finite-element mesh
+ */
+ void integrateJacobianAssembled(pylith::topology::Jacobian* jacobian,
+ const double t,
+ pylith::topology::SolutionFields* const fields);
+
+ /** 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);
+
+ /** Constrain solution space based on friction.
+ *
+ * @param fields Solution fields.
+ * @param t Current time.
+ * @param jacobian Sparse matrix for system Jacobian.
+ */
+ void constrainSolnSpace(pylith::topology::SolutionFields* const fields,
+ const double t,
+ const pylith::topology::Jacobian& jacobian);
+
+ /** Verify configuration is acceptable.
+ *
+ * @param mesh Finite-element mesh
+ */
+ void verifyConfiguration(const pylith::topology::Mesh& mesh) const;
+
+ /** Get vertex field associated with integrator.
+ *
+ * @param name Name of cell field.
+ * @param fields Solution fields.
+ * @returns Vertex field.
+ */
+ const pylith::topology::Field<pylith::topology::SubMesh>&
+ vertexField(const char* name,
+ const pylith::topology::SolutionFields* fields =0);
+
+ /** Get cell field associated with integrator.
+ *
+ * @param name Name of cell field.
+ * @param fields Solution fields.
+ * @returns Cell field.
+ */
+ const pylith::topology::Field<pylith::topology::SubMesh>&
+ cellField(const char* name,
+ const pylith::topology::SolutionFields* fields =0);
+
+ /** Cohesive cells use Lagrange multiplier constraints?
+ *
+ * @returns True if implementation using Lagrange multiplier
+ * constraints, false otherwise.
+ */
+ bool useLagrangeConstraints(void) const;
+
+ /** Get fields associated with fault.
+ *
+ * @returns Fields associated with fault.
+ */
+ const pylith::topology::Fields<pylith::topology::Field<pylith::topology::SubMesh> >*
+ fields(void) const;
+
+ }; // class FaultCohesiveDynL
+
+ } // faults
+} // pylith
+
+
+// End of file
Modified: short/3D/PyLith/branches/pylith-friction/modulesrc/faults/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-friction/modulesrc/faults/Makefile.am 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/modulesrc/faults/Makefile.am 2009-10-09 22:17:00 UTC (rev 15792)
@@ -28,6 +28,7 @@
Fault.i \
FaultCohesive.i \
FaultCohesiveDyn.i \
+ FaultCohesiveDynL.i \
FaultCohesiveKin.i \
../topology/SubMesh.i \
../feassemble/Quadrature.i \
Modified: short/3D/PyLith/branches/pylith-friction/modulesrc/faults/faults.i
===================================================================
--- short/3D/PyLith/branches/pylith-friction/modulesrc/faults/faults.i 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/modulesrc/faults/faults.i 2009-10-09 22:17:00 UTC (rev 15792)
@@ -24,6 +24,7 @@
#include "pylith/faults/Fault.hh"
#include "pylith/faults/FaultCohesive.hh"
#include "pylith/faults/FaultCohesiveDyn.hh"
+#include "pylith/faults/FaultCohesiveDynL.hh"
#include "pylith/faults/FaultCohesiveKin.hh"
#include "pylith/topology/SubMesh.hh"
@@ -71,6 +72,7 @@
%include "Fault.i"
%include "FaultCohesive.i"
%include "FaultCohesiveDyn.i"
+%include "FaultCohesiveDynL.i"
%include "FaultCohesiveKin.i"
// End of file
Modified: short/3D/PyLith/branches/pylith-friction/playpen/friction/twoquad4-axial.cfg
===================================================================
--- short/3D/PyLith/branches/pylith-friction/playpen/friction/twoquad4-axial.cfg 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/playpen/friction/twoquad4-axial.cfg 2009-10-09 22:17:00 UTC (rev 15792)
@@ -145,7 +145,7 @@
# Provide information on the fault (interface).
[pylithapp.timedependent.interfaces]
-fault = pylith.faults.FaultCohesiveDyn
+fault = pylith.faults.FaultCohesiveDynL
# Define fault properties.
[pylithapp.timedependent.interfaces.fault]
@@ -209,8 +209,8 @@
# Give basename for VTK fault output.
[pylithapp.timedependent.interfaces.fault.output]
writer.filename = twoquad4-axial-fault.vtk
-cell_info_fields = []
-cell_data_fields = []
+vertex_info_fields = []
+vertex_data_fields = []
# Give basename for VTK output of state variables.
[pylithapp.timedependent.materials.material.output]
Modified: short/3D/PyLith/branches/pylith-friction/pylith/Makefile.am
===================================================================
--- short/3D/PyLith/branches/pylith-friction/pylith/Makefile.am 2009-10-09 21:34:33 UTC (rev 15791)
+++ short/3D/PyLith/branches/pylith-friction/pylith/Makefile.am 2009-10-09 22:17:00 UTC (rev 15792)
@@ -33,6 +33,7 @@
faults/FaultCohesive.py \
faults/FaultCohesiveKin.py \
faults/FaultCohesiveDyn.py \
+ faults/FaultCohesiveDynL.py \
faults/LiuCosSlipFn.py \
faults/SingleRupture.py \
faults/SlipTimeFn.py \
Added: short/3D/PyLith/branches/pylith-friction/pylith/faults/FaultCohesiveDynL.py
===================================================================
--- short/3D/PyLith/branches/pylith-friction/pylith/faults/FaultCohesiveDynL.py (rev 0)
+++ short/3D/PyLith/branches/pylith-friction/pylith/faults/FaultCohesiveDynL.py 2009-10-09 22:17:00 UTC (rev 15792)
@@ -0,0 +1,219 @@
+#!/usr/bin/env python
+#
+# ----------------------------------------------------------------------
+#
+# Brad T. Aagaard
+# U.S. Geological Survey
+#
+# <LicenseText>
+#
+# ----------------------------------------------------------------------
+#
+
+## @file pylith/faults/FaultCohesiveDynL.py
+##
+
+## @brief Python object for a fault surface with dynamic
+## (friction) fault implemented with cohesive elements.
+##
+## Factory: fault
+
+from FaultCohesive import FaultCohesive
+from pylith.feassemble.Integrator import Integrator
+from faults import FaultCohesiveDynL as ModuleFaultCohesiveDynL
+
+from pylith.utils.NullComponent import NullComponent
+
+# FaultCohesiveDynL class
+class FaultCohesiveDynL(FaultCohesive, Integrator, ModuleFaultCohesiveDynL):
+ """
+ Python object for a fault surface with kinematic (prescribed) slip
+ implemented with cohesive elements.
+
+ Inventory
+
+ @class Inventory
+ Python object for managing FaultCohesiveDynL facilities and properties.
+
+ \b Properties
+ @li None
+
+ \b Facilities
+ @li \b db_initial_tractions Spatial database for initial tractions.
+ @li \b output Output manager associated with fault data.
+
+ Factory: fault
+ """
+
+ # INVENTORY //////////////////////////////////////////////////////////
+
+ import pyre.inventory
+
+ db = pyre.inventory.facility("db_initial_tractions", family="spatial_database",
+ factory=NullComponent)
+ db.meta['tip'] = "Spatial database for initial tractions."
+
+ from pylith.meshio.OutputFaultKin import OutputFaultKin
+ output = pyre.inventory.facility("output", family="output_manager",
+ factory=OutputFaultKin)
+ output.meta['tip'] = "Output manager associated with fault data."
+
+
+ # PUBLIC METHODS /////////////////////////////////////////////////////
+
+ def __init__(self, name="faultcohesivedynl"):
+ """
+ Initialize configuration.
+ """
+ FaultCohesive.__init__(self, name)
+ Integrator.__init__(self)
+ self._loggingPrefix = "CoDy "
+
+ self.availableFields = \
+ {'vertex': \
+ {'info': ["normal_dir"],
+ 'data': ["slip",
+ "traction"]},
+ 'cell': \
+ {'info': [],
+ 'data': []}}
+ return
+
+
+ def preinitialize(self, mesh):
+ """
+ Do pre-initialization setup.
+ """
+ self._info.log("Pre-initializing fault '%s'." % self.label())
+ FaultCohesive.preinitialize(self, mesh)
+ Integrator.preinitialize(self, mesh)
+
+ ModuleFaultCohesiveDynL.quadrature(self, self.faultQuadrature)
+
+ if mesh.dimension() == 2:
+ self.availableFields['vertex']['info'] += ["strike_dir"]
+ elif mesh.dimension() == 3:
+ self.availableFields['vertex']['info'] += ["strike_dir",
+ "dip_dir"]
+
+ if not isinstance(self.inventory.db, NullComponent):
+ self.availableFields['vertex']['info'] += ["initial_traction"]
+ return
+
+
+ def verifyConfiguration(self):
+ """
+ Verify compatibility of configuration.
+ """
+ logEvent = "%sverify" % self._loggingPrefix
+ self._eventLogger.eventBegin(logEvent)
+
+ FaultCohesive.verifyConfiguration(self)
+ Integrator.verifyConfiguration(self)
+ ModuleFaultCohesiveDynL.verifyConfiguration(self, self.mesh)
+
+ self._eventLogger.eventEnd(logEvent)
+ return
+
+
+ def initialize(self, totalTime, numTimeSteps, normalizer):
+ """
+ Initialize cohesive elements.
+ """
+ logEvent = "%sinit" % self._loggingPrefix
+ self._eventLogger.eventBegin(logEvent)
+ self._info.log("Initializing fault '%s'." % self.label())
+
+ Integrator.initialize(self, totalTime, numTimeSteps, normalizer)
+
+ FaultCohesive.initialize(self, totalTime, numTimeSteps, normalizer)
+
+ self._eventLogger.eventEnd(logEvent)
+ return
+
+
+ def poststep(self, t, dt, totalTime, fields):
+ """
+ Hook for doing stuff after advancing time step.
+ """
+ logEvent = "%spoststep" % self._loggingPrefix
+ self._eventLogger.eventBegin(logEvent)
+
+ Integrator.poststep(self, t, dt, totalTime, fields)
+ FaultCohesive.poststep(self, t, dt, totalTime, fields)
+
+ self._eventLogger.eventEnd(logEvent)
+ return
+
+
+ def getVertexField(self, name, fields=None):
+ """
+ Get vertex field.
+ """
+ if None == fields:
+ field = ModuleFaultCohesiveDynL.vertexField(self, name)
+ else:
+ field = ModuleFaultCohesiveDynL.vertexField(self, name, fields)
+ return field
+
+
+ def getCellField(self, name, fields=None):
+ """
+ Get cell field.
+ """
+ if None == fields:
+ field = ModuleFaultCohesiveDynL.cellField(self, name)
+ else:
+ field = ModuleFaultCohesiveDynL.cellField(self, name, fields)
+ return field
+
+
+ def finalize(self):
+ """
+ Cleanup.
+ """
+ FaultCohesive.finalize(self)
+ Integrator.finalize(self)
+ return
+
+
+ # PRIVATE METHODS ////////////////////////////////////////////////////
+
+ def _configure(self):
+ """
+ Setup members using inventory.
+ """
+ FaultCohesive._configure(self)
+ if not isinstance(self.inventory.db, NullComponent):
+ ModuleFaultCohesiveDynL.dbInitial(self, self.inventory.db)
+ self.output = self.inventory.output
+ return
+
+
+ def _createModuleObj(self):
+ """
+ Create handle to C++ FaultCohesiveDynL.
+ """
+ ModuleFaultCohesiveDynL.__init__(self)
+ return
+
+
+ def _modelMemoryUse(self):
+ """
+ Model memory allocation.
+ """
+ self.perfLogger.logFault("Fault", self)
+ self.perfLogger.logFields("Fault", self.fields())
+ return
+
+
+# FACTORIES ////////////////////////////////////////////////////////////
+
+def fault():
+ """
+ Factory associated with FaultCohesiveDynL.
+ """
+ return FaultCohesiveDynL()
+
+
+# End of file
More information about the CIG-COMMITS
mailing list