[cig-commits] r20077 - in short/3D/PyLith/branches/pylith-scecdynrup: . libsrc/pylith/problems

[brad at geodynamics.org]

Author: brad
Date: 2012-05-10 15:15:39 -0700 (Thu, 10 May 2012)
New Revision: 20077

Modified:
   short/3D/PyLith/branches/pylith-scecdynrup/TODO
   short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/Solver.cc
   short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/SolverLinear.cc
Log:
Merge from trunk (v1.7-trunk).

Modified: short/3D/PyLith/branches/pylith-scecdynrup/TODO
===================================================================
--- short/3D/PyLith/branches/pylith-scecdynrup/TODO	2012-05-10 22:12:34 UTC (rev 20076)
+++ short/3D/PyLith/branches/pylith-scecdynrup/TODO	2012-05-10 22:15:39 UTC (rev 20077)
@@ -19,7 +19,7 @@
 
 * Add h5py to PyLith binary
 
-* Use h5py for automated full-scale tests.
+* Use HDF5/h5py in automated full-scale tests.
 
 * Setup fieldsplit for SolverNonlinear [Matt]
 

Modified: short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/Solver.cc
===================================================================
--- short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/Solver.cc	2012-05-10 22:12:34 UTC (rev 20076)
+++ short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/Solver.cc	2012-05-10 22:15:39 UTC (rev 20077)
@@ -167,15 +167,21 @@
 
   bool separateComponents = formulation->splitFieldComponents();
   if (separateComponents) {
+    PetscMat* precon = new PetscMat[numSpaces];
+    for (int i=0; i < numSpaces; ++i) {
+      precon[i] = PETSC_NULL;
+    } // for
+    precon[numSpaces-1] = _jacobianPCFault;
     constructFieldSplit(solutionSection, PETSC_DETERMINE, PETSC_NULL, PETSC_NULL, 
-			sieveMesh->getFactory()->getGlobalOrder(sieveMesh, "default", solutionSection), PETSC_NULL,
-			solution.vector(), *pc);
+			sieveMesh->getFactory()->getGlobalOrder(sieveMesh, "default", solutionSection), precon, PETSC_NULL, solution.vector(), *pc);
+    delete[] precon; precon = 0;
   } else {
     int numFields[2] = {spaceDim, (numSpaces > spaceDim) ? 1 : 0};
     MatNullSpace nullsp[2] = {PETSC_NULL, PETSC_NULL};
+    PetscMat precon[2] = {PETSC_NULL, _jacobianPCFault};
     int* fields = new int[numSpaces];
     
-    /* Create rigid body null space */
+    // Create rigid body null space.
     const ALE::Obj<RealSection>& coordinatesSection = sieveMesh->getRealSection("coordinates");
     assert(!coordinatesSection.isNull());
     const ALE::Obj<SieveMesh::label_sequence>& vertices = sieveMesh->depthStratum(0);
@@ -183,7 +189,7 @@
     const SieveMesh::label_sequence::iterator verticesBegin = vertices->begin();
     const SieveMesh::label_sequence::iterator verticesEnd = vertices->end();
     PetscInt dim = spaceDim;
-    Vec mode[6];
+    PetscVec mode[6];
 
     if (dim > 1) {
       PetscInt n;
@@ -241,14 +247,14 @@
       for(int i=0; i< m; ++i) {err = VecDestroy(&mode[i]);CHECK_PETSC_ERROR(err);}
     } // if
 
-    for (PetscInt f=0; f < numSpaces; ++f) {
+    for (int f=0; f < numSpaces; ++f) {
       fields[f] = f;
     } // for
     constructFieldSplit(solutionSection, 2, numFields, fields,
-                        sieveMesh->getFactory()->getGlobalOrder(sieveMesh, "default", solutionSection), nullsp,
+                        sieveMesh->getFactory()->getGlobalOrder(sieveMesh, "default", solutionSection), precon, nullsp,
                         solution.vector(), *pc);
     err = MatNullSpaceDestroy(&nullsp[0]);CHECK_PETSC_ERROR(err);
-    delete [] fields;
+    delete[] fields;
   } // if/else
 
   if (formulation->splitFields() && 

Modified: short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/SolverLinear.cc
===================================================================
--- short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/SolverLinear.cc	2012-05-10 22:12:34 UTC (rev 20076)
+++ short/3D/PyLith/branches/pylith-scecdynrup/libsrc/pylith/problems/SolverLinear.cc	2012-05-10 22:15:39 UTC (rev 20077)
@@ -125,7 +125,7 @@
   const ALE::Obj<SieveMesh>& sieveMesh = solution->mesh().sieveMesh();
   assert(!sieveMesh.isNull());
 
-#if 1 // :TODO: Does this get replaced by the use of the MatShell like SNES?
+#if 0 // OBSOLETE? Replaced by stuff in Solver?
   if (solutionSection->getNumSpaces() > sieveMesh->getDimension() &&
       _jacobianPCFault) {