]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Removed: Rogue code that used SLES (no longer supported by deal.II)
authorToby D. Young <tyoung@ippt.pan.pl>
Fri, 23 Jul 2010 18:24:21 +0000 (18:24 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Fri, 23 Jul 2010 18:24:21 +0000 (18:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@21563 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/source/petsc_solver.cc

index adeb904bf4d494ce15073689d72f98c6cf4276fd..8f5af44aa60ddd089fb5c4e83b8701d7a00e4fc3 100644 (file)
@@ -69,68 +69,6 @@ namespace PETScWrappers
   {
     int ierr;
 
-                                     // PETSc unfortunately changed its
-                                     // interfaces completely and incompatibly
-                                     // between 2.1.6 and 2.2.0. This forces
-                                     // us to work around this by offering two
-                                     // completely different
-                                     // implementations. sigh...
-#if DEAL_II_PETSC_VERSION_LT(2,2,0)
-
-                                     // first create a solver object
-    SLES sles;
-    ierr = SLESCreate (mpi_communicator, &sles);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-                                     // set the matrices involved. the
-                                     // last argument is irrelevant here,
-                                     // since we use the solver only once
-                                     // anyway
-    ierr = SLESSetOperators (sles, A, preconditioner,
-                             SAME_PRECONDITIONER);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-                                     // let derived classes set the solver
-                                     // type, and the preconditioning object
-                                     // set the type of preconditioner
-    KSP ksp;
-    ierr = SLESGetKSP (sles, &ksp);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-    set_solver_type (ksp);
-
-    PC pc;
-    ierr = SLESGetPC (sles, &pc);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-    preconditioner.set_preconditioner_type (pc);
-
-                                     // then a convergence monitor
-                                     // function. that function simply checks
-                                     // with the solver_control object we have
-                                     // in this object for convergence
-    KSPSetConvergenceTest (ksp, &convergence_test,
-                           reinterpret_cast<void *>(&solver_control));
-
-                                     // then do the real work: set up solver
-                                     // internal data and solve the
-                                     // system. this could be joined in one
-                                     // operation, but it is recommended this
-                                     // way to be able to separate statistic
-                                     // output if requested
-    int iterations = 0;
-    ierr = SLESSetUp (sles, b, x);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-    ierr = SLESSolve (sles, b, x, &iterations);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-#else
-                                     // and here is the PETSc post 2.1.6
-                                     // version: SLES objects have been
-                                     // completely abandoned, so use KSP and
-                                     // PC only
-
                                      // first create a solver object if this
                                      // is necessary
     if (solver_data.get() == 0)
@@ -177,28 +115,13 @@ namespace PETScWrappers
 
                                      // then do the real work: set up solver
                                      // internal data and solve the
-                                     // system. unfortunately, the call
-                                     // sequence is different between PETSc
-                                     // versions 2.2.0 and 2.2.1
-#if DEAL_II_PETSC_VERSION_LT(2,2,1)
-    ierr = KSPSetRhs(solver_data->ksp,b);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-    ierr = KSPSetSolution(solver_data->ksp, x);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-    ierr = KSPSetUp (solver_data->ksp);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-    ierr = KSPSolve (solver_data->ksp);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-#  else
+                                     // system. 
     ierr = KSPSetUp (solver_data->ksp);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
     ierr = KSPSolve (solver_data->ksp, b, x);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
-#  endif
 
-#endif
                                      // do not destroy solver object
 //    solver_data.reset ();
 
@@ -211,7 +134,6 @@ namespace PETScWrappers
   }
 
 
-
   void
   SolverBase::reset()
   {
@@ -219,7 +141,6 @@ namespace PETScWrappers
   }
 
 
-
   SolverControl &
   SolverBase::control() const
   {
@@ -227,7 +148,6 @@ namespace PETScWrappers
   }
 
 
-
   int
   SolverBase::convergence_test (KSP                 /*ksp*/,
 #ifdef PETSC_USE_64BIT_INDICES

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.