MFEM  v4.6.0
Finite element discretization library
handle.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_HANDLE_HPP
13 #define MFEM_HANDLE_HPP
14 
15 #include "../config/config.hpp"
16 #include "operator.hpp"
17 #ifdef MFEM_USE_MPI
18 #include "hypre.hpp"
19 #endif
20 
21 namespace mfem
22 {
23 
24 /// Pointer to an Operator of a specified type
25 /** This class provides a common interface for global, matrix-type operators to
26  be used in bilinear forms, gradients of nonlinear forms, static condensation,
27  hybridization, etc. The following backends are currently supported:
28  - HYPRE parallel sparse matrix (Hypre_ParCSR)
29  - PETSC globally assembled parallel sparse matrix (PETSC_MATAIJ)
30  - PETSC parallel matrix assembled on each processor (PETSC_MATIS)
31  See also Operator::Type.
32 */
34 {
35 protected:
36  static const char not_supported_msg[];
37 
40  bool own_oper;
41 
43 
44  template <typename OpType>
45  void pSet(OpType *A, bool own_A = true)
46  {
47  oper = A;
48  type_id = A->GetType();
49  own_oper = own_A;
50  }
51 
52 public:
53  /** @brief Create an OperatorHandle with type id = Operator::MFEM_SPARSEMAT
54  without allocating the actual matrix. */
56  : oper(NULL), type_id(Operator::MFEM_SPARSEMAT), own_oper(false) { }
57 
58  /** @brief Create a OperatorHandle with a specified type id, @a tid, without
59  allocating the actual matrix. */
61  : oper(NULL), type_id(CheckType(tid)), own_oper(false) { }
62 
63  /// Create an OperatorHandle for the given OpType pointer, @a A.
64  /** Presently, OpType can be SparseMatrix, HypreParMatrix, or PetscParMatrix.
65 
66  The operator ownership flag is set to the value of @a own_A.
67 
68  It is expected that @a A points to a valid object. */
69  template <typename OpType>
70  explicit OperatorHandle(OpType *A, bool own_A = true) { pSet(A, own_A); }
71 
72  /// Shallow copy. The ownership flag of the target is set to false.
74  oper(other.oper), type_id(other.type_id), own_oper(false)
75  { }
76 
77  ~OperatorHandle() { if (own_oper) { delete oper; } }
78 
79  /// Shallow copy. The ownership flag of the target is set to false.
81  {
82  Clear(); oper = master.oper; type_id = master.type_id; own_oper = false;
83  return *this;
84  }
85 
86  /// Access the underlying Operator pointer.
87  Operator *Ptr() const { return oper; }
88 
89  /// Support the use of -> to call methods of the underlying Operator.
90  Operator *operator->() const { return oper; }
91 
92  /// Access the underlying Operator.
93  Operator &operator*() { return *oper; }
94 
95  /// Access the underlying Operator.
96  const Operator &operator*() const { return *oper; }
97 
98  /// Get the currently set operator type id.
99  Operator::Type Type() const { return type_id; }
100 
101  /** @brief Return the Operator pointer statically cast to a specified OpType.
102  Similar to the method Get(). */
103  template <typename OpType>
104  OpType *As() const { return static_cast<OpType*>(oper); }
105 
106  /// Return the Operator pointer dynamically cast to a specified OpType.
107  template <typename OpType>
108  OpType *Is() const { return dynamic_cast<OpType*>(oper); }
109 
110  /// Return the Operator pointer statically cast to a given OpType.
111  /** Similar to the method As(), however the template type OpType can be
112  derived automatically from the argument @a A. */
113  template <typename OpType>
114  void Get(OpType *&A) const { A = static_cast<OpType*>(oper); }
115 
116  /// Return true if the OperatorHandle owns the held Operator.
117  bool OwnsOperator() const { return own_oper; }
118 
119  /// Set the ownership flag for the held Operator.
120  void SetOperatorOwner(bool own = true) { own_oper = own; }
121 
122  /** @brief Clear the OperatorHandle, deleting the held Operator (if owned),
123  while leaving the type id unchanged. */
124  void Clear()
125  {
126  if (own_oper) { delete oper; }
127  oper = NULL;
128  own_oper = false;
129  }
130 
131  /// Invoke Clear() and set a new type id.
133  {
134  Clear();
135  type_id = CheckType(tid);
136  }
137 
138  /// Reset the OperatorHandle to the given OpType pointer, @a A.
139  /** Presently, OpType can be SparseMatrix, HypreParMatrix, or PetscParMatrix.
140 
141  The operator ownership flag is set to the value of @a own_A.
142 
143  It is expected that @a A points to a valid object. */
144  template <typename OpType>
145  void Reset(OpType *A, bool own_A = true)
146  {
147  if (own_oper) { delete oper; }
148  pSet(A, own_A);
149  }
150 
151 #ifdef MFEM_USE_MPI
152  /** @brief Reset the OperatorHandle to hold a parallel square block-diagonal
153  matrix using the currently set type id. */
154  /** The operator ownership flag is set to true. */
155  void MakeSquareBlockDiag(MPI_Comm comm, HYPRE_BigInt glob_size,
156  HYPRE_BigInt *row_starts, SparseMatrix *diag);
157 
158  /** @brief Reset the OperatorHandle to hold a parallel rectangular
159  block-diagonal matrix using the currently set type id. */
160  /** The operator ownership flag is set to true. */
161  void MakeRectangularBlockDiag(MPI_Comm comm, HYPRE_BigInt glob_num_rows,
162  HYPRE_BigInt glob_num_cols, HYPRE_BigInt *row_starts,
163  HYPRE_BigInt *col_starts, SparseMatrix *diag);
164 #endif // MFEM_USE_MPI
165 
166  /// Reset the OperatorHandle to hold the product @a P^t @a A @a P.
167  /** The type id of the result is determined by that of @a A and @a P. The
168  operator ownership flag is set to true. */
170 
171  /** @brief Reset the OperatorHandle to hold the product R @a A @a P, where
172  R = @a Rt^t. */
173  /** The type id of the result is determined by that of @a Rt, @a A, and @a P.
174  The operator ownership flag is set to true. */
176 
177  /// Convert the given OperatorHandle @a A to the currently set type id.
178  /** The operator ownership flag is set to false if the object held by @a A
179  will be held by this object as well, e.g. when the source and destination
180  types are the same; otherwise it is set to true. */
181  void ConvertFrom(OperatorHandle &A);
182 
183  /// Convert the given OpType pointer, @a A, to the currently set type id.
184  /** This method creates a temporary OperatorHandle for @a A and invokes
185  ConvertFrom(OperatorHandle &) with it. */
186  template <typename OpType>
187  void ConvertFrom(OpType *A)
188  {
189  OperatorHandle Ah(A, false);
190  ConvertFrom(Ah);
191  }
192 
193  /** @brief Reset the OperatorHandle to be the eliminated part of @a A after
194  elimination of the essential dofs @a ess_dof_list. */
195  void EliminateRowsCols(OperatorHandle &A, const Array<int> &ess_dof_list);
196 
197  /// Eliminate the rows corresponding to the essential dofs @a ess_dof_list
198  void EliminateRows(const Array<int> &ess_dof_list);
199 
200  /// Eliminate columns corresponding to the essential dofs @a ess_dof_list
201  void EliminateCols(const Array<int> &ess_dof_list);
202 
203  /// Eliminate essential dofs from the solution @a X into the r.h.s. @a B.
204  /** The argument @a A_e is expected to be the result of the method
205  EliminateRowsCols(). */
206  void EliminateBC(const OperatorHandle &A_e, const Array<int> &ess_dof_list,
207  const Vector &X, Vector &B) const;
208 };
209 
210 
211 /// Add an alternative name for OperatorHandle -- OperatorPtr.
213 
214 } // namespace mfem
215 
216 #endif
OperatorHandle()
Create an OperatorHandle with type id = Operator::MFEM_SPARSEMAT without allocating the actual matrix...
Definition: handle.hpp:55
void ConvertFrom(OperatorHandle &A)
Convert the given OperatorHandle A to the currently set type id.
Definition: handle.cpp:200
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
OpType * Is() const
Return the Operator pointer dynamically cast to a specified OpType.
Definition: handle.hpp:108
void Get(OpType *&A) const
Return the Operator pointer statically cast to a given OpType.
Definition: handle.hpp:114
void EliminateRows(const Array< int > &ess_dof_list)
Eliminate the rows corresponding to the essential dofs ess_dof_list.
Definition: handle.cpp:303
const Operator & operator*() const
Access the underlying Operator.
Definition: handle.hpp:96
Operator * oper
Definition: handle.hpp:38
void EliminateCols(const Array< int > &ess_dof_list)
Eliminate columns corresponding to the essential dofs ess_dof_list.
Definition: handle.cpp:321
static const char not_supported_msg[]
Definition: handle.hpp:36
Operator & operator*()
Access the underlying Operator.
Definition: handle.hpp:93
void EliminateBC(const OperatorHandle &A_e, const Array< int > &ess_dof_list, const Vector &X, Vector &B) const
Eliminate essential dofs from the solution X into the r.h.s. B.
Definition: handle.cpp:340
Data type sparse matrix.
Definition: sparsemat.hpp:50
void EliminateRowsCols(OperatorHandle &A, const Array< int > &ess_dof_list)
Reset the OperatorHandle to be the eliminated part of A after elimination of the essential dofs ess_d...
Definition: handle.cpp:252
Operator::Type CheckType(Operator::Type tid)
Definition: handle.cpp:32
void MakeRectangularBlockDiag(MPI_Comm comm, HYPRE_BigInt glob_num_rows, HYPRE_BigInt glob_num_cols, HYPRE_BigInt *row_starts, HYPRE_BigInt *col_starts, SparseMatrix *diag)
Reset the OperatorHandle to hold a parallel rectangular block-diagonal matrix using the currently set...
Definition: handle.cpp:91
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:283
void MakeSquareBlockDiag(MPI_Comm comm, HYPRE_BigInt glob_size, HYPRE_BigInt *row_starts, SparseMatrix *diag)
Reset the OperatorHandle to hold a parallel square block-diagonal matrix using the currently set type...
Definition: handle.cpp:60
OperatorHandle OperatorPtr
Add an alternative name for OperatorHandle – OperatorPtr.
Definition: handle.hpp:212
OperatorHandle(Operator::Type tid)
Create a OperatorHandle with a specified type id, tid, without allocating the actual matrix...
Definition: handle.hpp:60
HYPRE_Int HYPRE_BigInt
void ConvertFrom(OpType *A)
Convert the given OpType pointer, A, to the currently set type id.
Definition: handle.hpp:187
void Clear()
Clear the OperatorHandle, deleting the held Operator (if owned), while leaving the type id unchanged...
Definition: handle.hpp:124
OperatorHandle(const OperatorHandle &other)
Shallow copy. The ownership flag of the target is set to false.
Definition: handle.hpp:73
bool OwnsOperator() const
Return true if the OperatorHandle owns the held Operator.
Definition: handle.hpp:117
Operator::Type Type() const
Get the currently set operator type id.
Definition: handle.hpp:99
void SetOperatorOwner(bool own=true)
Set the ownership flag for the held Operator.
Definition: handle.hpp:120
OpType * As() const
Return the Operator pointer statically cast to a specified OpType. Similar to the method Get()...
Definition: handle.hpp:104
void pSet(OpType *A, bool own_A=true)
Definition: handle.hpp:45
Operator::Type type_id
Definition: handle.hpp:39
Operator * Ptr() const
Access the underlying Operator pointer.
Definition: handle.hpp:87
Vector data type.
Definition: vector.hpp:58
OperatorHandle & operator=(const OperatorHandle &master)
Shallow copy. The ownership flag of the target is set to false.
Definition: handle.hpp:80
Operator * operator->() const
Support the use of -> to call methods of the underlying Operator.
Definition: handle.hpp:90
Abstract operator.
Definition: operator.hpp:24
void MakePtAP(OperatorHandle &A, OperatorHandle &P)
Reset the OperatorHandle to hold the product P^t A P.
Definition: handle.cpp:123
void SetType(Operator::Type tid)
Invoke Clear() and set a new type id.
Definition: handle.hpp:132
void MakeRAP(OperatorHandle &Rt, OperatorHandle &A, OperatorHandle &P)
Reset the OperatorHandle to hold the product R A P, where R = Rt^t.
Definition: handle.cpp:161
OperatorHandle(OpType *A, bool own_A=true)
Create an OperatorHandle for the given OpType pointer, A.
Definition: handle.hpp:70
void Reset(OpType *A, bool own_A=true)
Reset the OperatorHandle to the given OpType pointer, A.
Definition: handle.hpp:145