[cig-commits] r7337 - in short/3D/PyLith/trunk: . libsrc/faults libsrc/feassemble libsrc/topology modulesrc/faults modulesrc/topology pylith/faults pylith/problems pylith/topology unittests/libtests/faults unittests/libtests/faults/data unittests/libtests/topology unittests/pytests/faults unittests/pytests/faults/data unittests/pytests/topology

brad at geodynamics.org brad at geodynamics.org
Wed Jun 20 19:39:27 PDT 2007


Author: brad
Date: 2007-06-20 19:39:26 -0700 (Wed, 20 Jun 2007)
New Revision: 7337

Added:
   short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_1d.spatialdb
   short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_2d.spatialdb
   short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_3d.spatialdb
   short/3D/PyLith/trunk/unittests/pytests/faults/data/bulkprops_2d.spatialdb
Modified:
   short/3D/PyLith/trunk/TODO
   short/3D/PyLith/trunk/libsrc/faults/Fault.hh
   short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.cc
   short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.hh
   short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.cc
   short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh
   short/3D/PyLith/trunk/libsrc/feassemble/ElasticityExplicit.cc
   short/3D/PyLith/trunk/libsrc/topology/FieldsManager.cc
   short/3D/PyLith/trunk/libsrc/topology/FieldsManager.hh
   short/3D/PyLith/trunk/modulesrc/faults/faults.pyxe.src
   short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src
   short/3D/PyLith/trunk/pylith/faults/Fault.py
   short/3D/PyLith/trunk/pylith/problems/Explicit.py
   short/3D/PyLith/trunk/pylith/problems/Formulation.py
   short/3D/PyLith/trunk/pylith/problems/Implicit.py
   short/3D/PyLith/trunk/pylith/topology/FieldsManager.py
   short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesive.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesiveKin.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.cc
   short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.hh
   short/3D/PyLith/trunk/unittests/libtests/faults/data/Makefile.am
   short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.cc
   short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.hh
   short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py
   short/3D/PyLith/trunk/unittests/pytests/faults/data/Makefile.am
   short/3D/PyLith/trunk/unittests/pytests/topology/TestFieldsManager.py
Log:
Fixed integrateResidual() for kinematic fault. Need to add in internal forces. Implemented scaling of Lagrange multiplier constraint entries to improve conditioning of Jacobian matrix. Cleaned up Implicit.py (moved integration of Jacobian to Formulation). Cleaned up time stamping some more. Reverted FieldsManager to have solutionField(). Removed some debugging output.

Modified: short/3D/PyLith/trunk/TODO
===================================================================
--- short/3D/PyLith/trunk/TODO	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/TODO	2007-06-21 02:39:26 UTC (rev 7337)
@@ -2,28 +2,19 @@
 CURRENT ISSUES
 ======================================================================
 
-  1. Fix scaling and time stepping for fault implementation.
+  2. Debug mesh distribution. [MATT]
 
-    BruneSlipFn
-      Add slipIncr()
-      [Needs unit test]
-    EqKinSrc
-      slipIncr()
-      [Needs unit test]
-    FaultCohesiveKin
-      Fix scaling.
-      C++ unit tests with solution increment
-        integrateJacobian
-        integrateResidual
+  3. Debug Maxwell viscoelastic model. [CHARLES]
 
-  2. Debug mesh distribution.
-     [MATT]
-
-Charles is debugging Maxwell viscoelastic model.
+  5. Add check to make sure every material supplied by user exists in
+  the mesh.
   
-  3. Create list of missing features (EqSim/PyLith 0.8) and new
-  features.
+  6. Need PETSc settings for absolute convergence (i.e., don't want to
+  max out iterations when there is no increment in displacement).
 
+----------------------------------------------------------------------
+List of missing features (EqSim/PyLith 0.8) and new features.
+
     Missing features
       EqSim
         Absorbing boundaries
@@ -38,6 +29,7 @@
         Multiple cell types
         1-D and 2-D simulations
         Importing of meshes from CUBIT
+	Exporting to VTK files
         Availability as open-source with documentation
       PyLith 0.8
         Dislocation-based fault implementation
@@ -48,6 +40,12 @@
 MAIN PRIORITIES (Brad)
 ======================================================================
 
+add check in fields manager when returning solution (solution name
+shouldn't be empty)
+
+add ability to have comments in meshio and spatial data files (if easy
+to do, would use '#' as delimiter).
+
 Create meshes for benchmarks
 
   strike-slip
@@ -57,9 +55,6 @@
     tet (LaGriT)
     hex (CUBIT)
 
-add ability to have comments in meshio and spatial data files (if easy
-to do, would use '#' as delimiter).
-
   FaultCohesiveKin
     multiple cohesive cells for hex8 mesh
 

Modified: short/3D/PyLith/trunk/libsrc/faults/Fault.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/Fault.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/faults/Fault.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -38,6 +38,10 @@
   namespace geocoords {
     class CoordSys;
   } // geocoords
+
+  namespace spatialdb {
+    class SpatialDB; // USES SpatialDB
+  } // spatialdb
 } // spatialdata
 
 /// C++ abstract base class for Fault object.
@@ -97,12 +101,15 @@
    * @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).
+   * @param matDB Database of bulk elastic properties for fault region
+   *   (used to improve conditioning of Jacobian matrix)
    */
   virtual
   void initialize(const ALE::Obj<ALE::Mesh>& mesh,
 		  const spatialdata::geocoords::CoordSys* cs,
 		  const double_array& upDir,
-		  const double_array& normalDir) = 0;
+		  const double_array& normalDir,
+		  spatialdata::spatialdb::SpatialDB* matDB) = 0;
 
   // PROTECTED METHODS //////////////////////////////////////////////////
 protected :

Modified: short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -40,7 +40,8 @@
 pylith::faults::FaultCohesiveDyn::initialize(const ALE::Obj<ALE::Mesh>& mesh,
 					     const spatialdata::geocoords::CoordSys* cs,
 					     const double_array& upDir,
-					     const double_array& normalDir)
+					     const double_array& normalDir,
+					     spatialdata::spatialdb::SpatialDB* matDB)
 { // initialize
   throw std::logic_error("FaultCohesiveDyn::initialize() not implemented.");
 } // initialize

Modified: short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveDyn.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -62,11 +62,14 @@
    * @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).
