]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
A patch by Vijay S. Mahadevan that allo0ws override of PETSc default solver parameters
authoryoung <young@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 22 Dec 2011 12:05:47 +0000 (12:05 +0000)
committeryoung <young@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 22 Dec 2011 12:05:47 +0000 (12:05 +0000)
git-svn-id: https://svn.dealii.org/trunk@24849 0785d39b-7218-0410-832d-ea1e28bc413d

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

index c76863d55505d742a4ed061223f0940295a8aabb..9b4b7eb761639ba7e6bce6ddeef993032ce1d947 100644 (file)
@@ -61,13 +61,33 @@ namespace PETScWrappers
  *   [1]PETSC ERROR: KSPSetUp() line 195 in src/ksp/ksp/interface/itfunc.c
  * @endverbatim
  *
- * This error, on which one can spend a very long time figuring out what
- * exactly goes wrong, results from not specifying an MPI communicator. Note
- * that the communicator @em must match that of the matrix and all vectors in
- * the linear system which we want to solve. Aggravating the situation is the
- * fact that the default argument to the solver classes, @p PETSC_COMM_SELF,
- * is the appropriate argument for the sequential case (which is why it is the
- * default argument), so this error only shows up in parallel mode.
+ * This error, on which one can spend a very long time figuring out
+ * what exactly goes wrong, results from not specifying an MPI
+ * communicator. Note that the communicator @em must match that of the
+ * matrix and all vectors in the linear system which we want to
+ * solve. Aggravating the situation is the fact that the default
+ * argument to the solver classes, @p PETSC_COMM_SELF, is the
+ * appropriate argument for the sequential case (which is why it is
+ * the default argument), so this error only shows up in parallel
+ * mode.
+ *
+ * Optionally, the user can create a solver derived from the
+ * SolverBase class and can set the default arguments necessary to
+ * solve the linear system of equations with SolverControl. These
+ * default options can be overridden by specifying command line
+ * arguments of the form @verbatim -ksp_*@verbatim. For example,
+ * @verbatim -ksp_monitor_true_residual@verbatim prints out true
+ * residual norm (unpreconditioned) at each iteration and @verbatim
+ * -ksp_view@verbatim provides information about the linear solver and
+ * the preconditioner used in the current context. The type of the
+ * solver can also be changed during runtime by specifying -ksp_type
+ * {richardson, cg, gmres, fgmres, ..} to dynamically test the optimal
+ * solver along with a suitable preconditioner set using -pc_type
+ * {jacobi, bjacobi, ilu, lu, ..}. There are several other command
+ * line options available to modify the behavior of the PETSc linear
+ * solver and can be obtained from the <a
+ * href="http://www.mcs.anl.gov/petsc">documentation and manual
+ * pages</a>.
  *
  * @ingroup PETScWrappers
  * @author Wolfgang Bangerth, 2004
@@ -126,6 +146,15 @@ namespace PETScWrappers
       virtual void reset();
                                        
       
+                                      /**
+                                        * Sets a prefix name for the solver 
+                                        * object. Useful when customizing the
+                                        * PETSc KSP object with command-line
+                                        * options.
+                                        */
+      void set_prefix(const std::string &prefix);
+                                       
+      
                                        /**
                                         * Access to object that controls
                                         * convergence.
@@ -169,6 +198,17 @@ namespace PETScWrappers
       virtual void set_solver_type (KSP &ksp) const = 0;
 
     private:
+                                       /**
+                                        * Solver prefix name to qualify options
+                                        * specific to the PETSc KSP object in the
+                                        * current context.
+                                        * Note: A hyphen (-) must NOT be given 
+                                        * at the beginning of the prefix name. 
+                                        * The first character of all runtime 
+                                        * options is AUTOMATICALLY the hyphen.
+                                        */
+      std::string prefix_name;
+
                                        /**
                                         * A function that is used in PETSc as
                                         * a callback to check on
index c35282df56546da40d4c069267bf8dcbe948f534..b4271900804d9450a7bb21f5d05089eab0a396b2 100644 (file)
@@ -115,6 +115,15 @@ namespace PETScWrappers
 
       }
 
+                                     // set the command line option prefix name
+    ierr = KSPSetOptionsPrefix(solver_data->ksp, prefix_name.c_str());
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+    
+                                     // set the command line options provided
+                                     // by the user to override the defaults
+    ierr = KSPSetFromOptions (solver_data->ksp);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+    
                                      // then do the real work: set up solver
                                      // internal data and solve the
                                      // system. 
@@ -136,6 +145,13 @@ namespace PETScWrappers
   }
 
 
+  void
+  SolverBase::set_prefix(const std::string &prefix)
+  {
+    prefix_name = prefix ;
+  }
+
+
   void
   SolverBase::reset()
   {

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.