MFEM  v3.2
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, 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 "../linalg/hypre.hpp"
21 #include "pfespace.hpp"
22 #include "pgridfunc.hpp"
23 #include "bilinearform.hpp"
24 
25 namespace mfem
26 {
27 
30 {
31 protected:
33  mutable ParGridFunction X, Y; // used in TrueAddMult
34 
36 
38 
39  // called when (mat == NULL && fbfi.Size() > 0)
40  void pAllocMat();
41 
42  void AssembleSharedFaces(int skip_zeros = 1);
43 
44 public:
46  : BilinearForm(pf), pfes(pf), p_mat(NULL), p_mat_e(NULL)
47  { keep_nbr_block = false; }
48 
50  : BilinearForm(pf, bf), pfes(pf), p_mat(NULL), p_mat_e(NULL)
51  { keep_nbr_block = false; }
52 
57  void KeepNbrBlock(bool knb = true) { keep_nbr_block = knb; }
58 
60  void Assemble(int skip_zeros = 1);
61 
64 
67 
70 
74  void ParallelEliminateEssentialBC(const Array<int> &bdr_attr_is_ess,
75  HypreParMatrix &A,
76  const HypreParVector &X,
77  HypreParVector &B) const;
78 
88  HypreParMatrix &A) const;
89 
98  HypreParMatrix &A) const
99  { return A.EliminateRowsCols(tdofs_list); }
100 
102  void TrueAddMult(const Vector &x, Vector &y, const double a = 1.0) const;
103 
105  ParFiniteElementSpace *ParFESpace() const { return pfes; }
106 
109  { return static_cond ? static_cond->GetParTraceFESpace() : NULL; }
110 
132  void FormLinearSystem(Array<int> &ess_tdof_list, Vector &x, Vector &b,
133  HypreParMatrix &A, Vector &X, Vector &B,
134  int copy_interior = 0);
135 
139  void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
140 
141  virtual void Update(FiniteElementSpace *nfes = NULL);
142 
143  virtual ~ParBilinearForm() { delete p_mat_e; delete p_mat; }
144 };
145 
148 {
149 protected:
152  mutable ParGridFunction X, Y; // used in TrueAddMult
153 
154 public:
157  : MixedBilinearForm(trial_fes, test_fes)
158  {
161  }
162 
165 
167  void TrueAddMult(const Vector &x, Vector &y, const double a = 1.0) const;
168 
169  virtual ~ParMixedBilinearForm() { }
170 };
171 
175 {
176 protected:
179 
180 public:
182  ParFiniteElementSpace *rfes)
183  : DiscreteLinearOperator(dfes, rfes) { domain_fes=dfes; range_fes=rfes; }
184 
187 
190  void GetParBlocks(Array2D<HypreParMatrix *> &blocks) const;
191 
193 };
194 
195 }
196 
197 #endif // MFEM_USE_MPI
198 
199 #endif
HypreParMatrix * p_mat_e
ParDiscreteLinearOperator(ParFiniteElementSpace *dfes, ParFiniteElementSpace *rfes)
void EliminateRowsCols(const Array< int > &rows_cols, const HypreParVector &X, HypreParVector &B)
Definition: hypre.cpp:1260
ParFiniteElementSpace * SCParFESpace() const
Return the parallel trace FE space associated with static condensation.
void KeepNbrBlock(bool knb=true)
ParFiniteElementSpace * GetParTraceFESpace()
Return a pointer to the parallel reduced/trace FE space.
Definition: staticcond.hpp:111
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.
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)
HypreParMatrix * p_mat
HypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
ParFiniteElementSpace * pfes
Data type sparse matrix.
Definition: sparsemat.hpp:38
StaticCondensation * static_cond
ParFiniteElementSpace * ParFESpace() const
Return the parallel FE space associated with the ParBilinearForm.
SparseMatrix * mat
Sparse matrix to be associated with the form.
HypreParMatrix * ParallelEliminateTDofs(const Array< int > &tdofs_list, HypreParMatrix &A) const
void Assemble(int skip_zeros=1)
Assemble the local matrix.
ParFiniteElementSpace * test_pfes
HypreParMatrix * ParallelAssemble() const
Returns the matrix &quot;assembled&quot; on the true dofs.
SparseMatrix * mat_e
Matrix used to eliminate b.c.
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:58
ParFiniteElementSpace * domain_fes
HypreParMatrix * ParallelAssembleElim()
Returns the eliminated matrix assembled on the true dofs, i.e. P^t A_e P.
void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
ParBilinearForm(ParFiniteElementSpace *pf)
FiniteElementSpace * test_fes
Class for parallel bilinear form.
ParBilinearForm(ParFiniteElementSpace *pf, ParBilinearForm *bf)
ParFiniteElementSpace * range_fes
void FormLinearSystem(Array< int > &ess_tdof_list, Vector &x, Vector &b, HypreParMatrix &A, Vector &X, Vector &B, int copy_interior=0)
void ParallelEliminateEssentialBC(const Array< int > &bdr_attr_is_ess, HypreParMatrix &A, const HypreParVector &X, HypreParVector &B) const
Class for parallel bilinear form.
ParFiniteElementSpace * trial_pfes
FiniteElementSpace * trial_fes
Vector data type.
Definition: vector.hpp:33
void GetParBlocks(Array2D< HypreParMatrix * > &blocks) const
ParMixedBilinearForm(ParFiniteElementSpace *trial_fes, ParFiniteElementSpace *test_fes)
Class for parallel grid function.
Definition: pgridfunc.hpp:31
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:150
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.