MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
pbilinearform.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, 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_PBILINEARFORM
13 #define MFEM_PBILINEARFORM
14 
15 #include "../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include <mpi.h>
20 #include "pfespace.hpp"
21 #include "pgridfunc.hpp"
22 #include "bilinearform.hpp"
23 
24 namespace mfem
25 {
26 
27 /// Class for parallel bilinear form
29 {
31 protected:
32  ParFiniteElementSpace *pfes; ///< Points to the same object as #fes
33 
34  /// Auxiliary objects used in TrueAddMult().
36  mutable Vector Ytmp;
37 
39 
41 
42  // Allocate mat - called when (mat == NULL && fbfi.Size() > 0)
43  void pAllocMat();
44 
45  void AssembleSharedFaces(int skip_zeros = 1);
46 
47 private:
48  /// Copy construction is not supported; body is undefined.
50 
51  /// Copy assignment is not supported; body is undefined.
52  ParBilinearForm &operator=(const ParBilinearForm &);
53 
54 public:
55  /// Creates parallel bilinear form associated with the FE space @a *pf.
56  /** The pointer @a pf is not owned by the newly constructed object. */
58  : BilinearForm(pf), pfes(pf),
60  { keep_nbr_block = false; }
61 
62  /** @brief Create a ParBilinearForm on the ParFiniteElementSpace @a *pf,
63  using the same integrators as the ParBilinearForm @a *bf.
64 
65  The pointer @a pf is not owned by the newly constructed object.
66 
67  The integrators in @a bf are copied as pointers and they are not owned by
68  the newly constructed ParBilinearForm. */
70  : BilinearForm(pf, bf), pfes(pf),
72  { keep_nbr_block = false; }
73 
74  /** When set to true and the ParBilinearForm has interior face integrators,
75  the local SparseMatrix will include the rows (in addition to the columns)
76  corresponding to face-neighbor dofs. The default behavior is to disregard
77  those rows. Must be called before the first Assemble call. */
78  void KeepNbrBlock(bool knb = true) { keep_nbr_block = knb; }
79 
80  /** @brief Set the operator type id for the parallel matrix/operator when
81  using AssemblyLevel::LEGACY. */
82  /** If using static condensation or hybridization, call this method *after*
83  enabling it. */
85  {
86  p_mat.SetType(tid); p_mat_e.SetType(tid);
89  }
90 
91  /// Assemble the local matrix
92  void Assemble(int skip_zeros = 1);
93 
94  /** @brief Assemble the diagonal of the bilinear form into @a diag. Note that
95  @a diag is a true-dof Vector.
96 
97  When the AssemblyLevel is not LEGACY, and the mesh is nonconforming,
98  this method returns |P^T| d_l, where d_l is the local diagonal of the
99  form before applying parallel/conforming assembly, P^T is the transpose
100  of the parallel/conforming prolongation, and |.| denotes the entry-wise
101  absolute value. In general, this is just an approximation of the exact
102  diagonal for this case. */
103  virtual void AssembleDiagonal(Vector &diag) const;
104 
105  /// Returns the matrix assembled on the true dofs, i.e. P^t A P.
106  /** The returned matrix has to be deleted by the caller. */
108 
109  /// Returns the eliminated matrix assembled on the true dofs, i.e. P^t A_e P.
110  /** The returned matrix has to be deleted by the caller. */
112 
113  /// Return the matrix @a m assembled on the true dofs, i.e. P^t A P.
114  /** The returned matrix has to be deleted by the caller. */
116 
117  /** @brief Compute parallel RAP operator and store it in @a A as a HypreParMatrix.
118 
119  @param[in] loc_A The rank-local `SparseMatrix`.
120  @param[out] A The `OperatorHandle` containing the global `HypreParMatrix`.
121  @param[in] steal_loc_A Have the `HypreParMatrix` in @a A take ownership of
122  the memory objects in @a loc_A.
123  */
124  void ParallelRAP(SparseMatrix &loc_A,
125  OperatorHandle &A,
126  bool steal_loc_A = false);
127 
128  /** @brief Returns the matrix assembled on the true dofs, i.e.
129  @a A = P^t A_local P, in the format (type id) specified by @a A. */
131 
132  /** Returns the eliminated matrix assembled on the true dofs, i.e.
133  @a A_elim = P^t A_elim_local P in the format (type id) specified by @a A.
134  */
136  { ParallelAssemble(A_elim, mat_e); }
137 
138  /** Returns the matrix @a A_local assembled on the true dofs, i.e.
139  @a A = P^t A_local P in the format (type id) specified by @a A. */
140  void ParallelAssemble(OperatorHandle &A, SparseMatrix *A_local);
141 
142  /// Eliminate essential boundary DOFs from a parallel assembled system.
143  /** The array @a bdr_attr_is_ess marks boundary attributes that constitute
144  the essential part of the boundary. */
145  void ParallelEliminateEssentialBC(const Array<int> &bdr_attr_is_ess,
146  HypreParMatrix &A,
147  const HypreParVector &X,
148  HypreParVector &B) const;
149 
150  /// Eliminate essential boundary DOFs from a parallel assembled matrix @a A.
151  /** The array @a bdr_attr_is_ess marks boundary attributes that constitute
152  the essential part of the boundary. The eliminated part is stored in a
153  matrix A_elim such that A_original = A_new + A_elim. Returns a pointer to
154  the newly allocated matrix A_elim which should be deleted by the caller.
155  The matrices @a A and A_elim can be used to eliminate boundary conditions
156  in multiple right-hand sides, by calling the function EliminateBC() (from
157  hypre.hpp). */
159  HypreParMatrix &A) const;
160 
161  /// Eliminate essential true DOFs from a parallel assembled matrix @a A.
162  /** Given a list of essential true dofs and the parallel assembled matrix
163  @a A, eliminate the true dofs from the matrix, storing the eliminated
164  part in a matrix A_elim such that A_original = A_new + A_elim. Returns a
165  pointer to the newly allocated matrix A_elim which should be deleted by
166  the caller. The matrices @a A and A_elim can be used to eliminate
167  boundary conditions in multiple right-hand sides, by calling the function
168  EliminateBC() (from hypre.hpp). */
170  HypreParMatrix &A) const
171  { return A.EliminateRowsCols(tdofs_list); }
172 
173  /** @brief Compute @a y += @a a (P^t A P) @a x, where @a x and @a y are
174  vectors on the true dofs. */
175  void TrueAddMult(const Vector &x, Vector &y, const double a = 1.0) const;
176 
177  /// Return the parallel FE space associated with the ParBilinearForm.
178  ParFiniteElementSpace *ParFESpace() const { return pfes; }
179 
180  /// Return the parallel trace FE space associated with static condensation.
182  { return static_cond ? static_cond->GetParTraceFESpace() : NULL; }
183 
184  /// Get the parallel finite element space prolongation matrix
185  virtual const Operator *GetProlongation() const
186  { return pfes->GetProlongationMatrix(); }
187  /// Get the transpose of GetRestriction, useful for matrix-free RAP
188  virtual const Operator *GetRestrictionTranspose() const
190  /// Get the parallel finite element space restriction matrix
191  virtual const Operator *GetRestriction() const
192  { return pfes->GetRestrictionMatrix(); }
193 
196 
197  virtual void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x,
198  Vector &b, OperatorHandle &A, Vector &X,
199  Vector &B, int copy_interior = 0);
200 
201  virtual void FormSystemMatrix(const Array<int> &ess_tdof_list,
202  OperatorHandle &A);
203 
204  /** Call this method after solving a linear system constructed using the
205  FormLinearSystem method to recover the solution as a ParGridFunction-size
206  vector in x. Use the same arguments as in the FormLinearSystem call. */
207  virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
208 
209  virtual void Update(FiniteElementSpace *nfes = NULL);
210 
211  void EliminateVDofsInRHS(const Array<int> &vdofs, const Vector &x, Vector &b);
212 
213  virtual ~ParBilinearForm() { }
214 };
215 
216 /// Class for parallel bilinear form using different test and trial FE spaces.
218 {
219 protected:
220  /// Points to the same object as #trial_fes
222  /// Points to the same object as #test_fes
224  /// Auxiliary objects used in TrueAddMult().
226 
227  /// Matrix and eliminated matrix
229 
230 private:
231  /// Copy construction is not supported; body is undefined.
233 
234  /// Copy assignment is not supported; body is undefined.
235  ParMixedBilinearForm &operator=(const ParMixedBilinearForm &);
236 
237 public:
238  /** @brief Construct a ParMixedBilinearForm on the given FiniteElementSpace%s
239  @a trial_fes and @a test_fes. */
240  /** The pointers @a trial_fes and @a test_fes are not owned by the newly
241  constructed object. */
244  : MixedBilinearForm(trial_fes, test_fes),
246  {
249  }
250 
251  /** @brief Create a ParMixedBilinearForm on the given FiniteElementSpace%s
252  @a trial_fes and @a test_fes, using the same integrators as the
253  ParMixedBilinearForm @a mbf.
254 
255  The pointers @a trial_fes and @a test_fes are not owned by the newly
256  constructed object.
257 
258  The integrators in @a mbf are copied as pointers and they are not owned
259  by the newly constructed ParMixedBilinearForm. */
262  ParMixedBilinearForm * mbf)
263  : MixedBilinearForm(trial_fes, test_fes, mbf),
265  {
268  }
269 
270  /// Returns the matrix assembled on the true dofs, i.e. P_test^t A P_trial.
272 
273  /** @brief Returns the matrix assembled on the true dofs, i.e.
274  @a A = P_test^t A_local P_trial, in the format (type id) specified by
275  @a A. */
277 
280 
281  /** @brief Return in @a A a parallel (on truedofs) version of this operator.
282 
283  This returns the same operator as FormRectangularLinearSystem(), but does
284  without the transformations of the right-hand side. */
285  virtual void FormRectangularSystemMatrix(const Array<int> &trial_tdof_list,
286  const Array<int> &test_tdof_list,
287  OperatorHandle &A);
288 
289  /** @brief Form the parallel linear system A X = B, corresponding to this mixed
290  bilinear form and the linear form @a b(.).
291 
292  Return in @a A a *reference* to the system matrix that is column-constrained.
293  The reference will be invalidated when SetOperatorType(), Update(), or the
294  destructor is called. */
295  virtual void FormRectangularLinearSystem(const Array<int> &trial_tdof_list,
296  const Array<int> &test_tdof_list, Vector &x,
297  Vector &b, OperatorHandle &A, Vector &X,
298  Vector &B);
299 
300  /// Compute y += a (P^t A P) x, where x and y are vectors on the true dofs
301  void TrueAddMult(const Vector &x, Vector &y, const double a = 1.0) const;
302 
303  virtual ~ParMixedBilinearForm() { }
304 };
305 
306 /** The parallel matrix representation a linear operator between parallel finite
307  element spaces */
309 {
310 protected:
311  /// Points to the same object as #trial_fes
313  /// Points to the same object as #test_fes
315 
316 private:
317  /// Copy construction is not supported; body is undefined.
319 
320  /// Copy assignment is not supported; body is undefined.
322 
323 public:
324  /** @brief Construct a ParDiscreteLinearOperator on the given
325  FiniteElementSpace%s @a dfes (domain FE space) and @a rfes (range FE
326  space). */
327  /** The pointers @a dfes and @a rfes are not owned by the newly constructed
328  object. */
330  ParFiniteElementSpace *rfes)
331  : DiscreteLinearOperator(dfes, rfes) { domain_fes=dfes; range_fes=rfes; }
332 
333  /// Returns the matrix "assembled" on the true dofs
335 
336  /** @brief Returns the matrix assembled on the true dofs, i.e.
337  @a A = R_test A_local P_trial, in the format (type id) specified by
338  @a A. */
340 
341  /** Extract the parallel blocks corresponding to the vector dimensions of the
342  domain and range parallel finite element spaces */
343  void GetParBlocks(Array2D<HypreParMatrix *> &blocks) const;
344 
346 
347  /** @brief Return in @a A a parallel (on truedofs) version of this operator. */
349 
351 };
352 
353 }
354 
355 #endif // MFEM_USE_MPI
356 
357 #endif
ParDiscreteLinearOperator(ParFiniteElementSpace *dfes, ParFiniteElementSpace *rfes)
Construct a ParDiscreteLinearOperator on the given FiniteElementSpaces dfes (domain FE space) and rfe...
virtual void FormRectangularSystemMatrix(OperatorHandle &A)
Return in A a parallel (on truedofs) version of this operator.
void SetOperatorType(Operator::Type tid)
Set the operator type id for the parallel reduced matrix/operator.
Definition: staticcond.hpp:174
void EliminateRowsCols(const Array< int > &rows_cols, const HypreParVector &X, HypreParVector &B)
Definition: hypre.cpp:2272
ParFiniteElementSpace * SCParFESpace() const
Return the parallel trace FE space associated with static condensation.
virtual void AssembleDiagonal(Vector &diag) const
Assemble the diagonal of the bilinear form into diag. Note that diag is a true-dof Vector...
virtual const Operator * GetRestrictionTranspose() const
Get the transpose of GetRestriction, useful for matrix-free RAP.
virtual void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
Form the linear system A X = B, corresponding to this bilinear form and the linear form b(...
void KeepNbrBlock(bool knb=true)
virtual void FormRectangularSystemMatrix(const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, OperatorHandle &A)
Return in A that is column-constrained.
virtual const Operator * GetProlongationMatrix() const
The returned Operator is owned by the FiniteElementSpace.
Definition: pfespace.cpp:1153
ParGridFunction Xaux
Auxiliary objects used in TrueAddMult().
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
ParFiniteElementSpace * GetParTraceFESpace()
Return a pointer to the parallel reduced/trace FE space.
Definition: staticcond.hpp:113
void ParallelRAP(SparseMatrix &loc_A, OperatorHandle &A, bool steal_loc_A=false)
Compute parallel RAP operator and store it in A as a HypreParMatrix.
virtual void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
Form the linear system A X = B, corresponding to this bilinear form and the linear form b(...
void TrueAddMult(const Vector &x, Vector &y, const double a=1.0) const
Compute y += a (P^t A P) x, where x and y are vectors on the true dofs.
OperatorHandle p_mat_e
Abstract parallel finite element space.
Definition: pfespace.hpp:28
void AssembleSharedFaces(int skip_zeros=1)
HypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
virtual void Update(FiniteElementSpace *nfes=NULL)
Update the FiniteElementSpace and delete all data associated with the old one.
HypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P_test^t A P_trial.
void SetOperatorType(Operator::Type tid)
Set the operator type id for the parallel hybridized matrix/operator.
void SetOperatorType(Operator::Type tid)
Set the operator type id for the parallel matrix/operator when using AssemblyLevel::LEGACY.
ParFiniteElementSpace * pfes
Points to the same object as fes.
void ParallelAssembleElim(OperatorHandle &A_elim)
Data type sparse matrix.
Definition: sparsemat.hpp:50
StaticCondensation * static_cond
Owned.
ParFiniteElementSpace * ParFESpace() const
Return the parallel FE space associated with the ParBilinearForm.
SparseMatrix * mat
Sparse matrix to be associated with the form. Owned.
virtual void FormRectangularLinearSystem(const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B)
Form the linear system A X = B, corresponding to this mixed bilinear form and the linear form b(...
Data and methods for fully-assembled bilinear forms.
double b
Definition: lissajous.cpp:42
HypreParMatrix * ParallelEliminateTDofs(const Array< int > &tdofs_list, HypreParMatrix &A) const
Eliminate essential true DOFs from a parallel assembled matrix A.
void EliminateVDofsInRHS(const Array< int > &vdofs, const Vector &x, Vector &b)
void Assemble(int skip_zeros=1)
Assemble the local matrix.
virtual void FormRectangularSystemMatrix(const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, OperatorHandle &A)
Return in A a parallel (on truedofs) version of this operator.
virtual void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Form the linear system matrix A, see FormLinearSystem() for details.
ParFiniteElementSpace * test_pfes
Points to the same object as test_fes.
Dynamic 2D array using row-major layout.
Definition: array.hpp:351
HypreParMatrix * ParallelAssemble() const
Returns the matrix &quot;assembled&quot; on the true dofs.
SparseMatrix * mat_e
Sparse Matrix used to store the eliminations from the b.c. Owned. .
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:161
ParFiniteElementSpace * domain_fes
Points to the same object as trial_fes.
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:256
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:96
HypreParMatrix * ParallelAssembleElim()
Returns the eliminated matrix assembled on the true dofs, i.e. P^t A_e P.
ParGridFunction Xaux
Auxiliary objects used in TrueAddMult().
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
ParBilinearForm(ParFiniteElementSpace *pf)
Creates parallel bilinear form associated with the FE space *pf.
virtual const Operator * GetProlongation() const
Get the parallel finite element space prolongation matrix.
FiniteElementSpace * test_fes
Not owned.
double a
Definition: lissajous.cpp:41
virtual void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Form the linear system matrix A, see FormLinearSystem() for details.
A &quot;square matrix&quot; operator for the associated FE space and BLFIntegrators The sum of all the BLFInteg...
Class for parallel bilinear form using different test and trial FE spaces.
ParBilinearForm(ParFiniteElementSpace *pf, ParBilinearForm *bf)
Create a ParBilinearForm on the ParFiniteElementSpace *pf, using the same integrators as the ParBilin...
virtual const SparseMatrix * GetRestrictionMatrix() const
Get the R matrix which restricts a local dof vector to true dof vector.
Definition: pfespace.hpp:389
ParFiniteElementSpace * range_fes
Points to the same object as test_fes.
void ParallelEliminateEssentialBC(const Array< int > &bdr_attr_is_ess, HypreParMatrix &A, const HypreParVector &X, HypreParVector &B) const
Eliminate essential boundary DOFs from a parallel assembled system.
virtual void FormRectangularLinearSystem(const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B)
Form the parallel linear system A X = B, corresponding to this mixed bilinear form and the linear for...
ParMixedBilinearForm(ParFiniteElementSpace *trial_fes, ParFiniteElementSpace *test_fes, ParMixedBilinearForm *mbf)
Create a ParMixedBilinearForm on the given FiniteElementSpaces trial_fes and test_fes, using the same integrators as the ParMixedBilinearForm mbf.
Class for parallel bilinear form.
virtual const Operator * GetRestriction() const
Get the parallel finite element space restriction matrix.
ParFiniteElementSpace * trial_pfes
Points to the same object as trial_fes.
FiniteElementSpace * trial_fes
Not owned.
virtual const Operator * GetRestrictionTransposeOperator() const
Return logical transpose of restriction matrix, but in non-assembled optimized matrix-free form...
Definition: pfespace.cpp:1217
OperatorHandle p_mat
Vector data type.
Definition: vector.hpp:60
void GetParBlocks(Array2D< HypreParMatrix * > &blocks) const
ID for class HypreParMatrix.
Definition: operator.hpp:260
Hybridization * hybridization
Owned.
ParMixedBilinearForm(ParFiniteElementSpace *trial_fes, ParFiniteElementSpace *test_fes)
Construct a ParMixedBilinearForm on the given FiniteElementSpaces trial_fes and test_fes.
Class for parallel grid function.
Definition: pgridfunc.hpp:32
Array< int > vdofs
Abstract operator.
Definition: operator.hpp:24
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:343
void TrueAddMult(const Vector &x, Vector &y, const double a=1.0) const
Compute y += a (P^t A P) x, where x and y are vectors on the true dofs.
OperatorHandle p_mat
Matrix and eliminated matrix.
ParGridFunction Yaux
void SetType(Operator::Type tid)
Invoke Clear() and set a new type id.
Definition: handle.hpp:132
void ParallelAssemble(OperatorHandle &A)
Returns the matrix assembled on the true dofs, i.e. A = P^t A_local P, in the format (type id) specif...