[cig-commits] r19267 - short/3D/PyLith/branches/v1.6-stable/libsrc/pylith/problems

brad at geodynamics.org brad at geodynamics.org
Mon Dec 5 13:13:35 PST 2011


Author: brad
Date: 2011-12-05 13:13:35 -0800 (Mon, 05 Dec 2011)
New Revision: 19267

Modified:
   short/3D/PyLith/branches/v1.6-stable/libsrc/pylith/problems/SolverNonlinear.cc
Log:
Cleanup of nonlinear solver code. Remove redundant call to constrain solution space (already called when computing residual). Remove commented out code.

Modified: short/3D/PyLith/branches/v1.6-stable/libsrc/pylith/problems/SolverNonlinear.cc
===================================================================
--- short/3D/PyLith/branches/v1.6-stable/libsrc/pylith/problems/SolverNonlinear.cc	2011-12-05 21:12:13 UTC (rev 19266)
+++ short/3D/PyLith/branches/v1.6-stable/libsrc/pylith/problems/SolverNonlinear.cc	2011-12-05 21:13:35 UTC (rev 19267)
@@ -118,43 +118,6 @@
 { // solve
   assert(0 != solution);
 
-#if 0
-  // Update KSP operators with custom preconditioner if necessary.
-  const ALE::Obj<RealSection>& solutionSection = solution->section();
-  assert(!solutionSection.isNull());
-  const ALE::Obj<SieveMesh>& sieveMesh = solution->mesh().sieveMesh();
-   assert(!sieveMesh.isNull());
-  if (solutionSection->getNumSpaces() > sieveMesh->getDimension() &&
-      0 != _jacobianPCFault) {
-    PetscKSP ksp = 0;
-    PetscPC pc = 0;
-    PetscKSP *ksps = 0;
-    PetscMat A = 0;
-    PetscInt num = 0;
-
-    PetscErrorCode err = 0;
-    err = SNESGetKSP(_snes, &ksp); CHECK_PETSC_ERROR(err);
-    err = KSPSetUp(ksp); CHECK_PETSC_ERROR(err);
-    err = KSPGetPC(ksp, &pc); CHECK_PETSC_ERROR(err);
-    err = PCFieldSplitGetSubKSP(pc, &num, &ksps); CHECK_PETSC_ERROR(err);
-    assert(solutionSection->getNumSpaces() == num);
-
-#if 0 // debugging
-    std::cout << "Preconditioner Matrix" << std::endl;
-    MatView(_jacobianPCFault, PETSC_VIEWER_STDOUT_WORLD);
-#endif
-
-
-    MatStructure flag;
-    err = KSPGetOperators(ksps[num-1], &A, 
-			  PETSC_NULL, &flag); CHECK_PETSC_ERROR(err);
-    err = PetscObjectReference((PetscObject) A); CHECK_PETSC_ERROR(err);
-    err = KSPSetOperators(ksps[num-1], A, _jacobianPCFault, 
-			  flag); CHECK_PETSC_ERROR(err);
-    err = PetscFree(ksps); CHECK_PETSC_ERROR(err);
-  } // if
-#endif
-
   const int solveEvent = _logger->eventId("SoNl solve");
   const int scatterEvent = _logger->eventId("SoNl scatter");
   _logger->eventBegin(solveEvent);
@@ -365,7 +328,7 @@
       Formulation* formulation = (Formulation*) lsctx;
       assert(formulation);
       formulation->printState(&w, &g, &x, &y);
-#if 0 // TEMPORARY
+#if 0 // Original PETSc code
       *flag = PETSC_FALSE; // DIVERGED_LINE_SEARCH
 #else
       std::cerr << "WARNING: Line search diverged ... continuing nonlinear iterations anyway in hopes that solution will converge anyway."
@@ -382,24 +345,7 @@
     if (a == 0.0) {
       lambdatemp = -initslope/(2.0*b);
     } else {
-
-      // MATT: Check this
-      // 
-      // Temporary fix by Brad to keep lambda in proper
-      // range. Necessary due to underflow and overflow of a, b, c,
-      // and d.
-#if 0 
       lambdatemp = (-b + PetscSqrtScalar(d))/(3.0*a);
-#else
-      if ((-b + PetscSqrtScalar(d) > 0.0 && a > 0.0) ||
-	  (-b + PetscSqrtScalar(d) < 0.0 && a < 0.0)) {
-	lambdatemp = (-b + PetscSqrtScalar(d))/(3.0*a);
-      } else {
-	lambdatemp = 0.05*lambda;
-      } // else
-#endif   
-
- 
     } // if/else
     lambdaprev = lambda;
     gnormprev  = *gnorm;
@@ -410,8 +356,6 @@
     else
       lambda = lambdatemp;
 
-
-
     ierr  = VecWAXPY(w,-lambda,y,x);CHKERRQ(ierr);
     if (snes->nfuncs >= snes->max_funcs) {
       ierr = PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);CHKERRQ(ierr);
@@ -462,12 +406,6 @@
   ierr = PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);CHKERRQ(ierr);
 
   // ======================================================================
-  // Code to constrain solution space.
-  assert(lsctx);
-  Formulation* formulation = (Formulation*) lsctx;
-  assert(formulation);
-  formulation->constrainSolnSpace(&w);
-
   ierr = SNESComputeFunction(snes,w,g);CHKERRQ(ierr);
   if (snes->domainerror) {
     ierr = PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);CHKERRQ(ierr);



More information about the CIG-COMMITS mailing list