MFEM  v4.1.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
pgridfunc.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, 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_PGRIDFUNC
13 #define MFEM_PGRIDFUNC
14 
15 #include "../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include "../general/globals.hpp"
20 #include "pfespace.hpp"
21 #include "gridfunc.hpp"
22 #include <iostream>
23 #include <limits>
24 
25 namespace mfem
26 {
27 
28 /// Compute a global Lp norm from the local Lp norms computed by each processor
29 double GlobalLpNorm(const double p, double loc_norm, MPI_Comm comm);
30 
31 /// Class for parallel grid function
33 {
34 protected:
35  ParFiniteElementSpace *pfes; ///< Points to the same object as #fes
36 
37  /** @brief Vector used to store data from face-neighbor processors,
38  initialized by ExchangeFaceNbrData(). */
40 
42  Array<int> &attr);
43 
44 public:
45  ParGridFunction() { pfes = NULL; }
46 
47  /// Copy constructor. The internal vector #face_nbr_data is not copied.
49  : GridFunction(orig), pfes(orig.pfes) { }
50 
52 
53  /// Construct a ParGridFunction using previously allocated array @a data.
54  /** The ParGridFunction does not assume ownership of @a data which is assumed
55  to be of size at least `pf->GetVSize()`. Similar to the GridFunction and
56  Vector constructors for externally allocated array, the pointer @a data
57  can be NULL. The data array can be replaced later using the method
58  SetData().
59  */
61  GridFunction(pf, data), pfes(pf) { }
62 
63  /// Construct a ParGridFunction using a GridFunction as external data.
64  /** The parallel space @a *pf and the space used by @a *gf should match. The
65  data from @a *gf is used as the local data of the ParGridFunction on each
66  processor. The ParGridFunction does not assume ownership of the data. */
68 
69  /** @brief Creates grid function on (all) dofs from a given vector on the
70  true dofs, i.e. P tv. */
72 
73  /** @brief Construct a local ParGridFunction from the given *global*
74  GridFunction. If @a partitioning is NULL (default), the data from @a gf
75  is NOT copied. */
76  ParGridFunction(ParMesh *pmesh, const GridFunction *gf,
77  const int *partitioning = NULL);
78 
79  /** @brief Construct a ParGridFunction on a given ParMesh, @a pmesh, reading
80  from an std::istream.
81 
82  In the process, a ParFiniteElementSpace and a FiniteElementCollection are
83  constructed. The new ParGridFunction assumes ownership of both. */
84  ParGridFunction(ParMesh *pmesh, std::istream &input);
85 
86  /// Copy assignment. Only the data of the base class Vector is copied.
87  /** It is assumed that this object and @a rhs use ParFiniteElementSpace%s
88  that have the same size.
89 
90  @note Defining this method overwrites the implicitly defined copy
91  assignment operator. */
93  { return operator=((const Vector &)rhs); }
94 
95  /// Assign constant values to the ParGridFunction data.
96  ParGridFunction &operator=(double value)
97  { GridFunction::operator=(value); return *this; }
98 
99  /// Copy the data from a Vector to the ParGridFunction data.
101  { GridFunction::operator=(v); return *this; }
102 
103  ParFiniteElementSpace *ParFESpace() const { return pfes; }
104 
105  virtual void Update();
106 
107  /// Associate a new FiniteElementSpace with the ParGridFunction.
108  /** The ParGridFunction is resized using the SetSize() method. The new space
109  @a f is expected to be a ParFiniteElementSpace. */
110  virtual void SetSpace(FiniteElementSpace *f);
111 
112  /// Associate a new parallel space with the ParGridFunction.
114 
115  using GridFunction::MakeRef;
116 
117  /** @brief Make the ParGridFunction reference external data on a new
118  FiniteElementSpace. */
119  /** This method changes the FiniteElementSpace associated with the
120  ParGridFunction and sets the pointer @a v as external data in the
121  ParGridFunction. The new space @a f is expected to be a
122  ParFiniteElementSpace. */
123  virtual void MakeRef(FiniteElementSpace *f, double *v);
124 
125  /** @brief Make the ParGridFunction reference external data on a new
126  ParFiniteElementSpace. */
127  /** This method changes the ParFiniteElementSpace associated with the
128  ParGridFunction and sets the pointer @a v as external data in the
129  ParGridFunction. */
130  void MakeRef(ParFiniteElementSpace *f, double *v);
131 
132  /** @brief Make the ParGridFunction reference external data on a new
133  FiniteElementSpace. */
134  /** This method changes the FiniteElementSpace associated with the
135  ParGridFunction and sets the data of the Vector @a v (plus the @a
136  v_offset) as external data in the ParGridFunction. The new space @a f is
137  expected to be a ParFiniteElementSpace.
138  @note This version of the method will also perform bounds checks when
139  the build option MFEM_DEBUG is enabled. */
140  virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset);
141 
142  /** @brief Make the ParGridFunction reference external data on a new
143  ParFiniteElementSpace. */
144  /** This method changes the ParFiniteElementSpace associated with the
145  ParGridFunction and sets the data of the Vector @a v (plus the
146  @a v_offset) as external data in the ParGridFunction.
147  @note This version of the method will also perform bounds checks when
148  the build option MFEM_DEBUG is enabled. */
149  void MakeRef(ParFiniteElementSpace *f, Vector &v, int v_offset);
150 
151  /** Set the grid function on (all) dofs from a given vector on the
152  true dofs, i.e. P tv. */
153  void Distribute(const Vector *tv);
154  void Distribute(const Vector &tv) { Distribute(&tv); }
155  void AddDistribute(double a, const Vector *tv);
156  void AddDistribute(double a, const Vector &tv) { AddDistribute(a, &tv); }
157 
158  /// Set the GridFunction from the given true-dof vector.
159  virtual void SetFromTrueDofs(const Vector &tv) { Distribute(tv); }
160 
161  /// Short semantic for Distribute()
163  { Distribute(&tv); return (*this); }
164 
166 
167  /// Returns the true dofs in a new HypreParVector
168  HypreParVector *GetTrueDofs() const;
169 
170  /// Returns the vector averaged on the true dofs.
171  void ParallelAverage(Vector &tv) const;
172 
173  /// Returns the vector averaged on the true dofs.
174  void ParallelAverage(HypreParVector &tv) const;
175 
176  /// Returns a new vector averaged on the true dofs.
178 
179  /// Returns the vector restricted to the true dofs.
180  void ParallelProject(Vector &tv) const;
181 
182  /// Returns the vector restricted to the true dofs.
183  void ParallelProject(HypreParVector &tv) const;
184 
185  /// Returns a new vector restricted to the true dofs.
187 
188  /// Returns the vector assembled on the true dofs.
189  void ParallelAssemble(Vector &tv) const;
190 
191  /// Returns the vector assembled on the true dofs.
192  void ParallelAssemble(HypreParVector &tv) const;
193 
194  /// Returns a new vector assembled on the true dofs.
196 
197  void ExchangeFaceNbrData();
199  const Vector &FaceNbrData() const { return face_nbr_data; }
200 
201  // Redefine to handle the case when i is a face-neighbor element
202  virtual double GetValue(int i, const IntegrationPoint &ip,
203  int vdim = 1) const;
205  { return GetValue(T.ElementNo, T.GetIntPoint()); }
206 
208  virtual void ProjectCoefficient(Coefficient &coeff);
209 
211  /** @brief Project a discontinuous vector coefficient as a grid function on
212  a continuous finite element space. The values in shared dofs are
213  determined from the element with maximal attribute. */
214  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff);
215 
216  virtual void ProjectDiscCoefficient(Coefficient &coeff, AvgType type);
217 
218  virtual void ProjectDiscCoefficient(VectorCoefficient &vcoeff, AvgType type);
219 
221 
222  // Only the values in the master are guaranteed to be correct!
224  Array<int> &attr)
225  { ProjectBdrCoefficient(NULL, &vcoeff, attr); }
226 
227  // Only the values in the master are guaranteed to be correct!
228  virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array<int> &attr)
229  { ProjectBdrCoefficient(coeff, NULL, attr); }
230 
231  // Only the values in the master are guaranteed to be correct!
232  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff,
233  Array<int> &bdr_attr);
234 
235  virtual double ComputeL1Error(Coefficient *exsol[],
236  const IntegrationRule *irs[] = NULL) const
237  {
239  *exsol, NULL, 1, NULL, irs), pfes->GetComm());
240  }
241 
242  virtual double ComputeL1Error(Coefficient &exsol,
243  const IntegrationRule *irs[] = NULL) const
244  { return ComputeLpError(1.0, exsol, NULL, irs); }
245 
246  virtual double ComputeL1Error(VectorCoefficient &exsol,
247  const IntegrationRule *irs[] = NULL) const
248  { return ComputeLpError(1.0, exsol, NULL, NULL, irs); }
249 
250  virtual double ComputeL2Error(Coefficient *exsol[],
251  const IntegrationRule *irs[] = NULL) const
252  {
253  return GlobalLpNorm(2.0, GridFunction::ComputeL2Error(exsol, irs),
254  pfes->GetComm());
255  }
256 
257  virtual double ComputeL2Error(Coefficient &exsol,
258  const IntegrationRule *irs[] = NULL) const
259  { return ComputeLpError(2.0, exsol, NULL, irs); }
260 
261  virtual double ComputeL2Error(VectorCoefficient &exsol,
262  const IntegrationRule *irs[] = NULL,
263  Array<int> *elems = NULL) const
264  {
265  return GlobalLpNorm(2.0, GridFunction::ComputeL2Error(exsol, irs, elems),
266  pfes->GetComm());
267  }
268 
269  virtual double ComputeMaxError(Coefficient *exsol[],
270  const IntegrationRule *irs[] = NULL) const
271  {
272  return GlobalLpNorm(infinity(),
273  GridFunction::ComputeMaxError(exsol, irs),
274  pfes->GetComm());
275  }
276 
277  virtual double ComputeMaxError(Coefficient &exsol,
278  const IntegrationRule *irs[] = NULL) const
279  {
280  return ComputeLpError(infinity(), exsol, NULL, irs);
281  }
282 
283  virtual double ComputeMaxError(VectorCoefficient &exsol,
284  const IntegrationRule *irs[] = NULL) const
285  {
286  return ComputeLpError(infinity(), exsol, NULL, NULL, irs);
287  }
288 
289  virtual double ComputeLpError(const double p, Coefficient &exsol,
290  Coefficient *weight = NULL,
291  const IntegrationRule *irs[] = NULL) const
292  {
294  p, exsol, weight, irs), pfes->GetComm());
295  }
296 
297  /** When given a vector weight, compute the pointwise (scalar) error as the
298  dot product of the vector error with the vector weight. Otherwise, the
299  scalar error is the l_2 norm of the vector error. */
300  virtual double ComputeLpError(const double p, VectorCoefficient &exsol,
301  Coefficient *weight = NULL,
302  VectorCoefficient *v_weight = NULL,
303  const IntegrationRule *irs[] = NULL) const
304  {
306  p, exsol, weight, v_weight, irs), pfes->GetComm());
307  }
308 
309  virtual void ComputeFlux(BilinearFormIntegrator &blfi,
310  GridFunction &flux,
311  bool wcoef = true, int subdomain = -1);
312 
313  /** Save the local portion of the ParGridFunction. It differs from the
314  serial GridFunction::Save in that it takes into account the signs of
315  the local dofs. */
316  virtual void Save(std::ostream &out) const;
317 
318  /// Merge the local grid functions
319  void SaveAsOne(std::ostream &out = mfem::out);
320 
321  virtual ~ParGridFunction() { }
322 };
323 
324 
325 /** Performs a global L2 projection (through a HypreBoomerAMG solve) of flux
326  from supplied discontinuous space into supplied smooth (continuous, or at
327  least conforming) space, and computes the Lp norms of the differences
328  between them on each element. This is one approach to handling conforming
329  and non-conforming elements in parallel. Returns the total error estimate. */
330 double L2ZZErrorEstimator(BilinearFormIntegrator &flux_integrator,
331  const ParGridFunction &x,
332  ParFiniteElementSpace &smooth_flux_fes,
333  ParFiniteElementSpace &flux_fes,
334  Vector &errors, int norm_p = 2, double solver_tol = 1e-12,
335  int solver_max_it = 200);
336 
337 }
338 
339 #endif // MFEM_USE_MPI
340 
341 #endif
virtual double ComputeLpError(const double p, VectorCoefficient &exsol, Coefficient *weight=NULL, VectorCoefficient *v_weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:300
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:27
Memory< double > data
Definition: vector.hpp:52
void AddDistribute(double a, const Vector &tv)
Definition: pgridfunc.hpp:156
ParGridFunction(const ParGridFunction &orig)
Copy constructor. The internal vector face_nbr_data is not copied.
Definition: pgridfunc.hpp:48
HypreParVector * ParallelAssemble() const
Returns a new vector assembled on the true dofs.
Definition: pgridfunc.cpp:198
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:318
HypreParVector * ParallelProject() const
Returns a new vector restricted to the true dofs.
Definition: pgridfunc.cpp:181
virtual ~ParGridFunction()
Definition: pgridfunc.hpp:321
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
Definition: gridfunc.cpp:1816
ParFiniteElementSpace * pfes
Points to the same object as fes.
Definition: pgridfunc.hpp:35
ParGridFunction & operator=(const HypreParVector &tv)
Short semantic for Distribute()
Definition: pgridfunc.hpp:162
void SaveAsOne(std::ostream &out=mfem::out)
Merge the local grid functions.
Definition: pgridfunc.cpp:501
ParGridFunction(ParFiniteElementSpace *pf)
Definition: pgridfunc.hpp:51
virtual double GetValue(int i, const IntegrationPoint &ip, int vdim=1) const
Definition: pgridfunc.cpp:263
virtual void Save(std::ostream &out) const
Definition: pgridfunc.cpp:485
virtual void ProjectDiscCoefficient(VectorCoefficient &coeff)
Project a discontinuous vector coefficient as a grid function on a continuous finite element space...
Definition: pgridfunc.cpp:322
virtual void SetSpace(FiniteElementSpace *f)
Associate a new FiniteElementSpace with the ParGridFunction.
Definition: pgridfunc.cpp:87
Abstract parallel finite element space.
Definition: pfespace.hpp:28
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: pgridfunc.cpp:301
const IntegrationPoint & GetIntPoint()
Definition: eltrans.hpp:58
void Distribute(const Vector &tv)
Definition: pgridfunc.hpp:154
virtual void ProjectBdrCoefficient(VectorCoefficient &vcoeff, Array< int > &attr)
Project a VectorCoefficient on the GridFunction, modifying only DOFs on the boundary associated with ...
Definition: pgridfunc.hpp:223
double GetValue(ElementTransformation &T)
Definition: pgridfunc.hpp:204
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the ParGridFunction reference external data on a new FiniteElementSpace.
Definition: pgridfunc.cpp:102
MPI_Comm GetComm() const
Definition: pfespace.hpp:239
ParGridFunction & operator=(double value)
Assign constant values to the ParGridFunction data.
Definition: pgridfunc.hpp:96
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
Definition: pgridfunc.hpp:159
const Vector & FaceNbrData() const
Definition: pgridfunc.hpp:199
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: pgridfunc.cpp:81
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:246
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
Definition: gridfunc.cpp:188
void GetTrueDofs(Vector &tv) const
Extract the true-dofs from the GridFunction. If all dofs are true, then tv will be set to point to th...
Definition: gridfunc.cpp:323
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:289
GridFunction & operator=(const GridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: gridfunc.hpp:106
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:2409
ParGridFunction & operator=(const Vector &v)
Copy the data from a Vector to the ParGridFunction data.
Definition: pgridfunc.hpp:100
Vector & FaceNbrData()
Definition: pgridfunc.hpp:198
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:73
virtual double ComputeL1Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:235
double GlobalLpNorm(const double p, double loc_norm, MPI_Comm comm)
Compute a global Lp norm from the local Lp norms computed by each processor.
Definition: pgridfunc.cpp:633
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:24
ParGridFunction & operator=(const ParGridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: pgridfunc.hpp:92
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:257
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Definition: pgridfunc.cpp:409
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:87
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
virtual void ComputeFlux(BilinearFormIntegrator &blfi, GridFunction &flux, bool wcoef=true, int subdomain=-1)
Definition: pgridfunc.cpp:669
virtual double ComputeMaxError(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:269
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:333
double L2ZZErrorEstimator(BilinearFormIntegrator &flux_integrator, const ParGridFunction &x, ParFiniteElementSpace &smooth_flux_fes, ParFiniteElementSpace &flux_fes, Vector &errors, int norm_p, double solver_tol, int solver_max_it)
Definition: pgridfunc.cpp:706
void Distribute(const Vector *tv)
Definition: pgridfunc.cpp:132
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:277
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:283
double a
Definition: lissajous.cpp:41
virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array< int > &attr)
Project a set of Coefficients on the components of the GridFunction, modifying only DOFs on the bound...
Definition: pgridfunc.hpp:228
Vector face_nbr_data
Vector used to store data from face-neighbor processors, initialized by ExchangeFaceNbrData().
Definition: pgridfunc.hpp:39
Class for integration point with weight.
Definition: intrules.hpp:25
virtual double ComputeL2Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:250
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
Definition: pgridfunc.cpp:447
HypreParVector * ParallelAverage() const
Returns a new vector averaged on the true dofs.
Definition: pgridfunc.cpp:164
HypreParVector * GetTrueDofs() const
Returns the true dofs in a new HypreParVector.
Definition: pgridfunc.cpp:143
void AddDistribute(double a, const Vector *tv)
Definition: pgridfunc.cpp:138
double infinity()
Define a shortcut for std::numeric_limits&lt;double&gt;::infinity()
Definition: vector.hpp:42
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: gridfunc.cpp:1685
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:242
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:2303
virtual double ComputeL2Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL, Array< int > *elems=NULL) const
Definition: pgridfunc.hpp:261
Vector data type.
Definition: vector.hpp:48
Class for parallel grid function.
Definition: pgridfunc.hpp:32
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
Class for parallel meshes.
Definition: pmesh.hpp:32
ParGridFunction(ParFiniteElementSpace *pf, double *data)
Construct a ParGridFunction using previously allocated array data.
Definition: pgridfunc.hpp:60
void ProjectBdrCoefficient(Coefficient &coeff, Array< int > &attr)
Project a Coefficient on the GridFunction, modifying only DOFs on the boundary associated with the bo...
Definition: gridfunc.hpp:287
ParFiniteElementSpace * ParFESpace() const
Definition: pgridfunc.hpp:103