]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implemented a Solver class for the KSPPREONLY type of PETSc, which
authoroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 27 Jul 2004 14:09:25 +0000 (14:09 +0000)
committeroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 27 Jul 2004 14:09:25 +0000 (14:09 +0000)
just applies the preconditioner. This is useful when the complete LU
decomposition is used as preconditioner.

git-svn-id: https://svn.dealii.org/trunk@9522 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/petsc_solver.h
deal.II/lac/source/petsc_solver.cc

index 983d05115b93cace8d7354f556b80444a2ac73d5..f59d4857f03537b8f12f24f69787db95a2c92c82 100644 (file)
@@ -1018,6 +1018,77 @@ namespace PETScWrappers
       virtual void set_solver_type (KSP &ksp) const;
   };
   
+
+/**
+ * An implementation of the solver interface using the PETSc PREONLY
+ * solver. Actually this is NOT a real solution algorithm. Its only
+ * purpose is to provide a solver object, when the preconditioner
+ * should be used as real solver. It is very useful in conjunction with
+ * the complete LU decomposition preconditioner <tt> PreconditionLU </tt>, 
+ * which in conjunction with this solver class becomes a direct solver.
+ *
+ * @ingroup PETScWrappers
+ * @author Wolfgang Bangerth, 2004, Oliver Kayser-Herold, 2004
+ */
+  class SolverPreOnly : public SolverBase
+  {
+    public:
+                                       /**
+                                        * Standardized data struct to
+                                        * pipe additional data to the
+                                        * solver.
+                                        */
+      struct AdditionalData
+      {};
+       
+                                       /**
+                                        * Constructor. In contrast to
+                                        * deal.II's own solvers, there is no
+                                        * need to give a vector memory
+                                        * object. However, PETSc solvers want
+                                        * to have an MPI communicator context
+                                        * over which computations are
+                                        * parallelized. By default,
+                                        * @p PETSC_COMM_SELF is used here,
+                                        * but you can change this. Note that
+                                        * for single processor (non-MPI)
+                                        * versions, this parameter does not
+                                        * have any effect.
+                                        *
+                                        * The last argument takes a structure
+                                        * with additional, solver dependent
+                                        * flags for tuning.
+                                        *
+                                        * Note that the communicator used here
+                                        * must match the communicator used in
+                                        * the system matrix, solution, and
+                                        * right hand side object of the solve
+                                        * to be done with this
+                                        * solver. Otherwise, PETSc will
+                                        * generate hard to track down errors,
+                                        * see the documentation of the
+                                        * SolverBase class.
+                                        */
+      SolverPreOnly (SolverControl        &cn,
+                    const MPI_Comm       &mpi_communicator = PETSC_COMM_SELF,
+                    const AdditionalData &data = AdditionalData());
+      
+    protected:
+                                       /**
+                                        * Store a copy of the flags for this
+                                        * particular solver.
+                                        */
+      const AdditionalData additional_data;
+
+                                       /**
+                                        * Function that takes a Krylov
+                                        * Subspace Solver context object, and
+                                        * sets the type of solver that is
+                                        * appropriate for this class.
+                                        */
+      virtual void set_solver_type (KSP &ksp) const;
+  };
+
 }
 
 #endif // DEAL_II_USE_PETSC
index 70487b688b2c2131e6bdfd8f06c54507a1701e4a..949982d06f7f9cbb1d29d425aef72e58b4b5105b 100644 (file)
@@ -103,11 +103,6 @@ namespace PETScWrappers
     
     preconditioner.set_preconditioner_type (pc);
   
-                                     // in the deal.II solvers, we always
-                                     // honor the initial guess in the
-                                     // solution vector. do so here as well:
-    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
-  
                                      // then a convergence monitor
                                      // function. that function simply checks
                                      // with the solver_control object we have
@@ -162,12 +157,6 @@ namespace PETScWrappers
 
         preconditioner.set_preconditioner_type (solver_data->pc);
 
-                                         // in the deal.II solvers, we always
-                                         // honor the initial guess in the
-                                         // solution vector. do so here as
-                                         // well:
-        KSPSetInitialGuessNonzero (solver_data->ksp, PETSC_TRUE);
-
                                          // then a convergence monitor
                                          // function. that function simply
                                          // checks with the solver_control
@@ -285,6 +274,11 @@ namespace PETScWrappers
                                      // set the damping factor from the data
     ierr = KSPRichardsonSetScale (ksp, additional_data.omega);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
   
 
@@ -309,6 +303,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYCHEV));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
 
   
@@ -333,6 +332,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPCG));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
   
 
@@ -357,6 +361,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPBICG));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
 
   
@@ -428,6 +437,11 @@ namespace PETScWrappers
        ierr = KSPSetPreconditionerSide(ksp, PC_RIGHT);
        AssertThrow (ierr == 0, ExcPETScError(ierr));
       }
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
   
 
@@ -452,6 +466,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPBCGS));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
 
   
@@ -476,6 +495,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPCGS));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
   
 
@@ -500,6 +524,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPTFQMR));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
 
   
@@ -524,6 +553,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPTCQMR));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
 
   
@@ -548,6 +582,11 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPCR));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
   }
 
   
@@ -572,9 +611,52 @@ namespace PETScWrappers
     int ierr;
     ierr = KSPSetType (ksp, const_cast<char *>(KSPLSQR));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // in the deal.II solvers, we always
+                                     // honor the initial guess in the
+                                     // solution vector. do so here as well:
+    KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
+  }
+
+/* ---------------------- SolverPreOnly ------------------------ */
+
+  SolverPreOnly::SolverPreOnly (SolverControl        &cn,
+                               const MPI_Comm       &mpi_communicator,
+                               const AdditionalData &data)
+                  :
+                  SolverBase (cn, mpi_communicator),
+                  additional_data (data)
+  {}
+
+
+  void
+  SolverPreOnly::set_solver_type (KSP &ksp) const
+  {
+                                     // set the type of solver. work around a
+                                     // problem in PETSc 2.1.6, where it asks
+                                     // for a char*, even though KSPXXXX is of
+                                     // type const char*
+    int ierr;
+    ierr = KSPSetType (ksp, const_cast<char *>(KSPPREONLY));
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                     // The KSPPREONLY solver of
+                                     // PETSc never calls the convergence
+                                     // monitor, which leads to failure 
+                                     // even when everything was ok.
+                                     // Therefore the SolverControl status 
+                                     // is set to some nice values, which
+                                     // guarantee a nice result at the end 
+                                     // of the solution process.
+    solver_control.check (1, 0.0);
+
+                                     // Using the PREONLY solver with
+                                     // a nonzero initial guess leads
+                                     // PETSc to produce some error messages.
+    KSPSetInitialGuessNonzero (ksp, PETSC_FALSE);
   }
 
-  
 }
 
 #else

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.