MFEM  v3.4
Finite element discretization library
pbilinearform.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
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 {
30 protected:
32  mutable ParGridFunction X, Y; // used in TrueAddMult
33 
35 
37 
38  // Allocate mat - called when (mat == NULL && fbfi.Size() > 0)
39  void pAllocMat();
40 
41  void AssembleSharedFaces(int skip_zeros = 1);
42 
43 public:
45  : BilinearForm(pf), pfes(pf),
47  { keep_nbr_block = false; }
48 
50  : BilinearForm(pf, bf), pfes(pf),
52  { keep_nbr_block = false; }
53 
54  /** When set to true and the ParBilinearForm has interior face integrators,
55  the local SparseMatrix will include the rows (in addition to the columns)
56  corresponding to face-neighbor dofs. The default behavior is to disregard
57  those rows. Must be called before the first Assemble call. */
58  void KeepNbrBlock(bool knb = true) { keep_nbr_block = knb; }
59 
60  /// Set the operator type id for the parallel matrix/operator.
61  /** If using static condensation or hybridization, call this method *after*
62  enabling it. */
64  {
65  p_mat.SetType(tid); p_mat_e.SetType(tid);
68  }
69 
70  /// Assemble the local matrix
71  void Assemble(int skip_zeros = 1);
72 
73  /// Returns the matrix assembled on the true dofs, i.e. P^t A P.
74  /** The returned matrix has to be deleted by the caller. */
76 
77  /// Returns the eliminated matrix assembled on the true dofs, i.e. P^t A_e P.
78  /** The returned matrix has to be deleted by the caller. */
80 
81  /// Return the matrix @a m assembled on the true dofs, i.e. P^t A P.
82  /** The returned matrix has to be deleted by the caller. */
84 
85  /** @brief Returns the matrix assembled on the true dofs, i.e.
86  @a A = P^t A_local P, in the format (type id) specified by @a A. */
88 
89  /** Returns the eliminated matrix assembled on the true dofs, i.e.
90  @a A_elim = P^t A_elim_local P in the format (type id) specified by @a A.
91  */
93  { ParallelAssemble(A_elim, mat_e); }
94 
95  /** Returns the matrix @a A_local assembled on the true dofs, i.e.
96  @a A = P^t A_local P in the format (type id) specified by @a A. */
97  void ParallelAssemble(OperatorHandle &A, SparseMatrix *A_local);
98 
99  /// Eliminate essential boundary DOFs from a parallel assembled system.
100  /** The array @a bdr_attr_is_ess marks boundary attributes that constitute
101  the essential part of the boundary. */
102  void ParallelEliminateEssentialBC(const Array<int> &bdr_attr_is_ess,
103  HypreParMatrix &A,
104  const HypreParVector &X,
105  HypreParVector &B) const;
106 
107  /// Eliminate essential boundary DOFs from a parallel assembled matrix @a A.
108  /** The array @a bdr_attr_is_ess marks boundary attributes that constitute
109  the essential part of the boundary. The eliminated part is stored in a
110  matrix A_elim such that A_original = A_new + A_elim. Returns a pointer to
111  the newly allocated matrix A_elim which should be deleted by the caller.
112  The matrices @a A and A_elim can be used to eliminate boundary conditions
113  in multiple right-hand sides, by calling the function EliminateBC() (from
114  hypre.hpp). */
116  HypreParMatrix &A) const;
117 
118  /// Eliminate essential true DOFs from a parallel assembled matrix @a A.
119  /** Given a list of essential true dofs and the parallel assembled matrix
120  @a A, eliminate the true dofs from the matrix, storing the eliminated
121  part in a matrix A_elim such that A_original = A_new + A_elim. Returns a
122  pointer to the newly allocated matrix A_elim which should be deleted by
123  the caller. The matrices @a A and A_elim can be used to eliminate
124  boundary conditions in multiple right-hand sides, by calling the function
125  EliminateBC() (from hypre.hpp). */
127  HypreParMatrix &A) const
128  { return A.EliminateRowsCols(tdofs_list); }
129 
130  /** @brief Compute @a y += @a a (P^t A P) @a x, where @a x and @a y are
131  vectors on the true dofs. */
132  void TrueAddMult(const Vector &x, Vector &y, const double a = 1.0) const;
133 
134  /// Return the parallel FE space associated with the ParBilinearForm.
135  ParFiniteElementSpace *ParFESpace() const { return pfes; }
136 
137  /// Return the parallel trace FE space associated with static condensation.
139  { return static_cond ? static_cond->GetParTraceFESpace() : NULL; }
140 
141  /// Get the parallel finite element space prolongation matrix
142  virtual const Operator *GetProlongation() const
143  { return pfes->GetProlongationMatrix(); }
144  /// Get the parallel finite element space restriction matrix
145  virtual const Operator *GetRestriction() const
146  { return pfes->GetRestrictionMatrix(); }
147 
148  /** Form the linear system A X = B, corresponding to the current bilinear
149  form and b(.), by applying any necessary transformations such as:
150  eliminating boundary conditions; applying conforming constraints for
151  non-conforming AMR; parallel assembly; static condensation;
152  hybridization.
153 
154  The ParGridFunction-size vector x must contain the essential b.c. The
155  ParBilinearForm and the ParLinearForm-size vector b must be assembled.
156 
157  The vector X is initialized with a suitable initial guess: when using
158  hybridization, the vector X is set to zero; otherwise, the essential
159  entries of X are set to the corresponding b.c. and all other entries are
160  set to zero (copy_interior == 0) or copied from x (copy_interior != 0).
161 
162  This method can be called multiple times (with the same ess_tdof_list
163  array) to initialize different right-hand sides and boundary condition
164  values.
165 
166  After solving the linear system, the finite element solution x can be
167  recovered by calling RecoverFEMSolution (with the same vectors X, b, and
168  x). */
169  void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x, Vector &b,
170  OperatorHandle &A, Vector &X, Vector &B,
171  int copy_interior = 0);
172 
173  /** Version of the method FormLinearSystem() where the system matrix is
174  returned in the variable @a A, of type OpType, holding a *reference* to
175  the system matrix (created with the method OpType::MakeRef()). The
176  reference will be invalidated when SetOperatorType(), Update(), or the
177  destructor is called. */
178  template <typename OpType>
179  void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x, Vector &b,
180  OpType &A, Vector &X, Vector &B,
181  int copy_interior = 0)
182  {
183  OperatorHandle Ah;
184  FormLinearSystem(ess_tdof_list, x, b, Ah, X, B, copy_interior);
185  OpType *A_ptr = Ah.Is<OpType>();
186  MFEM_VERIFY(A_ptr, "invalid OpType used");
187  A.MakeRef(*A_ptr);
188  }
189 
190  /// Form the linear system matrix @a A, see FormLinearSystem() for details.
191  void FormSystemMatrix(const Array<int> &ess_tdof_list, OperatorHandle &A);
192 
193  /** Version of the method FormSystemMatrix() where the system matrix is
194  returned in the variable @a A, of type OpType, holding a *reference* to
195  the system matrix (created with the method OpType::MakeRef()). The
196  reference will be invalidated when SetOperatorType(), Update(), or the
197  destructor is called. */
198  template <typename OpType>
199  void FormSystemMatrix(const Array<int> &ess_tdof_list, OpType &A)
200  {
201  OperatorHandle Ah;
202  FormSystemMatrix(ess_tdof_list, Ah);
203  OpType *A_ptr = Ah.Is<OpType>();
204  MFEM_VERIFY(A_ptr, "invalid OpType used");
205  A.MakeRef(*A_ptr);
206  }
207 
208  /** Call this method after solving a linear system constructed using the
209  FormLinearSystem method to recover the solution as a ParGridFunction-size
210  vector in x. Use the same arguments as in the FormLinearSystem call. */
211  virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
212 
213  virtual void Update(FiniteElementSpace *nfes = NULL);
214 
215  virtual ~ParBilinearForm() { }
216 };
217 
218 /// Class for parallel bilinear form using different test and trial FE spaces.
220 {
221 protected:
224  mutable ParGridFunction X, Y; // used in TrueAddMult
225 
226 public:
230  {
233  }
234 
235  /// Returns the matrix assembled on the true dofs, i.e. P_test^t A P_trial.
237 
238  /** @brief Returns the matrix assembled on the true dofs, i.e.
239  @a A = P_test^t A_local P_trial, in the format (type id) specified by
240  @a A. */
242 
243  /// Compute y += a (P^t A P) x, where x and y are vectors on the true dofs
244  void TrueAddMult(const Vector &x, Vector &y, const double a = 1.0) const;
245 
246  virtual ~ParMixedBilinearForm() { }
247 };
248 
249 /** The parallel matrix representation a linear operator between parallel finite
250  element spaces */
252 {
253 protected:
256 
257 public:
259  ParFiniteElementSpace *rfes)
260  : DiscreteLinearOperator(dfes, rfes) { domain_fes=dfes; range_fes=rfes; }
261 
262  /// Returns the matrix "assembled" on the true dofs
264 
265  /** Extract the parallel blocks corresponding to the vector dimensions of the
266  domain and range parallel finite element spaces */
267  void GetParBlocks(Array2D<HypreParMatrix *> &blocks) const;
268 
270 };
271 
272 }
273 
274 #endif // MFEM_USE_MPI
275 
276 #endif
ParDiscreteLinearOperator(ParFiniteElementSpace *dfes, ParFiniteElementSpace *rfes)
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:1368
ParFiniteElementSpace * ParFESpace() const
Return the parallel FE space associated with the ParBilinearForm.
virtual const Operator * GetProlongationMatrix() const
Definition: pfespace.cpp:742
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OpType &A, Vector &X, Vector &B, int copy_interior=0)
void KeepNbrBlock(bool knb=true)
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.
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:94
ParFiniteElementSpace * GetParTraceFESpace()
Return a pointer to the parallel reduced/trace FE space.
Definition: staticcond.hpp:113
HypreParMatrix * ParallelAssemble() const
Returns the matrix "assembled" on the true dofs.
OperatorHandle p_mat_e
void GetParBlocks(Array2D< HypreParMatrix *> &blocks) const
Abstract parallel finite element space.
Definition: pfespace.hpp:28
void AssembleSharedFaces(int skip_zeros=1)
virtual const SparseMatrix * GetRestrictionMatrix() const
Get the R matrix which restricts a local dof vector to true dof vector.
Definition: pfespace.hpp:322
HypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
virtual void Update(FiniteElementSpace *nfes=NULL)
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.
ParFiniteElementSpace * pfes
void ParallelAssembleElim(OperatorHandle &A_elim)
virtual const Operator * GetRestriction() const
Get the parallel finite element space restriction matrix.
Data type sparse matrix.
Definition: sparsemat.hpp:38
StaticCondensation * static_cond
virtual const Operator * GetProlongation() const
Get the parallel finite element space prolongation matrix.
SparseMatrix * mat
Sparse matrix to be associated with the form.
void Assemble(int skip_zeros=1)
Assemble the local matrix.
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.
ParFiniteElementSpace * test_pfes
ParFiniteElementSpace * SCParFESpace() const
Return the parallel trace FE space associated with static condensation.
Dynamic 2D array using row-major layout.
Definition: array.hpp:289
SparseMatrix * mat_e
Matrix used to eliminate b.c.
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:73
ParFiniteElementSpace * domain_fes
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:124
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:66
HypreParMatrix * ParallelEliminateTDofs(const Array< int > &tdofs_list, HypreParMatrix &A) const
Eliminate essential true DOFs from a parallel assembled matrix A.
HypreParMatrix * ParallelAssembleElim()
Returns the eliminated matrix assembled on the true dofs, i.e. P^t A_e P.
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 RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
ParBilinearForm(ParFiniteElementSpace *pf)
FiniteElementSpace * test_fes
void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Form the linear system matrix A, see FormLinearSystem() for details.
Class for parallel bilinear form using different test and trial FE spaces.
ParBilinearForm(ParFiniteElementSpace *pf, ParBilinearForm *bf)
ParFiniteElementSpace * range_fes
void FormSystemMatrix(const Array< int > &ess_tdof_list, OpType &A)
Class for parallel bilinear form.
ParFiniteElementSpace * trial_pfes
FiniteElementSpace * trial_fes
OperatorHandle p_mat
Vector data type.
Definition: vector.hpp:48
ID for class HypreParMatrix.
Definition: operator.hpp:128
Hybridization * hybridization
ParMixedBilinearForm(ParFiniteElementSpace *trial_fes, ParFiniteElementSpace *test_fes)
Class for parallel grid function.
Definition: pgridfunc.hpp:32
Abstract operator.
Definition: operator.hpp:21
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:175
void SetType(Operator::Type tid)
Invoke Clear() and set a new type id.
Definition: handle.hpp:118
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...