+   * @param matDB Database of bulk elastic properties for fault region
+   *   (used to improve conditioning of Jacobian matrix)
    */
   void initialize(const ALE::Obj<ALE::Mesh>& mesh,
 		  const spatialdata::geocoords::CoordSys* cs,
 		  const double_array& upDir,
-		  const double_array& normalDir);
+		  const double_array& normalDir,
+		  spatialdata::spatialdb::SpatialDB* matDB);
 
   /** Integrate contribution of cohesive cells to residual term.
    *

Modified: short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -23,6 +23,7 @@
 #include <petscmat.h> // USES PETSc Mat
 
 #include "spatialdata/geocoords/CoordSys.hh" // USES CoordSys
+#include "spatialdata/spatialdb/SpatialDB.hh" // USES CoordSys
 
 #include <Distribution.hh> // USES completeSection
 
@@ -59,7 +60,8 @@
 pylith::faults::FaultCohesiveKin::initialize(const ALE::Obj<ALE::Mesh>& mesh,
 					     const spatialdata::geocoords::CoordSys* cs,
 					     const double_array& upDir,
-					     const double_array& normalDir)
+					     const double_array& normalDir,
+					     spatialdata::spatialdb::SpatialDB* matDB)
 { // initialize
   assert(0 != _quadrature);
   assert(0 != _eqsrc);
@@ -312,9 +314,46 @@
       } // if
     } // for
   } // for
+
+  // Setup pseudo-stiffness of cohesive cells to improve conditioning
+  // of Jacobian matrix
+  _pseudoStiffness = new real_section_type(mesh->comm(), mesh->debug());
+  assert(!_pseudoStiffness.isNull());
+  for (std::set<Mesh::point_type>::const_iterator v_iter=vertConstraintBegin;
+       v_iter != vertConstraintEnd;
+       ++v_iter)
+    _pseudoStiffness->setFiberDimension(*v_iter, 1);
+  mesh->allocate(_pseudoStiffness);
+    
+  matDB->open();
+  const char* stiffnessVals[] = { "density", "vs" };
+  const int numStiffnessVals = 2;
+  matDB->queryVals(stiffnessVals, numStiffnessVals);
+  
+  double_array matprops(numStiffnessVals);
+  for (std::set<Mesh::point_type>::const_iterator v_iter=vertConstraintBegin;
+       v_iter != vertConstraintEnd;
+       ++v_iter) {
+    const real_section_type::value_type* vertexCoords = 
+      coordinates->restrictPoint(*v_iter);
+    int err = matDB->query(&matprops[0], numStiffnessVals, vertexCoords, 
+			   spaceDim, cs);
+    if (err) {
+      std::ostringstream msg;
+      msg << "Could not find material properties at (";
+      for (int i=0; i < spaceDim; ++i)
+	msg << "  " << vertexCoords[i];
+      msg << ") using spatial database " << matDB->label() << ".";
+      throw std::runtime_error(msg.str());
+    } // if
+    
+    const double density = matprops[0];
+    const double vs = matprops[1];
+    const double mu = density * vs*vs;
+    _pseudoStiffness->updatePoint(*v_iter, &mu);
+  } // for
 } // initialize
 
-const double fakeStiffness = 1.0;//e+9;
 // ----------------------------------------------------------------------
 // Integrate contribution of cohesive cells to residual term.
 void
@@ -328,7 +367,6 @@
   assert(0 != fields);
   assert(!mesh.isNull());
 
-#if 0
   // Cohesive cells with normal vertices i and j, and constraint
   // vertex k make 2 contributions to the residual:
   //
@@ -346,6 +384,7 @@
   double_array cellResidual(numCorners*spaceDim);
   double_array cellSoln(numCorners*spaceDim);
   double_array cellSlip(numConstraintVert*spaceDim);
+  double_array cellStiffness(numConstraintVert);
 
   // Get cohesive cells
   const ALE::Obj<ALE::Mesh::label_sequence>& cellsCohesive = 
@@ -357,12 +396,11 @@
     cellsCohesive->end();
 
   // Get section information
-  const ALE::Obj<real_section_type>& solution = fields->getReal("solution");
+  const ALE::Obj<real_section_type>& solution = fields->getSolution();
   assert(!solution.isNull());  
 
-  const bool _isSolnIncr = t > 0.0;
   ALE::Obj<real_section_type> slip;
-  if (!_isSolnIncr) {
+  if (!_useSolnIncr) {
     // Compute slip field at current time step
     assert(0 != _eqsrc);
     slip = _eqsrc->slip(t, _constraintVert);
@@ -374,9 +412,6 @@
     assert(!slip.isNull());
   } // else
   
-  solution->view("SOLUTION");
-  slip->view("SLIP");
-
   for (Mesh::label_sequence::iterator c_iter=cellsCohesiveBegin;
        c_iter != cellsCohesiveEnd;
        ++c_iter) {
@@ -389,6 +424,11 @@
     mesh->restrict(_orientation, *c_iter, &cellOrientation[0], 
 		   cellOrientation.size());
     
+    // Get pseudo stiffness at cells vertices (only valid at
+    // constraint vertices)
+    mesh->restrict(_pseudoStiffness, *c_iter, &cellStiffness[0], 
+		   cellStiffness.size());
+    
     // Get slip at cells vertices (only valid at constraint vertices)
     mesh->restrict(slip, *c_iter, &cellSlip[0], cellSlip.size());
 
@@ -406,11 +446,13 @@
       const int indexJ = iConstraint +   numConstraintVert;
       const int indexK = iConstraint + 2*numConstraintVert;
 
+      const double pseudoStiffness = cellStiffness[iConstraint];
+
       // Set slip values in residual vector; contributions are at DOF of
       // constraint vertices (k) of the cohesive cells
       for (int iDim=0; iDim < spaceDim; ++iDim)
 	cellResidual[indexK*spaceDim+iDim] = 
-	  cellSlip[iConstraint*spaceDim+iDim] * fakeStiffness;
+	  cellSlip[iConstraint*spaceDim+iDim] * pseudoStiffness;
       
       // Get orientation at constraint vertex
       const real_section_type::value_type* constraintOrient = 
@@ -421,41 +463,23 @@
 	for (int kDim=0; kDim < spaceDim; ++kDim)
 	  cellResidual[indexI*spaceDim+iDim] -=
 	    cellSoln[indexK*spaceDim+kDim] * 
-	    -constraintOrient[kDim*spaceDim+iDim];
+	    -constraintOrient[kDim*spaceDim+iDim] * pseudoStiffness;
 
       // Entries associated with constraint forces applied at node j
       for (int jDim=0; jDim < spaceDim; ++jDim)
 	for (int kDim=0; kDim < spaceDim; ++kDim)
 	  cellResidual[indexJ*spaceDim+jDim] -=
 	    cellSoln[indexK*spaceDim+kDim] * 
-	    constraintOrient[kDim*spaceDim+jDim];
+	    constraintOrient[kDim*spaceDim+jDim] * pseudoStiffness;
     } // for
+
     // Assemble cell contribution into field
     mesh->updateAdd(residual, *c_iter, &cellResidual[0]);
   } // for
-  PetscErrorCode err = PetscLogFlops(numConstraintVert*spaceDim*spaceDim*4 +
-				     numConstraintVert*spaceDim);
+  PetscErrorCode err = PetscLogFlops(numConstraintVert*spaceDim*spaceDim*6 +
+				     numConstraintVert*spaceDim*2);
   if (err)
     throw std::runtime_error("Logging PETSc flops failed.");
-#else
-  // Set slip values in residual vector; contributions are at DOF of
-  // constraint vertices (k) of the cohesive cells
-
-  typedef std::set<Mesh::point_type>::const_iterator vert_iterator;
-
-  assert(0 != _eqsrc);
-  const ALE::Obj<real_section_type>& slip = _eqsrc->slip(t, _constraintVert);
-  assert(!slip.isNull());
-  const vert_iterator vBegin = _constraintVert.begin();
-  const vert_iterator vEnd = _constraintVert.end();
-
-  for (vert_iterator v_iter=vBegin; v_iter != vEnd; ++v_iter) {
-    const int fiberDim = slip->getFiberDimension(*v_iter);
-    const real_section_type::value_type* values = slip->restrictPoint(*v_iter);
-    assert(fiberDim == residual->getFiberDimension(*v_iter));
-    residual->updatePoint(*v_iter, values);
-  } // for
-#endif
 } // integrateResidual
 
 // ----------------------------------------------------------------------
@@ -487,7 +511,7 @@
     cellsCohesive->end();
 
   // Get section information
-  const ALE::Obj<real_section_type>& solution = fields->getReal("solution");
+  const ALE::Obj<real_section_type>& solution = fields->getSolution();
   assert(!solution.isNull());  
 
   const int spaceDim = _quadrature->spaceDim();
@@ -498,6 +522,7 @@
   double_array cellMatrix(numCorners*spaceDim * numCorners*spaceDim);
   double_array cellOrientation(numConstraintVert*orientationSize);
   int_array cellConstraintCell(numConstraintVert);
+  double_array cellStiffness(numConstraintVert);
 
   for (Mesh::label_sequence::iterator c_iter=cellsCohesiveBegin;
        c_iter != cellsCohesiveEnd;
@@ -507,6 +532,11 @@
     mesh->restrict(_orientation, *c_iter, &cellOrientation[0], 
 		   cellOrientation.size());
 
+    // Get pseudo stiffness at cells vertices (only valid at
+    // constraint vertices)
+    mesh->restrict(_pseudoStiffness, *c_iter, &cellStiffness[0], 
+		   cellStiffness.size());
+    
     // Get constraint/cell pairings (only valid at constraint vertices)
     mesh->restrict(_constraintCell, *c_iter, &cellConstraintCell[0], 
 		   cellConstraintCell.size());
@@ -526,6 +556,8 @@
       const real_section_type::value_type* constraintOrient = 
 	&cellOrientation[iConstraint*orientationSize];
 
+      const double pseudoStiffness = cellStiffness[iConstraint];
+
       // Scale orientation information by pseudo-stiffness to bring
       // constraint forces in solution vector to the same order of
       // magnitude as the displacements to prevent ill-conditioning
@@ -536,7 +568,7 @@
 	  const int row = indexI*spaceDim+iDim;
 	  const int col = indexK*spaceDim+kDim;
 	  cellMatrix[row*numCorners*spaceDim+col] =
-	    -constraintOrient[kDim*spaceDim+iDim]*fakeStiffness;
+	    -constraintOrient[kDim*spaceDim+iDim]*pseudoStiffness;
 	  cellMatrix[col*numCorners*spaceDim+row] =
 	    cellMatrix[row*numCorners*spaceDim+col]; // symmetric
 	} // for
@@ -547,7 +579,7 @@
 	  const int row = indexJ*spaceDim+jDim;
 	  const int col = indexK*spaceDim+kDim;
 	  cellMatrix[row*numCorners*spaceDim+col] =
-	    constraintOrient[kDim*spaceDim+jDim]*fakeStiffness;
+	    constraintOrient[kDim*spaceDim+jDim]*pseudoStiffness;
 	  cellMatrix[col*numCorners*spaceDim+row] =
 	    cellMatrix[row*numCorners*spaceDim+col]; // symmetric
 	} // for

Modified: short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/faults/FaultCohesiveKin.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -21,9 +21,8 @@
  * fault. 
  *
  * The ordering of vertices in a cohesive cell is the vertices on the
- * POSITIVE/NEGATIVE (CHECK WHICH IT IS) side of the fault, the
- * corresponding entries on the other side of the fault, and then the
- * corresponding constraint vertices.
+ * one side of the fault, the corresponding entries on the other side
+ * of the fault, and then the corresponding constraint vertices.
  */
 
 #if !defined(pylith_faults_faultcohesivekin_hh)
@@ -76,11 +75,14 @@
    * @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).
