]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added interface to BDDC preconditioner
authorNicolas Barnafi <nabw91@gmail.com>
Wed, 6 Apr 2022 13:48:15 +0000 (15:48 +0200)
committerNicolas Barnafi <nabw91@gmail.com>
Wed, 6 Apr 2022 13:48:15 +0000 (15:48 +0200)
include/deal.II/lac/petsc_precondition.h
source/lac/petsc_precondition.cc

index b17c5e2a9f1d4421ff45790f9d6239c6d10c9b1f..678e4eab2fb328519d90afada1dfb402112a038f 100644 (file)
@@ -957,6 +957,110 @@ namespace PETScWrappers
     AdditionalData additional_data;
   };
 
+    /**
+     * A class that implements the interface to use the BDDC preconditioner from
+     * PETSc.
+     *
+     * @ingroup PETScWrappers
+     * @author Nicolas Barnafi, 2022
+     */
+    class PreconditionBDDC : public PreconditionBase
+    {
+    public:
+      /**
+       * Standardized data struct to pipe additional flags to the
+       * preconditioner.
+       */
+      struct AdditionalData
+      {
+        /**
+         * Constructor. Note that BDDC offers a lot more options to set
+         * than what is exposed here.
+         */
+        AdditionalData(const bool use_vertices_ = true,
+                       const bool use_edges_    = false,
+                       const bool use_faces_    = false,
+                       const bool symmetric_    = false,
+                       const unsigned int coords_cdim_ = 0,
+                       const types::global_dof_index coords_n_ = 0,
+                       const PetscReal *coords_data_ = nullptr);
+
+        /**
+         * This flag sets the use of degrees of freedom in the vertices of the
+         * subdomain as primal variables.
+         */
+        bool use_vertices;
+
+        /**
+         * This flag sets the use of degrees of freedom in the edges of the
+         * subdomain as primal variables.
+         */
+        bool use_edges;
+
+        /**
+         * This flag sets the use of degrees of freedom in the faces of the
+         * subdomain as primal variables.
+         */
+        bool use_faces;
+
+        /**
+         * Set whether the matrix is symmetric or not.
+         */
+        bool symmetric;
+
+        /**
+         * Support for H1 problems
+         */
+        unsigned int coords_cdim;
+        types::global_dof_index coords_n;
+        const PetscReal *coords_data; /* not referenced, consumed at initialize time */
+      };
+
+      /**
+       * Empty Constructor. You need to call initialize() before using this
+       * object.
+       */
+      PreconditionBDDC() = default;
+
+      /**
+       * Constructor. Take the matrix which is used to form the preconditioner,
+       * and additional flags if there are any.
+       */
+      PreconditionBDDC(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.
+       */
+      PreconditionBDDC(const MPI_Comm        communicator,
+                       const AdditionalData &additional_data = AdditionalData());
+
+      /**
+       * Initialize the preconditioner object and calculate all data that is
+       * necessary for applying it in a solver. This function is automatically
+       * called when calling the constructor with the same arguments and is only
+       * used if you create the preconditioner without arguments.
+       */
+      void
+      initialize(const MatrixBase &    matrix,
+                 const AdditionalData &additional_data = AdditionalData());
+
+    protected:
+      /**
+       * Store a copy of the flags for this particular preconditioner.
+       */
+      AdditionalData additional_data;
+
+      /**
+       * Initialize 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();
+    };
+
   /**
    * Alias for backwards-compatibility.
    * @deprecated Use PETScWrappers::PreconditionBase instead.
index cdcc9caac63ef507717c233ec04365217b708ab1..8cd0c0ca2c4011625fb2a91e6b68b3c292ad212b 100644 (file)
@@ -850,6 +850,104 @@ namespace PETScWrappers
     AssertThrow(ierr == 0, ExcPETScError(ierr));
   }
 
+  /* ----------------- PreconditionBDDC -------------------- */
+
+    PreconditionBDDC::AdditionalData::AdditionalData(const bool use_vertices_,
+                                                     const bool use_edges_,
+                                                     const bool use_faces_,
+                                                     const bool symmetric_,
+                                                     const unsigned int coords_cdim_,
+                                                     const types::global_dof_index coords_n_,
+                                                     const PetscReal *coords_data_)
+      : use_vertices(use_vertices_)
+      , use_edges(use_edges_)
+      , use_faces(use_faces_)
+      , symmetric(symmetric_)
+      , coords_cdim(coords_cdim_)
+      , coords_n(coords_n_)
+      , coords_data(coords_data_)
+    {}
+
+    PreconditionBDDC::PreconditionBDDC(const MPI_Comm        comm,
+                                       const AdditionalData &additional_data_)
+    {
+      additional_data = additional_data_;
+
+      PetscErrorCode ierr = PCCreate(comm, &pc);
+      AssertThrow(ierr == 0, ExcPETScError(ierr));
+
+      initialize();
+    }
+
+    PreconditionBDDC::PreconditionBDDC(const MatrixBase &    matrix,
+                                       const AdditionalData &additional_data)
+    {
+      initialize(matrix, additional_data);
+    }
+
+    void
+    PreconditionBDDC::initialize()
+    {
+      PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCBDDC));
+      AssertThrow(ierr == 0, ExcPETScError(ierr));
+
+      // TODO: Add BDDC options, in particular the additional primal nodes
+      std::stringstream ssStream;
+
+      if (additional_data.use_vertices)
+        set_option_value("-pc_bddc_use_vertices", "true");
+      else
+        set_option_value("-pc_bddc_use_vertices", "false");
+      if (additional_data.use_edges)
+        set_option_value("-pc_bddc_use_edges", "true");
+      else
+        set_option_value("-pc_bddc_use_edges", "false");
+      if (additional_data.use_faces)
+        set_option_value("-pc_bddc_use_faces", "true");
+      else
+        set_option_value("-pc_bddc_use_faces", "false");
+      if (additional_data.symmetric)
+        set_option_value("-pc_bddc_symmetric", "true");
+      else
+        set_option_value("-pc_bddc_symmetric", "false");
+      if (additional_data.coords_cdim)
+      {
+        set_option_value("-pc_bddc_corner_selection", "true");
+        ierr = PCSetCoordinates(pc,
+                                additional_data.coords_cdim,
+                                additional_data.coords_n,
+                                (PetscReal*)additional_data.coords_data);
+        AssertThrow(ierr == 0, ExcPETScError(ierr));
+      }
+      else
+      {
+        set_option_value("-pc_bddc_corner_selection", "false");
+        ierr = PCSetCoordinates(pc,0,0,NULL);
+        AssertThrow(ierr == 0, ExcPETScError(ierr));
+      }
+
+
+      ierr = PCSetFromOptions(pc);
+      AssertThrow(ierr == 0, ExcPETScError(ierr));
+    }
+
+    void
+    PreconditionBDDC::initialize(const MatrixBase &    matrix_,
+                                 const AdditionalData &additional_data_)
+    {
+      clear();
+
+      matrix          = static_cast<Mat>(matrix_);
+      additional_data = additional_data_;
+
+      create_pc();
+      initialize();
+
+      PetscErrorCode ierr = PCSetUp(pc);
+      AssertThrow(ierr == 0, ExcPETScError(ierr));
+    }
+
+
 } // namespace PETScWrappers
 
 DEAL_II_NAMESPACE_CLOSE

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.