MFEM  v3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
pfem_extras.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_PFEM_EXTRAS
13 #define MFEM_PFEM_EXTRAS
14 
15 #include "../../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include "mfem.hpp"
20 #include <cstddef>
21 
22 namespace mfem
23 {
24 
25 namespace miniapps
26 {
27 
33 {
34 public:
36  const int p, const int space_dim = 3, const int type = 0,
37  int vdim = 1, int order = Ordering::byNODES);
39 private:
40  const FiniteElementCollection *FEC_;
41 };
42 
48 {
49 public:
50  ND_ParFESpace(ParMesh *m, const int p, const int space_dim,
51  int vdim = 1, int order = Ordering::byNODES);
53 private:
54  const FiniteElementCollection *FEC_;
55 };
56 
62 {
63 public:
64  RT_ParFESpace(ParMesh *m, const int p, const int space_dim,
65  int vdim = 1, int order = Ordering::byNODES);
67 private:
68  const FiniteElementCollection *FEC_;
69 };
70 
76 {
77 public:
78  L2_ParFESpace(ParMesh *m, const int p, const int space_dim,
79  int vdim = 1, int order = Ordering::byNODES);
81 private:
82  const FiniteElementCollection *FEC_;
83 };
84 
86 {
87 public:
89 
91  HYPRE_Int Mult(HypreParVector &x, HypreParVector &y,
92  double alpha = 1.0, double beta = 0.0);
94  HYPRE_Int Mult(HYPRE_ParVector x, HYPRE_ParVector y,
95  double alpha = 1.0, double beta = 0.0);
96 
99  double alpha = 1.0, double beta = 0.0);
100 
102  void Mult(double a, const Vector &x, double b, Vector &y) const;
104  void MultTranspose(double a, const Vector &x, double b, Vector &y) const;
105 
107  void Mult(const Vector &x, Vector &y) const;
109  void MultTranspose(const Vector &x, Vector &y) const;
110 
111  void Update();
112 
113  const HypreParMatrix & GetMatrix() const { return *mat_; }
114 
116 
117 protected:
119 
120  void createMatrix() const;
121 
124 };
125 
127 {
128 public:
130  ParFiniteElementSpace *rfes);
131 };
132 
134 {
135 public:
137  ParFiniteElementSpace *rfes);
138 };
139 
141 {
142 public:
144  ParFiniteElementSpace *rfes);
145 };
146 
151 {
152 public:
154  ParFiniteElementSpace & HCurlFESpace,
156  virtual ~IrrotationalProjector();
157 
158  // Given a vector 'x' of Nedelec DoFs for an arbitrary vector field,
159  // compute the Nedelec DoFs of the irrotational portion, 'y', of
160  // this vector field. The resulting vector will satisfy Curl y = 0
161  // to machine precision.
162  virtual void Mult(const Vector &x, Vector &y) const;
163 
164  void Update();
165 
166 private:
167  ParFiniteElementSpace * H1FESpace_;
168  ParFiniteElementSpace * HCurlFESpace_;
169 
170  ParBilinearForm * s0_;
171  ParBilinearForm * m1_;
172 
173  HypreBoomerAMG * amg_;
174  HyprePCG * pcg_;
175  HypreParMatrix * S0_;
176  HypreParMatrix * M1_;
178  HypreParVector * gradYPot_;
179  HypreParVector * yPot_;
180  HypreParVector * xDiv_;
181 
182  // Array<int> dof_list_;
183  Array<int> ess_bdr_;
184 };
185 
190 {
191 public:
193  ParFiniteElementSpace & HCurlFESpace,
195  virtual ~DivergenceFreeProjector();
196 
197  // Given a vector 'x' of Nedelec DoFs for an arbitrary vector field,
198  // compute the Nedelec DoFs of the divergence free portion, 'y', of
199  // this vector field. The resulting vector will satisfy Div y = 0
200  // in a weak sense.
201  virtual void Mult(const Vector &x, Vector &y) const;
202 
203  void Update();
204 
205 private:
206  ParFiniteElementSpace * HCurlFESpace_;
207  HypreParVector * xIrr_;
208 };
209 
210 
214 void VisualizeField(socketstream &sock, const char *vishost, int visport,
215  ParGridFunction &gf, const char *title,
216  int x = 0, int y = 0, int w = 400, int h = 400);
217 
218 } // namespace miniapps
219 
220 } // namespace mfem
221 
222 #endif // MFEM_USE_MPI
223 #endif
RT_ParFESpace(ParMesh *m, const int p, const int space_dim, int vdim=1, int order=Ordering::byNODES)
Definition: pfem_extras.cpp:51
L2_ParFESpace(ParMesh *m, const int p, const int space_dim, int vdim=1, int order=Ordering::byNODES)
Definition: pfem_extras.cpp:63
IrrotationalProjector(ParFiniteElementSpace &H1FESpace, ParFiniteElementSpace &HCurlFESpace, ParDiscreteInterpolationOperator &Grad)
DivergenceFreeProjector(ParFiniteElementSpace &H1FESpace, ParFiniteElementSpace &HCurlFESpace, ParDiscreteInterpolationOperator &Grad)
HYPRE_Int MultTranspose(HypreParVector &x, HypreParVector &y, double alpha=1.0, double beta=0.0)
Computes y = alpha * A^t * x + beta * y.
Definition: pfem_extras.cpp:98
int vdim
Vector dimension (number of unknowns per degree of freedom).
Definition: fespace.hpp:69
Abstract parallel finite element space.
Definition: pfespace.hpp:28
ND_ParFESpace(ParMesh *m, const int p, const int space_dim, int vdim=1, int order=Ordering::byNODES)
Definition: pfem_extras.cpp:39
The BoomerAMG solver in hypre.
Definition: hypre.hpp:705
HYPRE_Int Mult(HypreParVector &x, HypreParVector &y, double alpha=1.0, double beta=0.0)
Computes y = alpha * A * x + beta * y.
Definition: pfem_extras.cpp:82
ParDiscreteCurlOperator(ParFiniteElementSpace *dfes, ParFiniteElementSpace *rfes)
void VisualizeField(socketstream &sock, const char *vishost, int visport, ParGridFunction &gf, const char *title, int x, int y, int w, int h)
ParDiscreteGradOperator(ParFiniteElementSpace *dfes, ParFiniteElementSpace *rfes)
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:58
PCG solver in hypre.
Definition: hypre.hpp:565
ParDiscreteDivOperator(ParFiniteElementSpace *dfes, ParFiniteElementSpace *rfes)
const HypreParMatrix & GetMatrix() const
Class for parallel bilinear form.
const double alpha
Definition: ex15.cpp:337
Vector data type.
Definition: vector.hpp:33
Class for parallel grid function.
Definition: pgridfunc.hpp:31
Abstract operator.
Definition: operator.hpp:21
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:150
virtual void Mult(const Vector &x, Vector &y) const
Operator application.
Class for parallel meshes.
Definition: pmesh.hpp:28
virtual void Mult(const Vector &x, Vector &y) const
Operator application.
H1_ParFESpace(ParMesh *m, const int p, const int space_dim=3, const int type=0, int vdim=1, int order=Ordering::byNODES)
Definition: pfem_extras.cpp:26