+   * @param matDB Database of bulk elastic properties for fault region
+   *   (used to improve conditioning of Jacobian matrix)
    */
   void initialize(const ALE::Obj<ALE::Mesh>& mesh,
 		  const spatialdata::geocoords::CoordSys* cs,
 		  const double_array& upDir,
-		  const double_array& normalDir);
+		  const double_array& normalDir,
+		  spatialdata::spatialdb::SpatialDB* matDB);
 
   /** Integrate contributions to residual term (r) for operator.
    *
@@ -131,6 +133,10 @@
 
   EqKinSrc* _eqsrc; ///< Kinematic earthquake source information
 
+  /// Pseudo-stiffness for scaling constraint information to improve
+  /// conditioning of Jacobian matrix
+  ALE::Obj<real_section_type> _pseudoStiffness;
+
   /// Orientation of fault surface at vertices (fiber dimension is
   /// nonzero only at constraint vertices)
   ALE::Obj<real_section_type> _orientation;
@@ -143,6 +149,8 @@
   /// each vertex).
   ALE::Obj<int_section_type> _constraintCell;
 
+
+
 }; // class FaultCohesiveKin
 
 #include "FaultCohesiveKin.icc" // inline methods

Modified: short/3D/PyLith/trunk/libsrc/feassemble/ElasticityExplicit.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/ElasticityExplicit.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/feassemble/ElasticityExplicit.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -75,7 +75,6 @@
   assert(0 != _material);
   if (!_needNewJacobian)
     _needNewJacobian = _material->needNewJacobian();
-  std::cout << "ElasticityExplicit.cc Need new jacobian: " << _needNewJacobian << std::endl;
   return _needNewJacobian;
 } // needNewJacobian
 

Modified: short/3D/PyLith/trunk/libsrc/topology/FieldsManager.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/topology/FieldsManager.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/topology/FieldsManager.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -22,10 +22,10 @@
 // Constructor
 pylith::topology::FieldsManager::FieldsManager(const ALE::Obj<Mesh>& mesh) :
   _mesh(mesh),
-  _outputName("")
+  _solutionName("")
 { // constructor
 } // constructor
-
+ 
 // ----------------------------------------------------------------------
 // Destructor
 pylith::topology::FieldsManager::~FieldsManager(void)
@@ -208,27 +208,27 @@
 } // copyLayout
 
 // ----------------------------------------------------------------------
-// Set name of solution field for output
+// Set name of solution field.
 void
-pylith::topology::FieldsManager::outputField(const char* name)
-{ // outputField
+pylith::topology::FieldsManager::solutionField(const char* name)
+{ // solutionField
   map_real_type::const_iterator iter = _real.find(name);
   if (iter == _real.end()) {
     std::ostringstream msg;
     msg << "Cannot use unknown field '" << name 
-	<< "' when setting name of solution field for output.";
+	<< "' when setting name of solution field.";
     throw std::runtime_error(msg.str());
   } // if
-  _outputName = name;
-} // outputField
+  _solutionName = name;
+} // solutionField
 
 // ----------------------------------------------------------------------
-// Get solution field for output.
+// Get solution field.
 const ALE::Obj<pylith::real_section_type>&
-pylith::topology::FieldsManager::getOutputSoln(void)
-{ // getOutputSoln
-  return getReal(_outputName.c_str());
-} // getOutputSoln
+pylith::topology::FieldsManager::getSolution(void)
+{ // getSolution
+  return getReal(_solutionName.c_str());
+} // getSolution
 
 // ----------------------------------------------------------------------
 // Create history manager for a subset of the managed fields.

Modified: short/3D/PyLith/trunk/libsrc/topology/FieldsManager.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/topology/FieldsManager.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/libsrc/topology/FieldsManager.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -14,7 +14,7 @@
  * @file pylith/topology/FieldsManager.hh
  *
  * @brief Object for managing fields associated with the fields
- * defined over a finite-element mesh.
+ * defined over a finite-element mesh. 
  */
 
 #if !defined(pylith_topology_fieldsmanager_hh)
@@ -93,17 +93,17 @@
    */
   void copyLayout(const ALE::Obj<real_section_type>& field);
 
-  /** Set name of solution field for output.
+  /** Set name of solution field.
    *
-   * @param name Name of output field.
+   * @param name Name of field that is the solution.
    */
-  void outputField(const char* name);
+  void solutionField(const char* name);
 
-  /** Get solution field for output.
+  /** Get solution field.
    *
-   * @returns Solution field for output.
+   * @returns Solution field.
    */
