Line data    Source code

       1             : // The libMesh Finite Element Library.
       2             : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
       3             : 
       4             : // This library is free software; you can redistribute it and/or
       5             : // modify it under the terms of the GNU Lesser General Public
       6             : // License as published by the Free Software Foundation; either
       7             : // version 2.1 of the License, or (at your option) any later version.
       8             : 
       9             : // This library is distributed in the hope that it will be useful,
      10             : // but WITHOUT ANY WARRANTY; without even the implied warranty of
      11             : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
      12             : // Lesser General Public License for more details.
      13             : 
      14             : // You should have received a copy of the GNU Lesser General Public
      15             : // License along with this library; if not, write to the Free Software
      16             : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
      17             : 
      18             : #ifndef LIBMESH_PETSC_SHELL_MATRIX_H
      19             : #define LIBMESH_PETSC_SHELL_MATRIX_H
      20             : 
      21             : 
      22             : #include "libmesh/libmesh_config.h"
      23             : 
      24             : #ifdef LIBMESH_HAVE_PETSC
      25             : 
      26             : // Local includes
      27             : #include "libmesh/libmesh_common.h"
      28             : #include "libmesh/reference_counted_object.h"
      29             : #include "libmesh/libmesh.h"
      30             : #include "libmesh/shell_matrix.h"
      31             : #include "libmesh/petsc_macro.h"
      32             : #include "libmesh/petsc_solver_exception.h"
      33             : #include "libmesh/petsc_vector.h"
      34             : #include "libmesh/libmesh_common.h"
      35             : #include "libmesh/wrapped_petsc.h"
      36             : 
      37             : // Petsc include files.
      38             : #ifdef I
      39             : # define LIBMESH_SAW_I
      40             : #endif
      41             : #include <petscmat.h>
      42             : #ifndef LIBMESH_SAW_I
      43             : # undef I // Avoid complex.h contamination
      44             : #endif
      45             : 
      46             : namespace libMesh
      47             : {
      48             : 
      49             : /**
      50             :  * Initialize a shell matrix object
      51             :  */
      52             : template <typename Obj>
      53             : void init_shell_mat(Obj & obj,
      54             :                     const numeric_index_type m,
      55             :                     const numeric_index_type n,
      56             :                     const numeric_index_type m_l,
      57             :                     const numeric_index_type n_l,
      58             :                     const numeric_index_type blocksize_in);
      59             : 
      60             : /**
      61             :  * Initialize a shell matrix object using information from the \p DofMap
      62             :  */
      63             : template <typename Obj>
      64             : void init_shell_mat(Obj & obj);
      65             : 
      66             : /**
      67             :  * This class allows to use a PETSc shell matrix.
      68             :  * All overridden virtual functions are documented in
      69             :  * shell_matrix.h.
      70             :  *
      71             :  * \author Fande Kong ([email protected])
      72             :  * \date 2019
      73             :  */
      74             : template <typename T>
      75             : class PetscShellMatrix : public ShellMatrix<T>
      76             : {
      77             : public:
      78             : 
      79             :   PetscShellMatrix (const Parallel::Communicator & comm_in);
      80             : 
      81             :   virtual ~PetscShellMatrix ();
      82             : 
      83             :   virtual numeric_index_type m () const override;
      84             : 
      85             :   virtual numeric_index_type n () const override;
      86             : 
      87             :   virtual numeric_index_type local_m () const;
      88             : 
      89             :   virtual numeric_index_type local_n () const;
      90             : 
      91             :   virtual void vector_mult (NumericVector<T> & dest,
      92             :                             const NumericVector<T> & arg) const override;
      93             : 
      94             :   virtual void vector_mult_add (NumericVector<T> & dest,
      95             :                                 const NumericVector<T> & arg) const override;
      96             : 
      97             :   virtual void get_diagonal (NumericVector<T> & dest) const override;
      98             : 
      99             :   virtual void clear () override;
     100             : 
     101             :   virtual void init () override;
     102             : 
     103             :   /**
     104             :    * \returns \p true if the matrix has been initialized,
     105             :    * \p false otherwise.
     106             :    */
     107             :   virtual bool initialized() const;
     108             : 
     109             :   /**
     110             :    * Returns a pointer to the underlying PETSc Mat object. Must call
     111             :    * init() before this.
     112             :    */
     113             :   Mat mat();
     114             : 
     115             : protected:
     116             : 
     117             :   /**
     118             :    * Petsc Shell Matrix
     119             :    */
     120             :   Mat _mat;
     121             : 
     122             :   bool _is_initialized;
     123             : 
     124             :   friend void init_shell_mat<PetscShellMatrix<T>>(PetscShellMatrix<T> & obj);
     125             :   friend void init_shell_mat<PetscShellMatrix<T>>(PetscShellMatrix<T> & obj,
     126             :                                                   const numeric_index_type m,
     127             :                                                   const numeric_index_type n,
     128             :                                                   const numeric_index_type m_l,
     129             :                                                   const numeric_index_type n_l,
     130             :                                                   const numeric_index_type blocksize_in);
     131             : };
     132             : 
     133             : 
     134             : //-----------------------------------------------------------------------
     135             : // PetscShellMatrix inline members
     136             : template <typename T>
     137             : inline
     138         132 : PetscShellMatrix<T>::PetscShellMatrix (const Parallel::Communicator & comm_in):
     139             :   ShellMatrix<T>(comm_in),
     140         132 :   _is_initialized(false)
     141           0 : {}
     142             : 
     143             : 
     144             : 
     145             : template <typename T>
     146             : inline
     147           0 : numeric_index_type PetscShellMatrix<T>::m () const
     148             : {
     149             :   PetscInt m;
     150             : 
     151           0 :   LibmeshPetscCall(MatGetSize(_mat, &m, nullptr));
     152             : 
     153           0 :   return m;
     154             : }
     155             : 
     156             : 
     157             : 
     158             : template <typename T>
     159             : inline
     160           0 : numeric_index_type PetscShellMatrix<T>::n () const
     161             : {
     162             :   PetscInt n;
     163             : 
     164           0 :   LibmeshPetscCall(MatGetSize(_mat, nullptr, &n));
     165             : 
     166           0 :   return n;
     167             : }
     168             : 
     169             : 
     170             : template <typename T>
     171             : inline
     172           0 : numeric_index_type PetscShellMatrix<T>::local_m () const
     173             : {
     174             :   PetscInt m;
     175             : 
     176           0 :   LibmeshPetscCall(MatGetLocalSize(_mat, &m, nullptr));
     177             : 
     178           0 :   return m;
     179             : }
     180             : 
     181             : 
     182             : 
     183             : template <typename T>
     184             : inline
     185           0 : numeric_index_type PetscShellMatrix<T>::local_n () const
     186             : {
     187             :   PetscInt n;
     188             : 
     189           0 :   LibmeshPetscCall(MatGetLocalSize(_mat, nullptr, &n));
     190             : 
     191           0 :   return n;
     192             : }
     193             : 
     194             : 
     195             : template <typename T>
     196             : inline
     197           0 : void PetscShellMatrix<T>::get_diagonal (NumericVector<T> & dest) const
     198             : {
     199             :   // Make sure the NumericVector passed in is really a PetscVector
     200           0 :   PetscVector<T> & petsc_dest = cast_ref<PetscVector<T> &>(dest);
     201             : 
     202           0 :   LibmeshPetscCall(MatGetDiagonal(_mat, petsc_dest.vec()));
     203           0 : }
     204             : 
     205             : 
     206             : 
     207             : } // namespace libMesh
     208             : 
     209             : #endif // LIBMESH_HAVE_PETSC
     210             : #endif // LIBMESH_SPARSE_SHELL_MATRIX_H