]> https://gitweb.dealii.org/ - dealii.git/commitdiff
allow creation of PETSc preconditioners without a matrix (Jacobi/BoomerAMG)
authorDenis Davydov <davydden@gmail.com>
Tue, 15 Dec 2015 10:32:18 +0000 (11:32 +0100)
committerDenis Davydov <davydden@gmail.com>
Tue, 29 Dec 2015 09:12:03 +0000 (10:12 +0100)
include/deal.II/lac/petsc_precondition.h
source/lac/petsc_precondition.cc

index 4750f06b0c6a7e09f17fa6b45c97bcf7dbd282f2..48008b4fcdbed2e005bc1232eb4d7f57b3f530f5 100644 (file)
@@ -137,12 +137,6 @@ namespace PETScWrappers
      */
     PreconditionJacobi ();
 
-    /**
-     * A constructor intended to be used with SLEPc objects.
-     */
-    PreconditionJacobi (const MPI_Comm communicator,
-                        const AdditionalData &additional_data = AdditionalData());
-
 
     /**
      * Constructor. Take the matrix which is used to form the preconditioner,
@@ -151,6 +145,14 @@ namespace PETScWrappers
     PreconditionJacobi (const MatrixBase     &matrix,
                         const AdditionalData &additional_data = AdditionalData());
 
+    /**
+     * Same as above but without setting a matrix to form the preconditioner.
+     * Intended to be used with SLEPc objects.
+     */
+    PreconditionJacobi (const MPI_Comm communicator,
+                        const AdditionalData &additional_data = AdditionalData());
+
+
     /**
      * Initializes the preconditioner object and calculate all data that is
      * necessary for applying it in a solver. This function is automatically
@@ -217,6 +219,14 @@ namespace PETScWrappers
     PreconditionBlockJacobi (const MatrixBase     &matrix,
                              const AdditionalData &additional_data = AdditionalData());
 
+    /**
+     * Same as above but without setting a matrix to form the preconditioner.
+     * Intended to be used with SLEPc objects.
+     */
+    PreconditionBlockJacobi (const MPI_Comm communicator,
+                             const AdditionalData &additional_data = AdditionalData());
+
+
     /**
      * Initializes the preconditioner object and calculate all data that is
      * necessary for applying it in a solver. This function is automatically
@@ -231,6 +241,14 @@ namespace PETScWrappers
      * Store a copy of the flags for this particular preconditioner.
      */
     AdditionalData additional_data;
+
+    /**
+     * Initializes the preconditioner object without knowing a particular matrix.
+     * This function sets up appropriate parameters to the underlying PETSc object
+     * after it has been created.
+     */
+    void initialize();
+
   };
 
 
@@ -713,6 +731,14 @@ namespace PETScWrappers
     PreconditionBoomerAMG (const MatrixBase     &matrix,
                            const AdditionalData &additional_data = AdditionalData());
 
+    /**
+     * Same as above but without setting a matrix to form the preconditioner.
+     * Intended to be used with SLEPc objects.
+     */
+    PreconditionBoomerAMG (const MPI_Comm communicator,
+                           const AdditionalData &additional_data = AdditionalData());
+
+
     /**
      * Initializes the preconditioner object and calculate all data that is
      * necessary for applying it in a solver. This function is automatically
@@ -727,6 +753,14 @@ namespace PETScWrappers
      * Store a copy of the flags for this particular preconditioner.
      */
     AdditionalData additional_data;
+
+    /**
+     * Initializes the preconditioner object without knowing a particular matrix.
+     * This function sets up appropriate parameters to the underlying PETSc object
+     * after it has been created.
+     */
+    void initialize();
+
   };
 
 
index b4ab5d23f41133f76426520a852d76a4db0419ec..8906e0c466468fb33de89ae9a67572bbef140f53 100644 (file)
@@ -150,6 +150,17 @@ namespace PETScWrappers
 
 
   /* ----------------- PreconditionBlockJacobi -------------------- */
+  PreconditionBlockJacobi::PreconditionBlockJacobi (const MPI_Comm comm,
+                                                    const AdditionalData &additional_data_)
+  {
+    additional_data = additional_data_;
+
+    int ierr = PCCreate(comm, &pc);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+    initialize();
+  }
+
 
   PreconditionBlockJacobi::PreconditionBlockJacobi ()
   {}
@@ -163,22 +174,28 @@ namespace PETScWrappers
   }
 
   void
-  PreconditionBlockJacobi::initialize (const MatrixBase     &matrix_,
-                                       const AdditionalData &additional_data_)
+  PreconditionBlockJacobi::initialize()
   {
-    matrix = static_cast<Mat>(matrix_);
-    additional_data = additional_data_;
-
-    create_pc();
-
     int ierr;
     ierr = PCSetType (pc, const_cast<char *>(PCBJACOBI));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
     ierr = PCSetFromOptions (pc);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+  }
 
-    ierr = PCSetUp (pc);
+
+  void
+  PreconditionBlockJacobi::initialize (const MatrixBase     &matrix_,
+                                       const AdditionalData &additional_data_)
+  {
+    matrix = static_cast<Mat>(matrix_);
+    additional_data = additional_data_;
+
+    create_pc();
+    initialize();
+
+    int ierr = PCSetUp (pc);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
   }
 
@@ -434,6 +451,24 @@ namespace PETScWrappers
   PreconditionBoomerAMG::PreconditionBoomerAMG ()
   {}
 
+  PreconditionBoomerAMG::PreconditionBoomerAMG (const MPI_Comm comm,
+                                                const AdditionalData &additional_data_)
+  {
+    additional_data = additional_data_;
+
+    int ierr = PCCreate(comm, &pc);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+#ifdef PETSC_HAVE_HYPRE
+    initialize();
+#else // PETSC_HAVE_HYPRE
+    (void)pc;
+    Assert (false,
+            ExcMessage ("Your PETSc installation does not include a copy of "
+                        "the hypre package necessary for this preconditioner."));
+#endif
+  }
+
 
   PreconditionBoomerAMG::PreconditionBoomerAMG (const MatrixBase     &matrix,
                                                 const AdditionalData &additional_data)
@@ -441,17 +476,9 @@ namespace PETScWrappers
     initialize(matrix, additional_data);
   }
 
-
   void
-  PreconditionBoomerAMG::initialize (const MatrixBase     &matrix_,
-                                     const AdditionalData &additional_data_)
+  PreconditionBoomerAMG::initialize ()
   {
-    matrix = static_cast<Mat>(matrix_);
-    additional_data = additional_data_;
-
-#ifdef PETSC_HAVE_HYPRE
-    create_pc();
-
     int ierr;
     ierr = PCSetType (pc, const_cast<char *>(PCHYPRE));
     AssertThrow (ierr == 0, ExcPETScError(ierr));
@@ -490,8 +517,20 @@ namespace PETScWrappers
 
     ierr = PCSetFromOptions (pc);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
+  }
 
-    ierr = PCSetUp (pc);
+  void
+  PreconditionBoomerAMG::initialize (const MatrixBase     &matrix_,
+                                     const AdditionalData &additional_data_)
+  {
+    matrix = static_cast<Mat>(matrix_);
+    additional_data = additional_data_;
+
+#ifdef PETSC_HAVE_HYPRE
+    create_pc();
+    initialize ();
+
+    int ierr = PCSetUp (pc);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
 #else // PETSC_HAVE_HYPRE

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.