-  const ALE::Obj<real_section_type>& getOutputSoln(void);
+  const ALE::Obj<real_section_type>& getSolution(void);
 
   /** Create history manager for a subset of the managed fields.
    *
@@ -149,8 +149,8 @@
   /// Map for fieldss stored as real fields
   map_real_type _real;
 
-  /// Name of field that is used to output the solution.
-  std::string _outputName;
+  /// Name of field that corresponds to the solution.
+  std::string _solutionName;
 
   /// History manager for a subset of the fields
   std::vector<std::string> _history;

Modified: short/3D/PyLith/trunk/modulesrc/faults/faults.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/faults/faults.pyxe.src	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/modulesrc/faults/faults.pyxe.src	2007-06-21 02:39:26 UTC (rev 7337)
@@ -119,19 +119,22 @@
     return
 
 
-  def initialize(self, mesh, cs, upDir, normalDir):
+  def initialize(self, mesh, cs, upDir, normalDir, matDB):
     """
     Initialize fault.
     """
     # create shim for method 'initialize'
-    #embed{ void Fault_initialize(void* objVptr, void* meshVptr, void* csVptr, double* upDirPtr, double* normalDirPtr)
+    #embed{ void Fault_initialize(void* objVptr, void* meshVptr, void* csVptr, double* upDirPtr, double* normalDirPtr, void* dbVptr)
     try {
       assert(0 != objVptr);
       assert(0 != meshVptr);
       assert(0 != csVptr);
+      assert(0 != dbVptr);
       ALE::Obj<ALE::Mesh>* mesh = (ALE::Obj<ALE::Mesh>*) meshVptr;
       spatialdata::geocoords::CoordSys* cs =
         (spatialdata::geocoords::CoordSys*) csVptr;
+      spatialdata::spatialdb::SpatialDB* matDB =
+        (spatialdata::spatialdb::SpatialDB*) dbVptr;
       pylith::double_array upDir(3);
       upDir[0] = upDirPtr[0];
       upDir[1] = upDirPtr[1];
@@ -140,7 +143,7 @@
       normalDir[0] = normalDirPtr[0];
       normalDir[1] = normalDirPtr[1];
       normalDir[2] = normalDirPtr[2];
-      ((pylith::faults::Fault*) objVptr)->initialize(*mesh, cs, upDir, normalDir);
+      ((pylith::faults::Fault*) objVptr)->initialize(*mesh, cs, upDir, normalDir, matDB);
       } catch (const std::exception& err) {
       PyErr_SetString(PyExc_RuntimeError,
                       const_cast<char*>(err.what()));
@@ -155,18 +158,22 @@
 
     if not mesh.name == "pylith_topology_Mesh":
       raise TypeError, \
-            "Argument must be extension module type " \
+            "Argument 'mesh' must be extension module type " \
             "'pylith::topology::Mesh'."
     if not cs.name == "spatialdata_geocoords_CoordSys":
       raise TypeError, \
-            "Argument must be extension module type " \
+            "Argument 'cs' must be extension module type " \
             "'spatialdata::geocoords::CoordSys'."
     if 3 != len(upDir):
       raise TypeError, \
-            "Argument must be a 3 vector (list)."
+            "Argument 'upDir' must be a 3 vector (list)."
     if 3 != len(normalDir):
       raise TypeError, \
-            "Argument must be a 3 vector (list)."
+            "Argument 'normalDir' must be a 3 vector (list)."
+    if not matDB.name == "spatialdata_spatialdb_SpatialDB":
+      raise TypeError, \
+            "Argument 'matDB' must be extension module type " \
+            "'spatialdata::spatialdb::SpatialDB'."
     cdef double upDirCpp[3]
     upDirCpp[0] = upDir[0]
     upDirCpp[1] = upDir[1]
@@ -176,7 +183,7 @@
     normalDirCpp[1] = normalDir[1]
     normalDirCpp[2] = normalDir[2]    
     Fault_initialize(self.thisptr, ptrFromHandle(mesh), ptrFromHandle(cs),
-                     upDirCpp, normalDirCpp)
+                     upDirCpp, normalDirCpp, ptrFromHandle(matDB))
     return
 
 

Modified: short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/modulesrc/topology/topology.pyxe.src	2007-06-21 02:39:26 UTC (rev 7337)
@@ -16,7 +16,7 @@
 #include "pylith/utils/petscfwd.h"
 #include <Generator.hh>
 #include <Distribution.hh>
-
+ 
 #include <stdexcept>
 #include <Python.h>
 #include <assert.h>
@@ -259,34 +259,8 @@
     fieldVptr = PyCObject_AsVoidPtr(field)
     ptr = Mesh_createMatrix(self.thisptr, fieldVptr)
     return PyCObject_FromVoidPtr(ptr, PetscMat_destructor)
-    
 
-  def view(self):
-    """
-    View the mesh.
-    """
-    # create shim for view
-    #embed{ void* Mesh_view(void* objVptr)
-    try {
-      ALE::Obj<ALE::Mesh>* mesh = (ALE::Obj<ALE::Mesh>*) objVptr;
-      assert(0 != mesh);
-      assert(!mesh->isNull());
-      (*mesh)->view("");
-    } 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
-    Mesh_view(self.thisptr)
-    return
 
-
   def _createHandle(self):
     """
     Wrap pointer to C++ object in PyCObject.
@@ -728,15 +702,15 @@
     return
 
 
-  def outputField(self, name):
+  def solutionField(self, name):
     """
-    Set name of field corresponding to solution field for output.
+    Set name of field corresponding to solution.
     """
-    # create shim for 'outputField'
-    #embed{ void FieldsManager_outputField(void* objVptr, char* name)
+    # create shim for solutionField
+    #embed{ void FieldsManager_solutionField(void* objVptr, char* name)
     try {
       assert(0 != objVptr);
-      ((pylith::topology::FieldsManager*) objVptr)->outputField(name);
+      ((pylith::topology::FieldsManager*) objVptr)->solutionField(name);
     } catch (const std::exception& err) {
       PyErr_SetString(PyExc_RuntimeError,
                       const_cast<char*>(err.what()));
@@ -748,20 +722,20 @@
                       "Caught unknown C++ exception.");
     } // try/catch
     #}embed
-    FieldsManager_outputField(self.thisptr, name)
+    FieldsManager_solutionField(self.thisptr, name)
 
 
-  def getOutputSoln(self):
+  def getSolution(self):
     """
-    Get field corresponding to solution for output.
+    Get field corresponding to solution.
     """
-    # create shim for 'getOutputSoln'
-    #embed{ void* FieldsManager_getOutputSoln(void* objVptr)
+    # create shim for getSolution
+    #embed{ void* FieldsManager_getSolution(void* objVptr)
     void* result = 0;
     try {
       assert(0 != objVptr);
       const ALE::Obj<pylith::real_section_type>& field =
-        ((pylith::topology::FieldsManager*) objVptr)->getOutputSoln();
+        ((pylith::topology::FieldsManager*) objVptr)->getSolution();
       result = (void*) &field;
     } catch (const std::exception& err) {
       PyErr_SetString(PyExc_RuntimeError,
@@ -776,7 +750,7 @@
     return result;
     #}embed
     cdef void* ptr
-    ptr = FieldsManager_getOutputSoln(self.thisptr)
+    ptr = FieldsManager_getSolution(self.thisptr)
     return PyCObject_FromVoidPtr(ptr, NULL)
 
 

Modified: short/3D/PyLith/trunk/pylith/faults/Fault.py
===================================================================
--- short/3D/PyLith/trunk/pylith/faults/Fault.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/pylith/faults/Fault.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -69,6 +69,8 @@
     ## @li \b normal_dir 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).
+    ## @li \b mat_db Spatial database for bulk material properties
+    ##   (used in improving conditioning of Jacobian matrix).
     ##
     ## \b Facilities
     ## @li \b quadrature Quadrature object for numerical integration
@@ -99,7 +101,13 @@
     quadrature = pyre.inventory.facility("quadrature", factory=Quadrature)
     quadrature.meta['tip'] = "Quadrature object for numerical integration."
 
+    from spatialdata.spatialdb.SimpleDB import SimpleDB
+    matDB = pyre.inventory.facility("mat_db", family="spatial_database",
+                                   factory=SimpleDB, args=["bulk materials"])
+    matDB.meta['tip'] = "Spatial database for bulk material properties " \
+                        "(used in improving conditioning of Jacobian matrix)."
 
+
   # PUBLIC METHODS /////////////////////////////////////////////////////
 
   def __init__(self, name="fault"):
@@ -128,10 +136,10 @@
     """
     Initialize fault.
     """
-    self._info.log("Initializing fault '%s'." % self.label)
     self._createCppHandle()
     
     self.quadrature.initialize()
+    self.matDB.initialize()
 
     faultDim = mesh.dimension() - 1
     if faultDim != self.quadrature.cell.cellDim:
@@ -145,7 +153,8 @@
     self.cppHandle.label = self.label
     self.cppHandle.quadrature = self.quadrature.cppHandle
     self.cppHandle.initialize(mesh.cppHandle, mesh.coordsys.cppHandle,
-                              self.upDir, self.normalDir)
+                              self.upDir, self.normalDir,
+                              self.matDB.cppHandle)
     return
 
 
@@ -161,6 +170,7 @@
     self.upDir = self.inventory.upDir
     self.normalDir = self.inventory.normalDir
     self.quadrature = self.inventory.quadrature
+    self.matDB = self.inventory.matDB
     return
 
   

Modified: short/3D/PyLith/trunk/pylith/problems/Explicit.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/Explicit.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/pylith/problems/Explicit.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -45,8 +45,8 @@
     Constructor.
     """
     Formulation.__init__(self, name)
-    self.outputField = {'name': "dispT",
-                        'label': "displacements"}
+    self.solnField = {'name': "dispT",
+                      'label': "displacements"}
     return
 
 
@@ -69,14 +69,14 @@
                            interfaceConditions, dimension, dt)
 
     self._info.log("Creating other fields and matrices.")
-    self.fields.addReal("dispT")
+    self.fields.addReal("dispTpdt")
     self.fields.addReal("dispTmdt")
     self.fields.addReal("residual")
-    self.fields.createHistory(["solution", "dispT", "dispTmdt"])    
-    self.fields.copyLayout("solution")
-    self.jacobian = mesh.createMatrix(self.fields.getReal("solution"))
+    self.fields.createHistory(["dispTpdt", "dispT", "dispTmdt"])    
+    self.fields.copyLayout("dispT")
+    self.jacobian = mesh.createMatrix(self.fields.getSolution())
 
-    self.solver.initialize(mesh, self.fields.getReal("solution"))
+    self.solver.initialize(mesh, self.fields.getSolution())
 
     # Solve for total displacement field
     for constraint in self.constraints:
@@ -86,7 +86,7 @@
     return
 
 
-  def startTime(self):
+  def startTime(self, dt):
     """
     Get time at which time stepping should start.
     """
@@ -108,7 +108,7 @@
     """
     Hook for doing stuff before advancing time step.
     """
-    dispTpdt = self.fields.getReal("solution")
+    dispTpdt = self.fields.getReal("dispTpdt")
     for constraint in self.constraints:
       constraint.setField(t+dt, dispTpdt)
 
@@ -140,11 +140,11 @@
       integrator.integrateResidual(residual, t, self.fields)
 
     self._info.log("Solving equations.")
-    self.solver.solve(self.fields.getReal("solution"), self.jacobian, residual)
+    self.solver.solve(self.fields.getReal("dispTpdt"), self.jacobian, residual)
     return
 
 
-  def poststep(self, t):
+  def poststep(self, t, dt):
     """
     Hook for doing stuff after advancing time step.
     """
@@ -154,7 +154,7 @@
     for integrator in self.integrators:
       integrator.updateState(t, self.fields.getReal("dispT"))
 
-    Formulation.poststep(self, t)
+    Formulation.poststep(self, t, dt)
     return
 
 

Modified: short/3D/PyLith/trunk/pylith/problems/Formulation.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/Formulation.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/pylith/problems/Formulation.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -67,7 +67,7 @@
     self.integrators = None
     self.constraints = None
     self.fields = None
-    self.outputField = None
+    self.solnField = None
     self._istep = 0
     return
 
@@ -148,23 +148,25 @@
       output.writeTopology()
 
     self._info.log("Creating solution field.")
-    self.fields.addReal("solution")
-    self.fields.setFiberDimension("solution", dimension)
+    solnName = self.solnField['name']
+    self.fields.addReal(solnName)
+    self.fields.solutionField(solnName)
+    self.fields.setFiberDimension(solnName, dimension)
     for constraint in self.constraints:
-      constraint.setConstraintSizes(self.fields.getReal("solution"))
-    self.fields.allocate("solution")
+      constraint.setConstraintSizes(self.fields.getSolution())
+    self.fields.allocate(solnName)
     for constraint in self.constraints:
-      constraint.setConstraints(self.fields.getReal("solution"))
+      constraint.setConstraints(self.fields.getSolution())
     return
 
 
-  def poststep(self, t):
+  def poststep(self, t, dt):
     """
     Hook for doing stuff after advancing time step.
     """
-    field = self.fields.getReal(self.outputField['name'])
+    field = self.fields.getSolution()
     for output in self.output.bin:
-      output.writeField(t, self._istep, field, self.outputField['label'])
+      output.writeField(t+dt, self._istep, field, self.solnField['label'])
     self._istep += 1
     return
 
@@ -194,6 +196,20 @@
     return
 
 
+  def _reformJacobian(self, t, dt):
+    """
+    Reform Jacobian matrix for operator.
+    """
+    self._info.log("Reforming Jacobian of operator.")
+    import pylith.utils.petsc as petsc
+    petsc.mat_setzero(self.jacobian)
+    for integrator in self.integrators:
+      integrator.timeStep(dt)
+      integrator.integrateJacobian(self.jacobian, t+dt, self.fields)
+    petsc.mat_assemble(self.jacobian)
+    return
+
+
 # FACTORIES ////////////////////////////////////////////////////////////
 
 def pde_formulation():

Modified: short/3D/PyLith/trunk/pylith/problems/Implicit.py
===================================================================
--- short/3D/PyLith/trunk/pylith/problems/Implicit.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/pylith/problems/Implicit.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -68,8 +68,8 @@
     Constructor.
     """
     Formulation.__init__(self, name)
-    self.outputField = {'name': "dispTBctpdt",
-                        'label': "displacements"}
+    self.solnField = {'name': "dispTBctpdt",
+                      'label': "displacements"}
     return
 
 
@@ -90,11 +90,11 @@
                            interfaceConditions, dimension, dt)
 
     self._info.log("Creating other fields and matrices.")
-    self.fields.addReal("dispTBctpdt")
+    self.fields.addReal("dispIncr")
     self.fields.addReal("residual")
-    self.fields.copyLayout("solution")
-    self.jacobian = mesh.createMatrix(self.fields.getReal("solution"))
-    self.solver.initialize(mesh, self.fields.getReal("solution"))
+    self.fields.copyLayout("dispTBctpdt")
+    self.jacobian = mesh.createMatrix(self.fields.getSolution())
+    self.solver.initialize(mesh, self.fields.getSolution())
 
     # Initial time step solves for total displacement field, not increment
     for constraint in self.constraints:
@@ -137,13 +137,7 @@
       if integrator.needNewJacobian():
         needNewJacobian = True
     if needNewJacobian:
-      self._info.log("Reforming Jacobian of operator.")
-      import pylith.utils.petsc as petsc
-      petsc.mat_setzero(self.jacobian)
-      for integrator in self.integrators:
-        integrator.timeStep(dt)
-        integrator.integrateJacobian(self.jacobian, t+dt, self.fields)
-      petsc.mat_assemble(self.jacobian)
+      self._reformJacobian()
     return
 
 
@@ -153,7 +147,7 @@
     """
     self._info.log("Integrating residual term in operator.")
     residual = self.fields.getReal("residual")
-    dispIncr = self.fields.getReal("solution")
+    dispIncr = self.fields.getReal("dispIncr")
     import pylith.topology.topology as bindings
     bindings.zeroRealSection(residual)
     bindings.zeroRealSection(dispIncr)
@@ -177,9 +171,9 @@
     # nonzero at unconstrained DOF) so that after poststep(),
     # dispTBctpdt contains the displacement field at time t+dt.
     import pylith.topology.topology as bindings
-    solution = self.fields.getReal("solution")
-    disp = self.fields.getReal("dispTBctpdt")
-    bindings.addRealSections(disp, disp, solution)
+    dispIncr = self.fields.getReal("dispIncr")
+    disp = self.fields.getSolution()
+    bindings.addRealSections(disp, disp, dispIncr)
 
     self._info.log("Updating integrators states.")
     for integrator in self.integrators:
@@ -194,8 +188,9 @@
         constraint.useSolnIncr(True)
       for integrator in self.integrators:
         integrator.useSolnIncr(True)
+      self._reformJacobian(t, dt)
 
-    Formulation.poststep(self, t+dt)
+    Formulation.poststep(self, t, dt)
     return
 
 

Modified: short/3D/PyLith/trunk/pylith/topology/FieldsManager.py
===================================================================
--- short/3D/PyLith/trunk/pylith/topology/FieldsManager.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/pylith/topology/FieldsManager.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -11,7 +11,7 @@
 #
 
 ## @file pylith/topology/FieldsManager.py
-##
+## 
 ## @brief Python manager for fields over a mesh.
 
 # FieldsManager class
@@ -89,21 +89,21 @@
     return self.cppHandle.copyLayoutFromSrc(field)
 
 
-  def outputField(self, name):
+  def solutionField(self, name):
     """
-    Set name of field corresponding to solution for output.
+    Set name of field corresponding to solution.
     """
     assert(None != self.cppHandle)
-    self.cppHandle.outputField(name)
+    self.cppHandle.solutionField(name)
     return
 
 
-  def getOutputSoln(self):
+  def getSolution(self):
     """
-    Get field corresponding to solution for output.
+    Get field corresponding to solution.
     """
     assert(None != self.cppHandle)
-    return self.cppHandle.getOutputSoln()
+    return self.cppHandle.getSolution()
 
 
   def createHistory(self, labels):

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesive.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesive.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesive.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -432,7 +432,7 @@
   ALE::Obj<ALE::Mesh> mesh;
   meshio::MeshIOAscii iohandler;
   iohandler.filename(data.filename);
-  iohandler.debug(true);
+  iohandler.debug(false);
   iohandler.interpolate(false);
   iohandler.read(&mesh);
 
@@ -478,7 +478,7 @@
 
   int iCell = 0;
   i = 0;
-  mesh->view(data.filename);
+  //mesh->view(data.filename);
   for(Mesh::label_sequence::iterator c_iter = cells->begin();
       c_iter != cells->end();
       ++c_iter) {

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesiveKin.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesiveKin.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/TestFaultCohesiveKin.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -163,6 +163,20 @@
     CPPUNIT_ASSERT_EQUAL(_data->constraintCells[iVertex], vertexCell[0]);
   } // for
 
+  // Check pseudoStiffness
+  iVertex = 0;
+  for (std::set<Mesh::point_type>::const_iterator v_iter=vertConstraintBegin;
+       v_iter != vertConstraintEnd;
+       ++v_iter, ++iVertex) {
+    const int fiberDim = fault._pseudoStiffness->getFiberDimension(*v_iter);
+    CPPUNIT_ASSERT_EQUAL(1, fiberDim);
+    const real_section_type::value_type* vertexStiffness = 
+      fault._pseudoStiffness->restrictPoint(*v_iter);
+    CPPUNIT_ASSERT(0 != vertexStiffness);
+    const double tolerance = 1.0e-06;
+    CPPUNIT_ASSERT_DOUBLES_EQUAL(_data->pseudoStiffness, vertexStiffness[0],
+				 tolerance);
+  } // for
 } // testInitialize
 
 // ----------------------------------------------------------------------
@@ -178,6 +192,7 @@
   topology::FieldsManager fields(mesh);
   fields.addReal("residual");
   fields.addReal("solution");
+  fields.solutionField("solution");
   
   const ALE::Obj<real_section_type>& residual = fields.getReal("residual");
   CPPUNIT_ASSERT(!residual.isNull());
@@ -221,10 +236,11 @@
       residual->restrictPoint(*v_iter);
     for (int i=0; i < fiberDimE; ++i) {
       const int index = iVertex*spaceDim+i;
-      if (valsE[index] > tolerance)
-	CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, vals[i]/valsE[index], tolerance);
+      const double valE = valsE[index] * _data->pseudoStiffness;
+      if (valE > tolerance)
+	CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, vals[i]/valE, tolerance);
       else
-	CPPUNIT_ASSERT_DOUBLES_EQUAL(valsE[index], vals[i], tolerance);
+	CPPUNIT_ASSERT_DOUBLES_EQUAL(valE, vals[i], tolerance);
     } // for
   } // for
 } // testIntegrateResidual
@@ -242,6 +258,7 @@
   topology::FieldsManager fields(mesh);
   fields.addReal("residual");
   fields.addReal("solution");
+  fields.solutionField("solution");
   
   const ALE::Obj<real_section_type>& residual = fields.getReal("residual");
   CPPUNIT_ASSERT(!residual.isNull());
@@ -307,17 +324,18 @@
   for (int iRow=0; iRow < nrows; ++iRow)
     for (int iCol=0; iCol < ncols; ++iCol) {
       const int index = ncols*iRow+iCol;
+      const double valE = valsE[index] * _data->pseudoStiffness;
 #if 0 // DEBUGGING
-      if (fabs(valsE[index]-vals[index]) > tolerance)
+      if (fabs(valE-vals[index]) > tolerance)
 	std::cout << "ERROR: iRow: " << iRow << ", iCol: " << iCol
-		  << "valE: " << valsE[index]
+		  << "valE: " << valE
 		  << ", val: " << vals[index]
 		  << std::endl;
 #endif // DEBUGGING
-      if (fabs(valsE[index]) > 1.0)
-	CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, vals[index]/valsE[index], tolerance);
+      if (fabs(valE) > 1.0)
+	CPPUNIT_ASSERT_DOUBLES_EQUAL(1.0, vals[index]/valE, tolerance);
       else
-	CPPUNIT_ASSERT_DOUBLES_EQUAL(valsE[index], vals[index], tolerance);
+	CPPUNIT_ASSERT_DOUBLES_EQUAL(valE, vals[index], tolerance);
     } // for
   MatDestroy(jDense);
   MatDestroy(jSparseAIJ);
@@ -387,7 +405,12 @@
     const double normalDirVals[] = { 1.0, 0.0, 0.0 };
     double_array normalDir(normalDirVals, 3);
 
-    fault->initialize(*mesh, &cs, upDir, normalDir); 
+    spatialdata::spatialdb::SimpleDB dbMatProp("material properties");
+    spatialdata::spatialdb::SimpleIOAscii ioMatProp;
+    ioMatProp.filename(_data->matPropsFilename);
+    dbMatProp.ioHandler(&ioMatProp);
+
+    fault->initialize(*mesh, &cs, upDir, normalDir, &dbMatProp); 
   } catch (const ALE::Exception& err) {
     throw std::runtime_error(err.msg());
   } // catch

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -30,12 +30,14 @@
   finalSlipFilename(0),
   slipTimeFilename(0),
   peakRateFilename(0),
+  matPropsFilename(0),
   fieldT(0),
   orientation(0),
   constraintVertices(0),
   constraintCells(0),
   valsResidual(0),
   valsJacobian(0),
+  pseudoStiffness(0),
   numConstraintVert(0)
 { // constructor
 } // constructor

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinData.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -56,6 +56,7 @@
   char* finalSlipFilename; ///< Name of db for final slip
   char* slipTimeFilename; ///< Name of db for slip time
   char* peakRateFilename; ///< Name of db for peak rate
+  char* matPropsFilename; ///< Name of db for bulk material properties
   //@}
 
   /// @name Input fields
@@ -70,6 +71,7 @@
   int* constraintCells; ///< Expected cells for constraint vertices
   double* valsResidual; ///< Expected values from residual calculation.
   double* valsJacobian; ///< Expected values from Jacobian calculation.
+  double pseudoStiffness; ///< Fake stiffness for conditioning
   int numConstraintVert; ///< Number of constraint vertices
   //@}
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -101,6 +101,9 @@
 const char* pylith::faults::CohesiveKinDataHex8::_peakRateFilename = 
   "data/hex8_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataHex8::_matPropsFilename = 
+  "data/bulkprops_3d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataHex8::_fieldT[] = {
   4.1, 6.1, 8.1,
   4.2, 6.2, 8.2,
@@ -115,13 +118,13 @@
   5.1, 7.1, 9.1,
   5.2, 7.2, 9.2,
   5.3, 7.3, 9.3,
-  5.4, 7.4, 9.4,
+  5.4, 7.4, 9.4, // 15
   5.5, 7.5, 9.5,
-  5.6, 7.6, 9.6,
+  5.6, 7.6, 9.6, // 17
   5.7, 7.7, 9.7,
-  5.8, 7.8, 9.8,
+  5.8, 7.8, 9.8, // 19
   5.9, 7.9, 9.9,
-  5.0, 7.0, 9.0,
+  5.0, 7.0, 9.0, // 21
 };
 
 const int pylith::faults::CohesiveKinDataHex8::_numConstraintVert = 4;
@@ -146,21 +149,21 @@
   0.0, 0.0, 0.0,
   0.0, 0.0, 0.0,
   0.0, 0.0, 0.0,
+  9.4, 5.4, 7.4, // 6
+  9.6, 5.6, 7.6, // 7
+  9.8, 5.8, 7.8, // 8
+  9.0, 5.0, 7.0, // 9
   0.0, 0.0, 0.0,
   0.0, 0.0, 0.0,
   0.0, 0.0, 0.0,
   0.0, 0.0, 0.0,
-  0.0, 0.0, 0.0,
-  0.0, 0.0, 0.0,
-  0.0, 0.0, 0.0,
-  0.0, 0.0, 0.0,
-  0.0, 0.0, 0.0,
+ -9.4,-5.4,-7.4, // 14
   1.07974939836, -0.32861938211, 0.04694562602, // 15 (constraint)
-  0.0, 0.0, 0.0,
+ -9.6,-5.6,-7.6, // 16
   1.00381374723, -0.33460458241, 0.08365114560, // 17 (constraint)
-  0.0, 0.0, 0.0,
+ -9.8,-5.8,-7.8, // 18
   0.90493237602, -0.32577565537, 0.10859188512, // 19 (constraint)
-  0.0, 0.0, 0.0,
+ -9.0,-5.0,-7.0, // 20
   0.78469841324, -0.30180708202, 0.12072283281, // 21 (constraint)
 };
 
@@ -1367,6 +1370,8 @@
   0.0, 0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataHex8::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataHex8::CohesiveKinDataHex8(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -1384,12 +1389,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataHex8.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -66,6 +66,9 @@
 const char* pylith::faults::CohesiveKinDataLine2::_peakRateFilename = 
   "data/line2_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataLine2::_matPropsFilename = 
+  "data/bulkprops_1d.spatialdb";
+
 // Don't expect these values to be used, so just use some values.
 const double pylith::faults::CohesiveKinDataLine2::_fieldT[] = {
   7.1,
@@ -90,11 +93,11 @@
 };
 
 const double pylith::faults::CohesiveKinDataLine2::_valsResidual[] = {
-  0.0,
-  0.0,
-  0.0,
-  0.0,
-  1.05168389458
+   0.0,
+   7.5,
+   0.0,
+  -7.5,
+   1.05168389458
 };
 
 const double pylith::faults::CohesiveKinDataLine2::_valsJacobian[] = {
@@ -105,6 +108,8 @@
   0.0, -1.0,  0.0, +1.0,  0.0,
 };
 
+const double pylith::faults::CohesiveKinDataLine2::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataLine2::CohesiveKinDataLine2(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -122,12 +127,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataLine2.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -89,6 +89,9 @@
 const char* pylith::faults::CohesiveKinDataQuad4::_peakRateFilename = 
   "data/quad4_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataQuad4::_matPropsFilename = 
+  "data/bulkprops_2d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataQuad4::_fieldT[] = {
   8.1, 9.1,
   8.2, 9.2,
@@ -97,9 +100,9 @@
   8.5, 9.5,
   8.6, 9.6,
   8.7, 9.7,
-  8.8, 9.8,
+  8.8, 9.8, // 9
   8.9, 9.9,
-  8.0, 9.0,
+  8.0, 9.0, // 11
 };
 
 
@@ -121,13 +124,13 @@
 const double pylith::faults::CohesiveKinDataQuad4::_valsResidual[] = {
   0.0,  0.0,
   0.0,  0.0,
+  9.8,  8.8, // 4
+  9.0,  8.0, // 5
   0.0,  0.0,
   0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
+ -9.8, -8.8, // 8
   0.989535448086, 0.0824612873405, // 9
-  0.0,  0.0,
+  -9.0, -8.0, // 10
   1.05057813143, 0.0456773100622, // 11
 };
 
@@ -334,6 +337,8 @@
   0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataQuad4::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataQuad4::CohesiveKinDataQuad4(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -351,12 +356,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -108,6 +108,9 @@
 const char* pylith::faults::CohesiveKinDataQuad4e::_peakRateFilename = 
   "data/quad4e_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataQuad4e::_matPropsFilename = 
+  "data/bulkprops_2d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataQuad4e::_fieldT[] = {
   3.1, 5.1,
   3.2, 5.2,
@@ -119,11 +122,11 @@
   3.8, 5.8,
   3.9, 5.9,
   3.0, 5.0,
-  4.1, 6.1,
+  4.1, 6.1, // 14
   4.2, 6.2,
-  4.3, 6.3,
+  4.3, 6.3, // 16
   4.4, 6.4,
-  4.5, 6.5,
+  4.5, 6.5, // 18
 };
 
 
@@ -146,18 +149,18 @@
 const double pylith::faults::CohesiveKinDataQuad4e::_valsResidual[] = {
   0.0,  0.0,
   0.0,  0.0,
+  6.1,  4.1, // 6
+  6.3,  4.3, // 7
   0.0,  0.0,
   0.0,  0.0,
   0.0,  0.0,
+  6.5,  4.5, // 11
   0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
+ -6.1, -4.1, // 13
   0.989535448086, 0.0824612873405, // 14
-  0.0,  0.0,
+ -6.3, -4.3, // 15
   1.05057813143, 0.0456773100622, // 16
-  0.0,  0.0,
+ -6.5, -4.5, // 17
   0.90435792846,  0.10852295130, // 18
 };
 
@@ -614,6 +617,8 @@
   0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataQuad4e::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataQuad4e::CohesiveKinDataQuad4e(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -631,12 +636,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataQuad4e.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -77,6 +77,9 @@
 const char* pylith::faults::CohesiveKinDataTet4::_peakRateFilename = 
   "data/tet4_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataTet4::_matPropsFilename = 
+  "data/bulkprops_3d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataTet4::_fieldT[] = {
   7.1, 8.1, 9.1,
   7.2, 8.2, 9.2,
@@ -84,11 +87,11 @@
   7.4, 8.4, 9.4,
   7.5, 8.5, 9.5,
   7.6, 8.6, 9.6,
-  7.7, 8.7, 9.7,
+  7.7, 8.7, 9.7, // 8
   7.8, 8.8, 9.8,
-  7.9, 8.9, 9.9,
+  7.9, 8.9, 9.9, // 10
   7.0, 8.0, 9.0,
-  7.1, 8.1, 9.1,
+  7.1, 8.1, 9.1, // 12
 };
 
 const int pylith::faults::CohesiveKinDataTet4::_numConstraintVert = 3;
@@ -109,15 +112,15 @@
 
 const double pylith::faults::CohesiveKinDataTet4::_valsResidual[] = {
   0.0,  0.0,  0.0,
+  9.7,  7.7,  8.7, // 3
+  9.9,  7.9,  8.9, // 4
+  9.1,  7.1,  8.1, // 5
   0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
+ -9.7, -7.7, -8.7, // 7
   1.07974939836, -0.32861938211, 0.04694562602, // 8
-  0.0,  0.0,  0.0,
+ -9.9, -7.9, -8.9, // 4
   1.00381374723, -0.33460458241, 0.08365114560, // 10
-  0.0,  0.0,  0.0,
+ -9.1, -7.1, -8.1, // 5
   0.90493237602, -0.32577565537, 0.10859188512, // 12
 };
 
@@ -487,6 +490,8 @@
   0.0, 0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataTet4::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataTet4::CohesiveKinDataTet4(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -504,12 +509,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -77,6 +77,9 @@
 const char* pylith::faults::CohesiveKinDataTet4e::_peakRateFilename = 
   "data/tet4e_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataTet4e::_matPropsFilename = 
+  "data/bulkprops_3d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataTet4e::_fieldT[] = {
   3.1, 5.1, 7.1,
   3.2, 5.2, 7.2,
@@ -85,13 +88,13 @@
   3.5, 5.5, 7.5,
   3.6, 5.6, 7.6,
   3.7, 5.7, 7.7,
-  3.8, 5.8, 7.8,
+  3.8, 5.8, 7.8, // 11
   3.9, 5.9, 7.9,
-  3.0, 5.0, 7.0,
+  3.0, 5.0, 7.0, // 13
   4.1, 6.1, 8.1,
-  4.2, 6.2, 8.2,
+  4.2, 6.2, 8.2, // 15
   4.3, 6.3, 8.3,
-  4.4, 6.4, 8.4,
+  4.4, 6.4, 8.4, // 17
 };
 
 const int pylith::faults::CohesiveKinDataTet4e::_numConstraintVert = 4;
@@ -113,18 +116,18 @@
 
 const double pylith::faults::CohesiveKinDataTet4e::_valsResidual[] = {
   0.0,  0.0,  0.0,
+  7.8,  3.8,  5.8, // 5
+  7.0,  3.0,  5.0, // 6
+  8.2,  4.2,  6.2, // 7
+  8.4,  4.4,  6.4, // 8
   0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
+ -7.8, -3.8, -5.8, // 10
   1.07974939836, -0.32861938211, 0.04694562602, // 11
-  0.0,  0.0,  0.0,
+ -7.0, -3.0, -5.0, // 12
   1.00381374723, -0.33460458241, 0.08365114560, // 13
-  0.0,  0.0,  0.0,
+ -8.2, -4.2, -6.2, // 14
   0.90493237602, -0.32577565537, 0.10859188512, // 15
-  0.0, 0.0, 0.0,
+ -8.4, -4.4, -6.4, // 16
   0.78469841324, -0.30180708202, 0.12072283281, // 17
 };
 
@@ -719,6 +722,8 @@
   0.0, 0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataTet4e::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataTet4e::CohesiveKinDataTet4e(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -736,12 +741,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4e.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -77,6 +77,9 @@
 const char* pylith::faults::CohesiveKinDataTet4f::_peakRateFilename = 
   "data/tet4_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataTet4f::_matPropsFilename = 
+  "data/bulkprops_3d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataTet4f::_fieldT[] = {
   7.1, 8.1, 9.1,
   7.2, 8.2, 9.2,
@@ -84,11 +87,11 @@
   7.4, 8.4, 9.4,
   7.5, 8.5, 9.5,
   7.6, 8.6, 9.6,
-  7.7, 8.7, 9.7,
+  7.7, 8.7, 9.7, // 8
   7.8, 8.8, 9.8,
-  7.9, 8.9, 9.9,
+  7.9, 8.9, 9.9, // 10
   7.0, 8.0, 9.0,
-  7.1, 8.1, 9.1,
+  7.1, 8.1, 9.1, // 12
 };
 
 const int pylith::faults::CohesiveKinDataTet4f::_numConstraintVert = 3;
@@ -109,15 +112,15 @@
 
 const double pylith::faults::CohesiveKinDataTet4f::_valsResidual[] = {
   0.0,  0.0,  0.0,
+  9.7,  7.7,  8.7, // 3
+  9.9,  7.9,  8.9, // 4
+  9.1,  7.1,  8.1, // 5
   0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
-  0.0,  0.0,  0.0,
+ -9.7, -7.7, -8.7, // 7
   1.07974939836, -0.32861938211, 0.04694562602, // 8
-  0.0,  0.0,  0.0,
+ -9.9, -7.9, -8.9, // 4
   1.00381374723, -0.33460458241, 0.08365114560, // 10
-  0.0,  0.0,  0.0,
+ -9.1, -7.1, -8.1, // 5
   0.90493237602, -0.32577565537, 0.10859188512, // 12
 };
 
@@ -487,6 +490,8 @@
   0.0, 0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataTet4f::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataTet4f::CohesiveKinDataTet4f(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -504,12 +509,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTet4f.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -92,15 +92,18 @@
 const char* pylith::faults::CohesiveKinDataTri3::_peakRateFilename = 
   "data/tri3_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataTri3::_matPropsFilename = 
+  "data/bulkprops_2d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataTri3::_fieldT[] = {
   8.1, 9.1,
   8.2, 9.2,
   8.3, 9.3,
   8.4, 9.4,
   8.5, 9.5,
-  8.6, 9.6,
+  8.6, 9.6, // 7
   8.7, 9.7,
-  8.8, 9.8,
+  8.8, 9.8, // 9
 };
 
 const int pylith::faults::CohesiveKinDataTri3::_numConstraintVert = 2;
@@ -120,12 +123,12 @@
 
 const double pylith::faults::CohesiveKinDataTri3::_valsResidual[] = {
   0.0,  0.0,
+  9.6,  8.6, // 3
+  9.8,  8.8, // 4
   0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
+ -9.6, -8.6, // 6
   1.05057813143, 0.0456773100622, // 7
-  0.0,  0.0,
+ -9.8, -8.8, // 8
   0.989535448086, 0.0824612873405, // 9
 };
 
@@ -260,6 +263,8 @@
   0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataTri3::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataTri3::CohesiveKinDataTri3(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -277,12 +282,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -107,6 +107,9 @@
 const char* pylith::faults::CohesiveKinDataTri3d::_peakRateFilename = 
   "data/tri3d_peakrate.spatialdb";
 
+const char* pylith::faults::CohesiveKinDataTri3d::_matPropsFilename = 
+  "data/bulkprops_2d.spatialdb";
+
 const double pylith::faults::CohesiveKinDataTri3d::_fieldT[] = {
   6.1, 8.1,
   6.2, 8.2,
@@ -115,11 +118,11 @@
   6.5, 8.5,
   6.6, 8.6,
   6.7, 8.7,
-  6.8, 8.8,
+  6.8, 8.8, // 11
   6.9, 8.9,
-  6.0, 8.0,
+  6.0, 8.0, // 13
   7.1, 9.1,
-  7.2, 9.2,
+  7.2, 9.2, // 15
 };
 
 const int pylith::faults::CohesiveKinDataTri3d::_numConstraintVert = 3;
@@ -141,16 +144,16 @@
 
 const double pylith::faults::CohesiveKinDataTri3d::_valsResidual[] = {
   0.0,  0.0,
+  1.4142135623730949,  11.030865786510143, // 5
+  8.0,  6.0, // 6
   0.0,  0.0,
+ -7.2,  9.2, // 8
   0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
-  0.0,  0.0,
+ -1.4142135623730949, -11.030865786510143, // 10
   1.05057813143, 0.0456773100622, // 11
-  0.0,  0.0,
+ -8.0, -6.0, // 12
   0.989535448086, 0.0824612873405, // 13
-  0.0,  0.0,
+  7.2, -9.2, // 14
   0.90435792846,  0.10852295130, // 15
 };
 
@@ -445,6 +448,8 @@
   0.0, 0.0,
 };
 
+const double pylith::faults::CohesiveKinDataTri3d::_pseudoStiffness = 2.4;
+
 pylith::faults::CohesiveKinDataTri3d::CohesiveKinDataTri3d(void)
 { // constructor
   meshFilename = const_cast<char*>(_meshFilename);
@@ -462,12 +467,14 @@
   finalSlipFilename = const_cast<char*>(_finalSlipFilename);
   slipTimeFilename = const_cast<char*>(_slipTimeFilename);
   peakRateFilename = const_cast<char*>(_peakRateFilename);
+  matPropsFilename = const_cast<char*>(_matPropsFilename);
   fieldT = const_cast<double*>(_fieldT);
   orientation = const_cast<double*>(_orientation);
   constraintVertices = const_cast<int*>(_constraintVertices);
   constraintCells = const_cast<int*>(_constraintCells);
   valsResidual = const_cast<double*>(_valsResidual);
   valsJacobian = const_cast<double*>(_valsJacobian);
+  pseudoStiffness = _pseudoStiffness;
   numConstraintVert = _numConstraintVert;  
 } // constructor
 

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/CohesiveKinDataTri3d.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -54,6 +54,7 @@
   static const char* _finalSlipFilename; ///< Name of db for final slip
   static const char* _slipTimeFilename; ///< Name of db for slip time
   static const char* _peakRateFilename; ///< Name of db for peak rate
+  static const char* _matPropsFilename; ///< Name of db for bulk mat properties.
   //@}
 
   static const double _fieldT[]; ///< Solution field at time t.
@@ -63,6 +64,7 @@
   static const int _constraintCells[]; ///< Expected cells for constraint vertices
   static const double _valsResidual[]; ///< Expected values from residual calculation.
   static const double _valsJacobian[]; ///< Expected values from Jacobian calculation.
+  static const double _pseudoStiffness; ///< Fake stiffness for conditioning
   static const int _numConstraintVert; ///< Number of constraint vertices
 
 };

Modified: short/3D/PyLith/trunk/unittests/libtests/faults/data/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/Makefile.am	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/Makefile.am	2007-06-21 02:39:26 UTC (rev 7337)
@@ -11,6 +11,9 @@
 #
 
 noinst_DATA = \
+	bulkprops_1d.spatialdb \
+	bulkprops_2d.spatialdb \
+	bulkprops_3d.spatialdb \
 	line2.mesh \
 	line2_finalslip.spatialdb \
 	line2_sliptime.spatialdb \

Added: short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_1d.spatialdb
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_1d.spatialdb	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_1d.spatialdb	2007-06-21 02:39:26 UTC (rev 7337)
@@ -0,0 +1,14 @@
+#SPATIAL.ascii 1
+SimpleDB {
+  num-values = 2
+  value-names =  density  vs
+  value-units =  none  none
+  num-locs = 1
+  data-dim = 0
+  space-dim = 1
+  cs-data = cartesian {
+    to-meters = 1.0
+    space-dim = 1
+  }
+}
+0.0    2.4  1.0

Added: short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_2d.spatialdb
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_2d.spatialdb	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_2d.spatialdb	2007-06-21 02:39:26 UTC (rev 7337)
@@ -0,0 +1,14 @@
+#SPATIAL.ascii 1
+SimpleDB {
+  num-values = 2
+  value-names =  density  vs
+  value-units =  none  none
+  num-locs = 1
+  data-dim = 0
+  space-dim = 2
+  cs-data = cartesian {
+    to-meters = 1.0
+    space-dim = 2
+  }
+}
+0.0  0.0    2.4  1.0

Added: short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_3d.spatialdb
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_3d.spatialdb	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/faults/data/bulkprops_3d.spatialdb	2007-06-21 02:39:26 UTC (rev 7337)
@@ -0,0 +1,14 @@
+#SPATIAL.ascii 1
+SimpleDB {
+  num-values = 2
+  value-names =  density  vs
+  value-units =  none  none
+  num-locs = 1
+  data-dim = 0
+  space-dim = 3
+  cs-data = cartesian {
+    to-meters = 1.0
+    space-dim = 3
+  }
+}
+0.0  0.0  0.0    2.4  1.0

Modified: short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.cc
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.cc	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.cc	2007-06-21 02:39:26 UTC (rev 7337)
@@ -31,7 +31,7 @@
   _initialize(&mesh);
   FieldsManager manager(mesh);
 } // testConstructor
-
+ 
 // ----------------------------------------------------------------------
 // Test addReal().
 void
@@ -253,25 +253,25 @@
 } // testCopyLayoutFromField
 
 // ----------------------------------------------------------------------
-// Test outputField().
+// Test solutionField().
 void
-pylith::topology::TestFieldsManager::testOutputField(void)
-{ // testOutputField
+pylith::topology::TestFieldsManager::testSolutionField(void)
+{ // testSolutionField
   ALE::Obj<Mesh> mesh;
   _initialize(&mesh);
   FieldsManager manager(mesh);
 
   const std::string& name = "my solution";
   manager.addReal(name.c_str());
-  manager.outputField(name.c_str());
-  CPPUNIT_ASSERT_EQUAL(name, manager._outputName);
+  manager.solutionField(name.c_str());
+  CPPUNIT_ASSERT_EQUAL(name, manager._solutionName);
 } // testSolutionField
 
 // ----------------------------------------------------------------------
-// Test getOutputSoln().
+// Test getSolution().
 void
-pylith::topology::TestFieldsManager::testGetOutputSoln(void)
-{ // testGetOutputSoln
+pylith::topology::TestFieldsManager::testGetSolution(void)
+{ // testGetSolution
   ALE::Obj<Mesh> mesh;
   _initialize(&mesh);
   FieldsManager manager(mesh);
@@ -293,11 +293,11 @@
   fieldB->setFiberDimension(vertices, fiberDimB);
   fieldC->setFiberDimension(vertices, fiberDimC);
 
-  manager.outputField(labels[1]);
-  const ALE::Obj<real_section_type>& solution = manager.getOutputSoln();
+  manager.solutionField(labels[1]);
+  const ALE::Obj<real_section_type>& solution = manager.getSolution();
   CPPUNIT_ASSERT_EQUAL(fiberDimB, 
 		       solution->getFiberDimension(*(vertices->begin())));
-} // testGetOutputSoln
+} // testGetSolution
 
 // ----------------------------------------------------------------------
 // Test createHistory().

Modified: short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.hh
===================================================================
--- short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.hh	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/libtests/topology/TestFieldsManager.hh	2007-06-21 02:39:26 UTC (rev 7337)
@@ -14,7 +14,7 @@
  * @file unittests/libtests/topology/TestFieldsManager.hh
  *
  * @brief C++ TestFieldsManager object.
- *
+ * 
  * C++ unit testing for FieldsManager.
  */
 
