[cig-commits] r7194 - in short/3D/PyLith/trunk: .
examples/templates/containers libsrc/faults libsrc/meshio
pylith/meshio pylith/problems
brad at geodynamics.org
brad at geodynamics.org
Wed Jun 13 08:11:20 PDT 2007
Author: brad
Date: 2007-06-13 08:11:19 -0700 (Wed, 13 Jun 2007)
New Revision: 7194
Added:
short/3D/PyLith/trunk/examples/templates/containers/faults.odb
short/3D/PyLith/trunk/examples/templates/containers/materials.odb
Modified:
short/3D/PyLith/trunk/TODO
short/3D/PyLith/trunk/examples/templates/containers/bc.odb
short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh
short/3D/PyLith/trunk/libsrc/meshio/SolutionIOVTK.cc
short/3D/PyLith/trunk/pylith/meshio/SolutionIO.py
short/3D/PyLith/trunk/pylith/problems/TimeDependent.py
Log:
Updated VTK output so that each time step is written to a different file. Added more templates for user-defined containers.
Modified: short/3D/PyLith/trunk/TODO
===================================================================
--- short/3D/PyLith/trunk/TODO 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/TODO 2007-06-13 15:11:19 UTC (rev 7194)
@@ -10,10 +10,6 @@
dislocation.cfg
-solutioniovtk
-
- need to write each time step to a different file.
-
1. Simple tests with analytical solutions
a. 2-D
@@ -40,10 +36,8 @@
i. Add checking of faultMesh [not currently used]
examples/templates
-
containers
- faults.odb
- materials.odb
+ verify containers work
Check trapping of errors in reading spatialdata files.
Modified: short/3D/PyLith/trunk/examples/templates/containers/bc.odb
===================================================================
--- short/3D/PyLith/trunk/examples/templates/containers/bc.odb 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/examples/templates/containers/bc.odb 2007-06-13 15:11:19 UTC (rev 7194)
@@ -20,7 +20,7 @@
## the boundary condition facility of the problem:
##
## [pylithapp.timedependent]
-## bc = pylith.bc.MyBC
+## bc = MyBC
##
## You can then set the properties of the facilities in this container
## just as you would any other Pyre component.
Added: short/3D/PyLith/trunk/examples/templates/containers/faults.odb
===================================================================
--- short/3D/PyLith/trunk/examples/templates/containers/faults.odb 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/examples/templates/containers/faults.odb 2007-06-13 15:11:19 UTC (rev 7194)
@@ -0,0 +1,131 @@
+#!/usr/bin/env python
+#
+# ----------------------------------------------------------------------
+#
+# Brad T. Aagaard
+# U.S. Geological Survey
+#
+# <LicenseText>
+#
+# ----------------------------------------------------------------------
+#
+
+## Template for user-defined faults container.
+##
+## Copy this file to your working directory and rename it as
+## appropriate. Note: the extension MUST rename ".odb" for PyLith to
+## find it.
+##
+## To use this container, in your .cfg file, bind this component to
+## the interfaces facility of the problem:
+##
+## [pylithapp.timedependent]
+## interfaces = MyFaults
+##
+## You can then set the properties of the facilities in this container
+## just as you would any other Pyre component.
+##
+## For example:
+##
+## [pylithapp.timedependent.interfaces.san_andreas]
+## label = San Andreas Fault
+## id = 2
+
+# Parent class
+from pylith.utils.ObjectBin import ObjectBin
+
+# Define new faults container class
+#
+# You can change the name of this class, but the name here MUST match
+# the one in the object_bin() function at the bottom of this file.
+class MyFaults(ObjectBin):
+ """
+ User-defined boundary conditions container.
+
+ Factory: object_bin
+ """
+
+ # INVENTORY //////////////////////////////////////////////////////////
+
+ class Inventory(ObjectBin.Inventory):
+ """
+ Python object for managing MyBC facilities and properties.
+ """
+
+ # Define the facilities in this faults container.
+ #
+ # You must import any Python objects that are bound to the
+ # facilities.
+ #
+ # Synopsis:
+ #
+ # facilityName = pyre.inventory.facility("facility_name",
+ # family="fault", factory=ClassNameOfComponent)
+ #
+ # where ClassNameOfComponent is the class name of the default
+ # component to bind to the facility 'facility_name'.
+
+ import pyre.inventory
+
+ from pylith.faults.FaultCohesiveKin import FaultCohesiveKin
+
+ sanandreas = pyre.inventory.facility("san_andreas", family="fault",
+ factory=CohesiveFaultKin)
+ sanandreas.meta['tip'] = "San Andreas fault."
+
+ sanjacinto = pyre.inventory.facility("san_jacinto", family="fault",
+ factory=CohesiveFaultKin)
+ sanjacinto.meta['tip'] = "San Jacinto fault."
+
+
+ # PUBLIC METHODS /////////////////////////////////////////////////////
+
+ # The 'name' argument on the next line defines the default name used
+ # to configure this component. Generally, we use the class name in
+ # lowercase.
+ def __init__(self, name="myfaults"):
+ """
+ Constructor.
+ """
+ ObjectBin.__init__(self, name)
+ return
+
+
+ # PRIVATE METHODS ////////////////////////////////////////////////////
+
+ def _configure(self):
+ """
+ Set attributes from inventory.
+ """
+ ObjectBin._configure(self)
+
+ # Define how the components are ordered in the container.
+ #
+ # Synopsis:
+ #
+ # self.bin = [self.inventory.facilityName1,
+ # self.inventory.facilityName2]
+ #
+ # Replace the names of the facilities as desired. The names MUST
+ # match the VARIABLES in the inventory, as opposed to the names
+ # used to bind the components to the facilities, which is the
+ # first arugment in the calls to pyre.inventory.facility()).
+ #
+ # DO NOT change the name 'self.bin'.
+ self.bin = [self.inventory.sanjacinto,
+ self.inventory.sanandreas]
+ return
+
+
+# FACTORIES ////////////////////////////////////////////////////////////
+
+def object_bin():
+ """
+ Factory associated with MyFaults.
+ """
+ # The class name for the container used above MUST match the name
+ # used on the next line.
+ return MyFaults()
+
+
+# End of file
Added: short/3D/PyLith/trunk/examples/templates/containers/materials.odb
===================================================================
--- short/3D/PyLith/trunk/examples/templates/containers/materials.odb 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/examples/templates/containers/materials.odb 2007-06-13 15:11:19 UTC (rev 7194)
@@ -0,0 +1,132 @@
+#!/usr/bin/env python
+#
+# ----------------------------------------------------------------------
+#
+# Brad T. Aagaard
+# U.S. Geological Survey
+#
+# <LicenseText>
+#
+# ----------------------------------------------------------------------
+#
+
+## Template for user-defined physical property container.
+##
+## Copy this file to your working directory and rename it as
+## appropriate. Note: the extension MUST rename ".odb" for PyLith to
+## find it.
+##
+## To use this container, in your .cfg file, bind this component to
+## the materials facility of the problem:
+##
+## [pylithapp.timedependent]
+## materials = MyMatBin
+##
+## You can then set the properties of the facilities in this container
+## just as you would any other Pyre component.
+##
+## For example:
+##
+## [pylithapp.timedepedent.materials.materialA]
+## label = elastic material
+## id = 1
+
+# Parent class
+from pylith.utils.ObjectBin import ObjectBin
+
+# Define new material container class
+#
+# You can change the name of this class, but the name here MUST match
+# the one in the object_bin() function at the bottom of this file.
+class MyMatBin(ObjectBin):
+ """
+ User-defined materials container.
+
+ Factory: object_bin
+ """
+
+ # INVENTORY //////////////////////////////////////////////////////////
+
+ class Inventory(ObjectBin.Inventory):
+ """
+ Python object for managing MyBC facilities and properties.
+ """
+
+ # Define the facilities in this material container.
+ #
+ # You must import any Python objects that are bound to the
+ # facilities.
+ #
+ # Synopsis:
+ #
+ # facilityName = pyre.inventory.facility("facility_name",
+ # family="material", factory=ClassNameOfComponent)
+ #
+ # where ClassNameOfComponent is the class name of the default
+ # component to bind to the facility 'facility_name'.
+
+ import pyre.inventory
+
+ from pylith.materials.ElasticIsotropic3D import ElasticIsotropic3D
+ from pylith.materials.MaxwellIsotropic3D import MaxwellIsotropic3D
+
+ elastic = pyre.inventory.facility("elastic", family="material",
+ factory=ElasticIsotropic3D)
+ elastic.meta['tip'] = "Elastic material."
+
+ viscoelastic = pyre.inventory.facility("viscoelastic", family="material",
+ factory=MaxwellIsotropic3D)
+ viscoelastic.meta['tip'] = "Viscoelastic material."
+
+
+ # PUBLIC METHODS /////////////////////////////////////////////////////
+
+ # The 'name' argument on the next line defines the default name used
+ # to configure this component. Generally, we use the class name in
+ # lowercase.
+ def __init__(self, name="mymatbin"):
+ """
+ Constructor.
+ """
+ ObjectBin.__init__(self, name)
+ return
+
+
+ # PRIVATE METHODS ////////////////////////////////////////////////////
+
+ def _configure(self):
+ """
+ Set attributes from inventory.
+ """
+ ObjectBin._configure(self)
+
+ # Define how the components are ordered in the container.
+ #
+ # Synopsis:
+ #
+ # self.bin = [self.inventory.facilityName1,
+ # self.inventory.facilityName2]
+ #
+ # Replace the names of the facilities as desired. The names MUST
+ # match the VARIABLES in the inventory, as opposed to the names
+ # used to bind the components to the facilities, which is the
+ # first arugment in the calls to pyre.inventory.facility()).
+ #
+ # DO NOT change the name 'self.bin'.
+ self.bin = [self.inventory.elastic,
+ self.inventory.viscoelastic]
+ return
+
+
+# FACTORIES ////////////////////////////////////////////////////////////
+
+def object_bin():
+ """
+ Factory associated with MyMatBin.
+ """
+ # The class name for the container used above MUST match the name
+ # used on the next line.
+ return MyMatBin()
+
+
+# End of file
Modified: short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh 2007-06-13 15:11:19 UTC (rev 7194)
@@ -159,6 +159,11 @@
/// Fault vertices associated with constraints
std::set<Mesh::point_type> _constraintVert;
+ /// Label of cell used to compute Jacobian for each vertex (must
+ /// prevent overlap so only allow 1 cell to contribute for each
+ /// vertex).
+ ALE::Obj<int_section_type> _vertexCells;
+
}; // class FaultCohesiveKin
#include "FaultCohesiveKin.icc" // inline methods
Modified: short/3D/PyLith/trunk/libsrc/meshio/SolutionIOVTK.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/meshio/SolutionIOVTK.cc 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/libsrc/meshio/SolutionIOVTK.cc 2007-06-13 15:11:19 UTC (rev 7194)
@@ -10,6 +10,10 @@
// ======================================================================
//
+/* For now, we write each time step to a different file, so we throw
+ * the whole implementation into writeField().
+ */
+
#include <portinfo>
#include "SolutionIOVTK.hh" // implementation of class methods
@@ -25,7 +29,7 @@
// ----------------------------------------------------------------------
// Constructor
pylith::meshio::SolutionIOVTK::SolutionIOVTK(void) :
- _viewer(NULL)
+ _viewer(0)
{ // constructor
} // constructor
@@ -33,7 +37,9 @@
// Destructor
pylith::meshio::SolutionIOVTK::~SolutionIOVTK(void)
{ // destructor
- if (_viewer) PetscViewerDestroy(_viewer); _viewer = NULL;
+ if (_viewer)
+ PetscViewerDestroy(_viewer);
+ _viewer = NULL;
} // destructor
// ----------------------------------------------------------------------
@@ -41,21 +47,21 @@
void
pylith::meshio::SolutionIOVTK::open(const ALE::Obj<ALE::Mesh>& mesh)
{ // open
- PetscErrorCode ierr;
+#if 0
+ PetscErrorCode err;
- ierr = PetscViewerCreate(mesh->comm(), &_viewer);
- ierr = PetscViewerSetType(_viewer, PETSC_VIEWER_ASCII);
- ierr = PetscViewerSetFormat(_viewer, PETSC_VIEWER_ASCII_VTK);
- ierr = PetscViewerFileSetName(_viewer, _filename.c_str());
- if (ierr) {
+ err = PetscViewerCreate(mesh->comm(), &_viewer);
+ err = PetscViewerSetType(_viewer, PETSC_VIEWER_ASCII);
+ err = PetscViewerSetFormat(_viewer, PETSC_VIEWER_ASCII_VTK);
+ err = PetscViewerFileSetName(_viewer, _filename.c_str());
+ if (err) {
std::ostringstream msg;
msg << "Could not open VTK file '" << _filename
<< "' for solution output.\n";
throw std::runtime_error(msg.str());
} // if
- // Write header
- // ADD STUFF HERE
+#endif
} // open
// ----------------------------------------------------------------------
@@ -63,7 +69,9 @@
void
pylith::meshio::SolutionIOVTK::close(void)
{ // close
- if (_viewer) PetscViewerDestroy(_viewer); _viewer = NULL;
+ if
+ (_viewer) PetscViewerDestroy(_viewer);
+ _viewer = NULL;
} // close
// ----------------------------------------------------------------------
@@ -72,15 +80,13 @@
pylith::meshio::SolutionIOVTK::writeTopology(const ALE::Obj<ALE::Mesh>& mesh,
const spatialdata::geocoords::CoordSys* csMesh)
{ // writeTopology
+ #if 0
try {
- PetscErrorCode ierr;
+ PetscErrorCode err = 0;
- ierr = VTKViewer::writeHeader(_viewer);
- ierr = VTKViewer::writeVertices(mesh, _viewer);
- ierr = VTKViewer::writeElements(mesh, _viewer);
- // Use spatialdata::geocoords::Converter::convert() to convert
- // coordinates of vertices from csMesh to _cs (postpone and wait
- // for more general implementation of SolutionIO?).
+ err = VTKViewer::writeHeader(_viewer);
+ err = VTKViewer::writeVertices(mesh, _viewer);
+ err = VTKViewer::writeElements(mesh, _viewer);
} catch (const std::exception& err) {
std::ostringstream msg;
msg << "Error while writing topology information to VTK file '"
@@ -92,6 +98,7 @@
<< _filename << "'.\n";
throw std::runtime_error(msg.str());
} // try/catch
+#endif
} // writeTopology
// ----------------------------------------------------------------------
@@ -103,12 +110,34 @@
const char* name,
const ALE::Obj<ALE::Mesh>& mesh)
{ // writeField
+
try {
- PetscErrorCode ierr;
+ PetscErrorCode err;
- // Ignore time for now
+ std::ostringstream buffer;
+ const int indexExt = _filename.find(".vtk");
+ buffer << std::string(_filename, 0, indexExt) << "_t" << t << ".vtk";
+
+ err = PetscViewerCreate(mesh->comm(), &_viewer);
+ err = PetscViewerSetType(_viewer, PETSC_VIEWER_ASCII);
+ err = PetscViewerSetFormat(_viewer, PETSC_VIEWER_ASCII_VTK);
+ err = PetscViewerFileSetName(_viewer, buffer.str().c_str());
+ if (err) {
+ std::ostringstream msg;
+ msg << "Could not open VTK file '" << buffer.str()
+ << "' for solution output.\n";
+ throw std::runtime_error(msg.str());
+ } // if
+
+ err = VTKViewer::writeHeader(_viewer);
+ err = VTKViewer::writeVertices(mesh, _viewer);
+ err = VTKViewer::writeElements(mesh, _viewer);
+
+ buffer.clear();
+ buffer << name << "_t" << t << std::endl;
+
field->view("");
- ierr = SectionView_Sieve_Ascii(mesh, field, name, _viewer);
+ err = SectionView_Sieve_Ascii(mesh, field, buffer.str().c_str(), _viewer);
} catch (const std::exception& err) {
std::ostringstream msg;
msg << "Error while writing field '" << name << "' at time "
Modified: short/3D/PyLith/trunk/pylith/meshio/SolutionIO.py
===================================================================
--- short/3D/PyLith/trunk/pylith/meshio/SolutionIO.py 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/pylith/meshio/SolutionIO.py 2007-06-13 15:11:19 UTC (rev 7194)
@@ -55,7 +55,7 @@
dt = pyre.inventory.dimensional("time_step", default=1.0*s)
dt.meta['tip'] = "Time step between solution output."
- skip = pyre.inventory.int("skip", default=1,
+ skip = pyre.inventory.int("skip", default=0,
validator=pyre.inventory.greaterEqual(0))
skip.meta['tip'] = "Number of time steps to skip between solution output."
Modified: short/3D/PyLith/trunk/pylith/problems/TimeDependent.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/TimeDependent.py 2007-06-13 13:48:17 UTC (rev 7193)
+++ short/3D/PyLith/trunk/pylith/problems/TimeDependent.py 2007-06-13 15:11:19 UTC (rev 7194)
@@ -111,6 +111,10 @@
from pyre.units.time import second
t = 0.0*second
+
+ # Do stuff for initial state
+ self._poststep(t)
+
while t.value <= self.totalTime.value:
self._info.log("Main time loop, t=%s" % t)
More information about the cig-commits
mailing list