[cig-commits] r20022 - in short/3D/PyLith/branches/v1.7-trunk: . libsrc/pylith libsrc/pylith/faults pylith/faults pylith/problems
brad at geodynamics.org
brad at geodynamics.org
Tue May 1 17:30:12 PDT 2012
Author: brad
Date: 2012-05-01 17:30:11 -0700 (Tue, 01 May 2012)
New Revision: 20022
Added:
short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.cc
short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.hh
Modified:
short/3D/PyLith/branches/v1.7-trunk/TODO
short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/Makefile.am
short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/Makefile.am
short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/faultsfwd.hh
short/3D/PyLith/branches/v1.7-trunk/pylith/faults/Fault.py
short/3D/PyLith/branches/v1.7-trunk/pylith/problems/GreensFns.py
Log:
Started work on FaultCohesiveImpulses for Green's functions.
Modified: short/3D/PyLith/branches/v1.7-trunk/TODO
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/TODO 2012-05-01 21:54:45 UTC (rev 20021)
+++ short/3D/PyLith/branches/v1.7-trunk/TODO 2012-05-02 00:30:11 UTC (rev 20022)
@@ -15,14 +15,15 @@
CURRENT ISSUES/PRIORITIES
======================================================================
-* FIATSimplex
- cell.shape -> cell.dimension
* FaultCohesiveDyn
Add switch for turning on/off zero tractions for fault opening
+ Add example
* Manual
+ - FIATSimplex
+ cell.shape -> cell.dimension
- Order of tensor components for Xdmf files
- Drucker Prager fit to yield surface
- Drucker Prage allow tensile yield
Modified: short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/Makefile.am
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/Makefile.am 2012-05-01 21:54:45 UTC (rev 20021)
+++ short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/Makefile.am 2012-05-02 00:30:11 UTC (rev 20022)
@@ -67,6 +67,7 @@
faults/FaultCohesiveLagrange.cc \
faults/FaultCohesiveKin.cc \
faults/FaultCohesiveDyn.cc \
+ faults/FaultCohesiveImpulses.cc \
faults/FaultCohesiveTract.cc \
feassemble/CellGeometry.cc \
feassemble/Constraint.cc \
Added: short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.cc
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.cc (rev 0)
+++ short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.cc 2012-05-02 00:30:11 UTC (rev 20022)
@@ -0,0 +1,423 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+// Brad T. Aagaard, U.S. Geological Survey
+// Charles A. Williams, GNS Science
+// Matthew G. Knepley, University of Chicago
+//
+// This code was developed as part of the Computational Infrastructure
+// for Geodynamics (http://geodynamics.org).
+//
+// Copyright (c) 2010-2012 University of California, Davis
+//
+// See COPYING for license information.
+//
+// ----------------------------------------------------------------------
+//
+
+#include <portinfo>
+
+#include "FaultCohesiveImpulses.hh" // implementation of object methods
+
+#include "CohesiveTopology.hh" // USES CohesiveTopology
+
+#include "pylith/feassemble/Quadrature.hh" // USES Quadrature
+#include "pylith/feassemble/CellGeometry.hh" // USES CellGeometry
+#include "pylith/topology/Mesh.hh" // USES Mesh
+#include "pylith/topology/SubMesh.hh" // USES SubMesh
+#include "pylith/topology/Field.hh" // USES Field
+#include "pylith/topology/Fields.hh" // USES Fields
+#include "pylith/topology/Jacobian.hh" // USES Jacobian
+#include "pylith/topology/SolutionFields.hh" // USES SolutionFields
+
+#include "spatialdata/geocoords/CoordSys.hh" // USES CoordSys
+#include "spatialdata/spatialdb/SpatialDB.hh" // USES SpatialDB
+
+#include <cmath> // USES pow(), sqrt()
+#include <strings.h> // USES strcasecmp()
+#include <cstring> // USES strlen()
+#include <cstdlib> // USES atoi()
+#include <cassert> // USES assert()
+#include <sstream> // USES std::ostringstream
+#include <stdexcept> // USES std::runtime_error
+
+//#define PRECOMPUTE_GEOMETRY
+//#define DETAILED_EVENT_LOGGING
+
+// ----------------------------------------------------------------------
+typedef pylith::topology::Mesh::SieveMesh SieveMesh;
+typedef pylith::topology::Mesh::RealSection RealSection;
+typedef pylith::topology::SubMesh::SieveMesh SieveSubMesh;
+
+// ----------------------------------------------------------------------
+// Default constructor.
+pylith::faults::FaultCohesiveImpulses::FaultCohesiveImpulses(void) :
+ _threshold(1.0e-6),
+ _dbImpulseAmp(0)
+{ // constructor
+} // constructor
+
+// ----------------------------------------------------------------------
+// Destructor.
+pylith::faults::FaultCohesiveImpulses::~FaultCohesiveImpulses(void)
+{ // destructor
+ deallocate();
+} // destructor
+
+// ----------------------------------------------------------------------
+// Deallocate PETSc and local data structures.
+void
+pylith::faults::FaultCohesiveImpulses::deallocate(void)
+{ // deallocate
+ FaultCohesiveLagrange::deallocate();
+
+ // :TODO: Use shared pointers for amplitudes of impulses
+ _dbImpulseAmp = 0;
+} // deallocate
+
+// ----------------------------------------------------------------------
+// Sets the spatial database for amplitudes of the impulses.
+void
+pylith::faults::FaultCohesiveImpulses::dbImpulseAmp(spatialdata::spatialdb::SpatialDB* db)
+{ // dbImpulseAmp
+ _dbImpulseAmp = db;
+} // dbImpulseAmp
+
+// ----------------------------------------------------------------------
+// Set indices of fault degrees of freedom associated with
+void
+pylith::faults::FaultCohesiveImpulses::impulseDOF(const int* flags,
+ const int size)
+{ // impulseDOF
+ if (size > 0)
+ assert(flags);
+
+ _impulseDOF.resize(size);
+ for (int i=0; i < size; ++i)
+ _impulseDOF[i] = flags[i];
+} // impulseDOF
+
+// ----------------------------------------------------------------------
+// Set threshold for nonzero impulse amplitude.
+void
+pylith::faults::FaultCohesiveImpulses::threshold(const double value)
+{ // threshold
+ if (value < 0) {
+ std::ostringstream msg;
+ msg << "Threshold (" << value << ") for nonzero amplitudes of impulses "
+ "must be nonnegative";
+ throw std::runtime_error(msg.str());
+ } // if
+
+ _threshold = value;
+} // threshold
+
+// ----------------------------------------------------------------------
+// Get number of impulses.
+int
+pylith::faults::FaultCohesiveImpulses::numImpulses(void) const
+{ // numImpulses
+ return _impulsePoints.size();
+} // numImpulses
+
+// ----------------------------------------------------------------------
+// Get number of components for impulses at each point.
+int
+pylith::faults::FaultCohesiveImpulses::numComponents(void) const
+{ // numComponents
+ return _impulseDOF.size();
+} // numComponents
+
+// ----------------------------------------------------------------------
+// Initialize fault. Determine orientation and setup boundary
+void
+pylith::faults::FaultCohesiveImpulses::initialize(const topology::Mesh& mesh,
+ const PylithScalar upDir[3])
+{ // initialize
+ assert(upDir);
+ assert(_quadrature);
+ assert(_normalizer);
+
+ FaultCohesiveLagrange::initialize(mesh, upDir);
+
+ // Setup impulses
+ _setupImpulseAmp();
+ _setupImpulseOrder();
+} // initialize
+
+// ----------------------------------------------------------------------
+// Integrate contribution of cohesive cells to residual term that do
+// not require assembly across cells, vertices, or processors.
+void
+pylith::faults::FaultCohesiveImpulses::integrateResidual(
+ const topology::Field<topology::Mesh>& residual,
+ const PylithScalar t,
+ topology::SolutionFields* const fields)
+{ // integrateResidual
+ assert(fields);
+ assert(_fields);
+ assert(_logger);
+
+ const int setupEvent = _logger->eventId("FaIR setup");
+ _logger->eventBegin(setupEvent);
+
+ topology::Field<topology::SubMesh>& dispRel = _fields->get("relative disp");
+ dispRel.zero();
+ // Set impulse corresponding to current time.
+ _setRelativeDisp(dispRel, int(t+0.1));
+
+ // Transform slip from local (fault) coordinate system to relative
+ // displacement field in global coordinate system
+ _faultToGlobal(&dispRel);
+
+ _logger->eventEnd(setupEvent);
+
+ FaultCohesiveLagrange::integrateResidual(residual, t, fields);
+
+} // integrateResidual
+
+// ----------------------------------------------------------------------
+// Get vertex field associated with integrator.
+const pylith::topology::Field<pylith::topology::SubMesh>&
+pylith::faults::FaultCohesiveImpulses::vertexField(const char* name,
+ const topology::SolutionFields* fields)
+{ // vertexField
+ assert(_faultMesh);
+ assert(_quadrature);
+ assert(_normalizer);
+ assert(_fields);
+
+ const int cohesiveDim = _faultMesh->dimension();
+ const int spaceDim = _quadrature->spaceDim();
+
+ PylithScalar scale = 0.0;
+ int fiberDim = 0;
+ if (0 == strcasecmp("slip", name)) {
+ const topology::Field<topology::SubMesh>& dispRel =
+ _fields->get("relative disp");
+ _allocateBufferVectorField();
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (vector)");
+ buffer.copy(dispRel);
+ buffer.label("slip");
+ _globalToFault(&buffer);
+ return buffer;
+
+ } else if (cohesiveDim > 0 && 0 == strcasecmp("strike_dir", name)) {
+ const ALE::Obj<RealSection>& orientationSection = _fields->get(
+ "orientation").section();
+ assert(!orientationSection.isNull());
+ const ALE::Obj<RealSection>& dirSection = orientationSection->getFibration(
+ 0);
+ assert(!dirSection.isNull());
+ _allocateBufferVectorField();
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (vector)");
+ buffer.copy(dirSection);
+ buffer.label("strike_dir");
+ buffer.scale(1.0);
+ return buffer;
+
+ } else if (2 == cohesiveDim && 0 == strcasecmp("dip_dir", name)) {
+ const ALE::Obj<RealSection>& orientationSection = _fields->get(
+ "orientation").section();
+ assert(!orientationSection.isNull());
+ const ALE::Obj<RealSection>& dirSection = orientationSection->getFibration(
+ 1);
+ _allocateBufferVectorField();
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (vector)");
+ buffer.copy(dirSection);
+ buffer.label("dip_dir");
+ buffer.scale(1.0);
+ return buffer;
+
+ } else if (0 == strcasecmp("normal_dir", name)) {
+ const ALE::Obj<RealSection>& orientationSection = _fields->get(
+ "orientation").section();
+ assert(!orientationSection.isNull());
+ const int space = (0 == cohesiveDim) ? 0 : (1 == cohesiveDim) ? 1 : 2;
+ const ALE::Obj<RealSection>& dirSection = orientationSection->getFibration(
+ space);
+ assert(!dirSection.isNull());
+ _allocateBufferVectorField();
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (vector)");
+ buffer.copy(dirSection);
+ buffer.label("normal_dir");
+ buffer.scale(1.0);
+ return buffer;
+
+ } else if (0 == strcasecmp("impulse_amplitude", name)) {
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (scalar)");
+ return buffer;
+
+ } else if (0 == strcasecmp("impulses", name)) {
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (scalar)");
+ return buffer;
+
+ } else if (0 == strcasecmp("traction_change", name)) {
+ assert(fields);
+ const topology::Field<topology::Mesh>& dispT = fields->get("disp(t)");
+ _allocateBufferVectorField();
+ topology::Field<topology::SubMesh>& buffer =
+ _fields->get("buffer (vector)");
+ _calcTractionsChange(&buffer, dispT);
+ return buffer;
+
+ } else {
+ std::ostringstream msg;
+ msg << "Request for unknown vertex field '" << name << "' for fault '"
+ << label() << "'.";
+ throw std::runtime_error(msg.str());
+ } // else
+
+
+ // Should never get here.
+ throw std::logic_error("Unknown field in FaultCohesiveImpulses::vertexField().");
+
+ // Satisfy return values
+ assert(_fields);
+ const topology::Field<topology::SubMesh>& buffer = _fields->get("buffer (vector)");
+ return buffer;
+} // vertexField
+
+
+// ----------------------------------------------------------------------
+// Setup amplitudes of impulses.
+void
+pylith::faults::FaultCohesiveImpulses::_setupImpulseAmp(void)
+{ // _setupImpulseAmp
+ // If no impulse amplitude specified, leave method
+ if (0 == _dbImpulseAmp)
+ return;
+
+ assert(_normalizer);
+ const PylithScalar lengthScale = _normalizer->lengthScale();
+
+ const int spaceDim = _quadrature->spaceDim();
+
+ // Create section to hold amplitudes of impulses.
+ _fields->add("impulse amplitude", "impulse_amplitude");
+ topology::Field<topology::SubMesh>& amplitude =
+ _fields->get("impulse amplitude");
+ topology::Field<topology::SubMesh>& dispRel = _fields->get("relative disp");
+ amplitude.cloneSection(dispRel);
+ amplitude.scale(lengthScale);
+
+ PylithScalar amplitudeVertex;
+ const ALE::Obj<RealSection>& amplitudeSection = amplitude.section();
+ assert(!amplitudeSection.isNull());
+
+ const spatialdata::geocoords::CoordSys* cs = _faultMesh->coordsys();
+ assert(cs);
+
+ const ALE::Obj<SieveSubMesh>& faultSieveMesh = _faultMesh->sieveMesh();
+ assert(!faultSieveMesh.isNull());
+
+ scalar_array coordsVertex(spaceDim);
+ const ALE::Obj<RealSection>& coordsSection =
+ faultSieveMesh->getRealSection("coordinates");
+ assert(!coordsSection.isNull());
+
+ assert(_dbImpulseAmp);
+ _dbImpulseAmp->open();
+ const char* valueNames[1] = { "slip" };
+ _dbImpulseAmp->queryVals(valueNames, 1);
+
+ const int numVertices = _cohesiveVertices.size();
+ for (int iVertex=0; iVertex < numVertices; ++iVertex) {
+ const int v_fault = _cohesiveVertices[iVertex].fault;
+
+ coordsSection->restrictPoint(v_fault, &coordsVertex[0], coordsVertex.size());
+ _normalizer->dimensionalize(&coordsVertex[0], coordsVertex.size(),
+ lengthScale);
+
+ amplitudeVertex = 0.0;
+ int err = _dbImpulseAmp->query(&litudeVertex, 1,
+ &coordsVertex[0], coordsVertex.size(), cs);
+ if (err) {
+ std::ostringstream msg;
+ msg << "Could not find amplitude for Green's function impulses at \n" << "(";
+ for (int i = 0; i < spaceDim; ++i)
+ msg << " " << coordsVertex[i];
+ msg << ") in fault " << label() << "\n"
+ << "using spatial database '" << _dbImpulseAmp->label() << "'.";
+ throw std::runtime_error(msg.str());
+ } // if
+ _normalizer->nondimensionalize(&litudeVertex, 1, lengthScale);
+
+ assert(1 == amplitudeSection->getFiberDimension(v_fault));
+ amplitudeSection->updatePoint(v_fault, &litudeVertex);
+ } // for
+
+ // Close properties database
+ _dbImpulseAmp->close();
+
+ amplitude.view("IMPULSE AMPLITUDE"); // DEBUGGING
+} // _setupImpulseAmp
+
+
+// ----------------------------------------------------------------------
+// Setup amplitudes of impulses.
+void
+pylith::faults::FaultCohesiveImpulses::_setupImpulseOrder(void)
+{ // _setupImpulseOrder
+ // If no impulse amplitude specified, leave method
+ if (0 == _dbImpulseAmp)
+ return;
+
+} // _setupImpulseOrder
+
+
+// ----------------------------------------------------------------------
+// Set relative displacemet associated with impulse.
+void
+pylith::faults::FaultCohesiveImpulses::_setRelativeDisp(const topology::Field<topology::SubMesh>& dispRel,
+ const int impulse)
+{ // _setRelativeDisp
+ assert(_fields);
+
+ // If no impulse amplitude specified, leave method
+ if (0 == _dbImpulseAmp)
+ return;
+
+ const spatialdata::geocoords::CoordSys* cs = _faultMesh->coordsys();
+ assert(cs);
+ const int spaceDim = cs->spaceDim();
+
+ const ALE::Obj<RealSection>& amplitudeSection = _fields->get("impulse amplitude").section();
+ assert(!amplitudeSection.isNull());
+
+ const ALE::Obj<RealSection>& dispRelSection = dispRel.section();
+ assert(!dispRelSection.isNull());
+
+ scalar_array dispRelVertex(spaceDim);
+ dispRelVertex = 0.0;
+
+ const srcs_type::const_iterator& impulseInfo = _impulsePoints.find(impulse);
+ if (impulseInfo != _impulsePoints.end()) {
+ const int iVertex = impulseInfo->second.indexCohesive;
+ const int v_fault = _cohesiveVertices[iVertex].fault;
+ const int v_lagrange = _cohesiveVertices[iVertex].lagrange;
+
+ // Get amplitude of slip impulse
+ assert(1 == amplitudeSection->getFiberDimension(v_fault));
+ const PylithScalar* amplitudeVertex = amplitudeSection->restrictPoint(v_fault);
+ assert(amplitudeVertex);
+
+ const int indexDOF = impulseInfo->second.indexDOF;
+ assert(indexDOF >= 0 && indexDOF < spaceDim);
+ dispRelVertex[indexDOF] = amplitudeVertex[0];
+
+ assert(dispRelVertex.size() == dispRelSection->getFiberDimension(v_lagrange));
+ dispRelSection->updatePoint(v_lagrange, &dispRelVertex[0]);
+ } // if
+
+} // _setRelativeDisp
+
+
+// End of file
Added: short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.hh
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.hh (rev 0)
+++ short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/FaultCohesiveImpulses.hh 2012-05-02 00:30:11 UTC (rev 20022)
@@ -0,0 +1,177 @@
+// -*- C++ -*-
+//
+// ----------------------------------------------------------------------
+//
+// Brad T. Aagaard, U.S. Geological Survey
+// Charles A. Williams, GNS Science
+// Matthew G. Knepley, University of Chicago
+//
+// This code was developed as part of the Computational Infrastructure
+// for Geodynamics (http://geodynamics.org).
+//
+// Copyright (c) 2010-2012 University of California, Davis
+//
+// See COPYING for license information.
+//
+// ----------------------------------------------------------------------
+//
+
+/** @file libsrc/faults/FaultCohesiveImpulses.hh
+ *
+ * @brief C++ implementation for Green's functions impulses
+ * implemented with cohesive elements.
+ */
+
+#if !defined(pylith_faults_faultcohesiveimpulses_hh)
+#define pylith_faults_faultcohesiveimpulses_hh
+
+// Include directives ---------------------------------------------------
+#include "FaultCohesiveLagrange.hh" // ISA FaultCohesive
+
+#include <map> // HASA std::map
+
+// FaultCohesiveImpulses -----------------------------------------------------
+/**
+ * @brief C++ implementation for Green's functions impulses
+ * implemented with cohesive elements.
+ */
+class pylith::faults::FaultCohesiveImpulses : public FaultCohesiveLagrange
+{ // class FaultCohesiveImpulses
+ friend class TestFaultCohesiveImpulses; // unit testing
+
+ // PUBLIC METHODS /////////////////////////////////////////////////////
+public :
+
+ /// Default constructor.
+ FaultCohesiveImpulses(void);
+
+ /// Destructor.
+ ~FaultCohesiveImpulses(void);
+
+ /// Deallocate PETSc and local data structures.
+ void deallocate(void);
+
+ /** Sets the spatial database for amplitudes of the impulses.
+ *
+ * @param db spatial database for amplitudes of impulses.
+ */
+ void dbImpulseAmp(spatialdata::spatialdb::SpatialDB* db);
+
+ /** Set indices of fault degrees of freedom associated with
+ * impulses.
+ *
+ * @param flags Array of indices for degrees of freedom.
+ * @param size Size of array
+ */
+ void impulseDOF(const int* flags,
+ const int size);
+
+ /** Set threshold for nonzero impulse amplitude.
+ *
+ * @param value Threshold for detecting nonzero amplitude.
+ */
+ void threshold(const double value);
+
+ /** Get number of impulses.
+ *
+ * Multiply by number of components to get total number of impulses.
+ *
+ * @returns Number of points with impulses.
+ */
+ int numImpulses(void) const;
+
+ /** Get number of components for impulses at each point.
+ *
+ * Multiply by number of components to get total number of impulses.
+ *
+ * @returns Number of points with impulses.
+ */
+ int numComponents(void) const;
+
+ /** 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; applies to fault surfaces in 2-D and 3-D).
+ */
+ void initialize(const topology::Mesh& mesh,
+ const PylithScalar upDir[3]);
+
+ /** 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 integrateResidual(const topology::Field<topology::Mesh>& residual,
+ const PylithScalar t,
+ topology::SolutionFields* const fields);
+
+ /** Get vertex field associated with integrator.
+ *
+ * @param name Name of cell field.
+ * @param fields Solution fields.
+ * @returns Vertex field.
+ */
+ const topology::Field<topology::SubMesh>&
+ vertexField(const char* name,
+ const topology::SolutionFields* fields =0);
+
+ // PRIVATE METHODS ////////////////////////////////////////////////////
+private :
+
+ /// Setup amplitudes of impulses.
+ void _setupImpulseAmp(void);
+
+ /// Setup impulse order.
+ void _setupImpulseOrder(void);
+
+ /** Set relative displacemet associated with impulse.
+ *
+ * @param dispRel Relative displacement field.
+ * @parm impulse Index of impulse.
+ */
+ void _setRelativeDisp(const topology::Field<topology::SubMesh>& dispRel,
+ const int impulse);
+
+ // PRIVATE TYPEDEFS ///////////////////////////////////////////////////
+private :
+
+ struct ImpulseInfoStruct {
+ int indexCohesive;
+ int indexDOF;
+ }; // ImpulseInfoStruct
+ typedef std::map<int, ImpulseInfoStruct> srcs_type;
+
+ // PRIVATE MEMBERS ////////////////////////////////////////////////////
+private :
+
+ /// Threshold for nonzero impulse amplitude.
+ PylithScalar _threshold;
+
+ /// Database for amplitudes of impulses.
+ spatialdata::spatialdb::SpatialDB* _dbImpulseAmp;
+
+ /// Map from impulse index to corresponding point
+ srcs_type _impulsePoints;
+
+ int_array _impulseDOF; ///< Degrees of freedom associated with impulses.
+
+ // NOT IMPLEMENTED ////////////////////////////////////////////////////
+private :
+
+ /// Not implemented
+ FaultCohesiveImpulses(const FaultCohesiveImpulses&);
+
+ /// Not implemented
+ const FaultCohesiveImpulses& operator=(const FaultCohesiveImpulses&);
+
+}; // class FaultCohesiveImpulses
+
+#endif // pylith_faults_faultcohesiveimpulses_hh
+
+
+// End of file
Modified: short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/Makefile.am
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/Makefile.am 2012-05-01 21:54:45 UTC (rev 20021)
+++ short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/Makefile.am 2012-05-02 00:30:11 UTC (rev 20022)
@@ -34,6 +34,7 @@
FaultCohesiveTract.hh \
FaultCohesiveDyn.hh \
FaultCohesiveKin.hh \
+ FaultCohesiveImpulses.hh \
LiuCosSlipFn.hh \
LiuCosSlipFn.icc \
SlipTimeFn.hh \
Modified: short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/faultsfwd.hh
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/faultsfwd.hh 2012-05-01 21:54:45 UTC (rev 20021)
+++ short/3D/PyLith/branches/v1.7-trunk/libsrc/pylith/faults/faultsfwd.hh 2012-05-02 00:30:11 UTC (rev 20022)
@@ -37,6 +37,7 @@
class FaultCohesiveLagrange;
class FaultCohesiveKin;
class FaultCohesiveDyn;
+ class FaultCohesiveImpulses;
class FaultCohesiveTract;
class EqKinSrc;
Modified: short/3D/PyLith/branches/v1.7-trunk/pylith/faults/Fault.py
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/pylith/faults/Fault.py 2012-05-01 21:54:45 UTC (rev 20021)
+++ short/3D/PyLith/branches/v1.7-trunk/pylith/faults/Fault.py 2012-05-02 00:30:11 UTC (rev 20022)
@@ -68,7 +68,7 @@
\b Properties
@li \b id Fault identifier
- @li \b name Name of fault
+ @li \b label Label identifier for fault.
@li \b up_dir Up-dip or up direction
(perpendicular to along-strike and not collinear with fault normal;
applies to fault surfaces in 2-D and 3-D).
Modified: short/3D/PyLith/branches/v1.7-trunk/pylith/problems/GreensFns.py
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/pylith/problems/GreensFns.py 2012-05-01 21:54:45 UTC (rev 20021)
+++ short/3D/PyLith/branches/v1.7-trunk/pylith/problems/GreensFns.py 2012-05-02 00:30:11 UTC (rev 20022)
@@ -44,7 +44,7 @@
## Python object for managing GreensFns facilities and properties.
##
## \b Properties
- ## None
+ ## @li \b faultId Id of fault on which to impose impulses.
##
## \b Facilities
## @li \b formulation Formulation for solving PDE.
@@ -52,6 +52,9 @@
import pyre.inventory
+ faultId = pyre.inventory.int("fault_id", default=100)
+ faultId.meta['tip'] = "Id of fault on which to impose impulses."
+
from Implicit import Implicit
formulation = pyre.inventory.facility("formulation",
family="pde_formulation",
@@ -90,6 +93,17 @@
self.mesh = mesh
self.formulation.preinitialize(mesh, self.materials, self.bc,
self.interfaces, self.gravityField)
+
+ # Find fault for impulses
+ found = False
+ for fault in self.interfaces.components():
+ if self.faultId == fault.id():
+ self.source = fault
+ found = True
+ break
+ if not found:
+ raise ValueError("Could not find fault interface with id '%d' for "
+ "Green's function impulses." % self.faultId)
return
@@ -99,6 +113,11 @@
"""
Problem.verifyConfiguration(self)
self.formulation.verifyConfiguration()
+
+ if not "numImpulses" in dir(self.source) or not "numComponents" in dir(self.source):
+ raise ValueError("Incompatible source for green's function impulses "
+ "with id '%d' and label '%s'." % \
+ (self.source.id(), self.source.label())
return
@@ -128,7 +147,7 @@
self._info.log("Computing Green's functions.")
self.checkpointTimer.toplevel = app # Set handle for saving state
- nimpulses = 0
+ nimpulses = self.source.numImpulses()*self.source.numComponents()
ipulse = 0;
while ipulse < nimpulses:
self._eventLogger.stagePush("Prestep")
@@ -192,6 +211,8 @@
Set members based using inventory.
"""
Problem._configure(self)
+
+ self.faultId = self.inventory.faultId
self.formulation = self.inventory.formulation
self.checkpointTimer = self.inventory.checkpointTimer
return
More information about the CIG-COMMITS
mailing list