@@ -48,8 +48,8 @@
   CPPUNIT_TEST( testAllocate );
   CPPUNIT_TEST( testCopyLayout );
   CPPUNIT_TEST( testCopyLayoutFromField );
-  CPPUNIT_TEST( testOutputField );
-  CPPUNIT_TEST( testGetOutputSoln );
+  CPPUNIT_TEST( testSolutionField );
+  CPPUNIT_TEST( testGetSolution );
   CPPUNIT_TEST( testCreateHistory );
   CPPUNIT_TEST( testShiftHistory );
   CPPUNIT_TEST( testGetHistoryItem );
@@ -82,11 +82,11 @@
   /// Test copyLayoutFromField().
   void testCopyLayoutFromField(void);
 
-  /// Test outputField().
-  void testOutputField(void);
+  /// Test solutionField().
+  void testSolutionField(void);
 
-  /// Test getOutputSoln().
-  void testGetOutputSoln(void);
+  /// Test getSolution().
+  void testGetSolution(void);
 
   /// Test createHistory().
   void testCreateHistory(void);

Modified: short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/pytests/faults/TestFaultCohesiveKin.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -265,6 +265,12 @@
     eqsrc = EqKinSrc()
     eqsrc.slipfn = slipfn
 
+    ioMatDB = SimpleIOAscii()
+    ioMatDB.filename = "data/bulkprops_2d.spatialdb"
+    dbMat = SimpleDB()
+    dbMat.iohandler = ioMatDB
+    dbMat.label = "bulk properties"
+    
     # Setup fault
     fault = FaultCohesiveKin()
     fault.id = 10
@@ -273,6 +279,7 @@
     fault.normalDir = [1, 0, 0]
     fault.quadrature = quadrature
     fault.eqsrc = eqsrc
+    fault.matDB = dbMat
     fault.timeStep(dt)
     fault.adjustTopology(mesh)
     fault.initialize(mesh)
@@ -283,6 +290,7 @@
     fields.addReal("residual")
     fields.addReal("solution")
     fields.addReal("disp")
+    fields.solutionField("solution")
     fields.setFiberDimension("residual", cs.spaceDim)
     fields.allocate("residual")
     fields.copyLayout("residual")

Modified: short/3D/PyLith/trunk/unittests/pytests/faults/data/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/faults/data/Makefile.am	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/pytests/faults/data/Makefile.am	2007-06-21 02:39:26 UTC (rev 7337)
@@ -11,6 +11,7 @@
 #
 
 noinst_DATA = \
+	bulkprops_2d.spatialdb \
 	tri3.mesh \
 	tri3_finalslip.spatialdb \
 	tri3_sliptime.spatialdb \

Added: short/3D/PyLith/trunk/unittests/pytests/faults/data/bulkprops_2d.spatialdb
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/faults/data/bulkprops_2d.spatialdb	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/pytests/faults/data/bulkprops_2d.spatialdb	2007-06-21 02:39:26 UTC (rev 7337)
@@ -0,0 +1,14 @@
+#SPATIAL.ascii 1
+SimpleDB {
+  num-values = 2
+  value-names =  density  vs
+  value-units =  none  none
+  num-locs = 1
+  data-dim = 0
+  space-dim = 2
+  cs-data = cartesian {
+    to-meters = 1.0
+    space-dim = 2
+  }
+}
+0.0  0.0    2.4  1.0

Modified: short/3D/PyLith/trunk/unittests/pytests/topology/TestFieldsManager.py
===================================================================
--- short/3D/PyLith/trunk/unittests/pytests/topology/TestFieldsManager.py	2007-06-21 01:36:57 UTC (rev 7336)
+++ short/3D/PyLith/trunk/unittests/pytests/topology/TestFieldsManager.py	2007-06-21 02:39:26 UTC (rev 7337)
@@ -11,7 +11,7 @@
 #
 
 ## @file unittests/pytests/topology/TestFieldsManager.py
-
+##
 ## @brief Unit testing of FieldsManager object.
 
 import unittest
@@ -179,11 +179,11 @@
     return
 
 
-  def test_outputField(self):
+  def test_solutionField(self):
     """
-    Test outputField().
+    Test solutionField().
 
-    WARNING: This is not a rigorous test of outputField() because we
+    WARNING: This is not a rigorous test of solutionField() because we
     don't verify the results.
     """
     mesh = self._initialize()
@@ -194,18 +194,18 @@
     for field in fields:
       manager.addReal(field)
 
-    manager.outputField("field B")
+    manager.solutionField("field B")
 
     # We should really add something here to check to make sure things
     # actually initialized correctly.
     return
 
 
-  def test_getOutputSoln(self):
+  def test_getSolution(self):
     """
-    Test getOutputSoln().
+    Test getSolution().
 
-    WARNING: This is not a rigorous test of getOutputSoln() because we
+    WARNING: This is not a rigorous test of getSolution() because we
     don't verify the results.
     """
     mesh = self._initialize()
@@ -216,8 +216,8 @@
     for field in fields:
       manager.addReal(field)
 
-    manager.outputField("field B")
-    solution = manager.getOutputSoln()
+    manager.solutionField("field B")
+    solution = manager.getSolution()
 
     # We should really add something here to check to make sure things
     # actually initialized correctly.



More information about the cig-commits mailing list