MFEM  v4.2.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 
41  /** @brief Vector used as an MPI buffer to send face-neighbor data
42  in ExchangeFaceNbrData() to neighboring processors. */
43  //TODO: Use temporary memory to avoid CUDA malloc allocation cost.
45 
47  Array<int> &attr);
48 
49 public:
50  ParGridFunction() { pfes = NULL; }
51 
52  /// Copy constructor. The internal vector #face_nbr_data is not copied.
54  : GridFunction(orig), pfes(orig.pfes) { }
55 
57 
58  /// Construct a ParGridFunction using previously allocated array @a data.
59  /** The ParGridFunction does not assume ownership of @a data which is assumed
60  to be of size at least `pf->GetVSize()`. Similar to the GridFunction and
61  Vector constructors for externally allocated array, the pointer @a data
62  can be NULL. The data array can be replaced later using the method
63  SetData().
64  */
66  GridFunction(pf, data), pfes(pf) { }
67 
68  /// Construct a ParGridFunction using a GridFunction as external data.
69  /** The parallel space @a *pf and the space used by @a *gf should match. The
70  data from @a *gf is used as the local data of the ParGridFunction on each
71  processor. The ParGridFunction does not assume ownership of the data. */
73 
74  /** @brief Creates grid function on (all) dofs from a given vector on the
75  true dofs, i.e. P tv. */
77 
78  /** @brief Construct a local ParGridFunction from the given *global*
79  GridFunction. If @a partitioning is NULL (default), the data from @a gf
80  is NOT copied. */
81  ParGridFunction(ParMesh *pmesh, const GridFunction *gf,
82  const int *partitioning = NULL);
83 
84  /** @brief Construct a ParGridFunction on a given ParMesh, @a pmesh, reading
85  from an std::istream.
86 
87  In the process, a ParFiniteElementSpace and a FiniteElementCollection are
88  constructed. The new ParGridFunction assumes ownership of both. */
89  ParGridFunction(ParMesh *pmesh, std::istream &input);
90 
91  /// Copy assignment. Only the data of the base class Vector is copied.
92  /** It is assumed that this object and @a rhs use ParFiniteElementSpace%s
93  that have the same size.
94 
95  @note Defining this method overwrites the implicitly defined copy
96  assignment operator. */
98  { return operator=((const Vector &)rhs); }
99 
100  /// Assign constant values to the ParGridFunction data.
101  ParGridFunction &operator=(double value)
102  { GridFunction::operator=(value); return *this; }
103 
104  /// Copy the data from a Vector to the ParGridFunction data.
106  { GridFunction::operator=(v); return *this; }
107 
108  ParFiniteElementSpace *ParFESpace() const { return pfes; }
109 
110  virtual void Update();
111 
112  /// Associate a new FiniteElementSpace with the ParGridFunction.
113  /** The ParGridFunction is resized using the SetSize() method. The new space
114  @a f is expected to be a ParFiniteElementSpace. */
115  virtual void SetSpace(FiniteElementSpace *f);
116 
117  /// Associate a new parallel space with the ParGridFunction.
119 
120  using GridFunction::MakeRef;
121 
122  /** @brief Make the ParGridFunction reference external data on a new
123  FiniteElementSpace. */
124  /** This method changes the FiniteElementSpace associated with the
125  ParGridFunction and sets the pointer @a v as external data in the
126  ParGridFunction. The new space @a f is expected to be a
127  ParFiniteElementSpace. */
128  virtual void MakeRef(FiniteElementSpace *f, double *v);
129 
130  /** @brief Make the ParGridFunction reference external data on a new
131  ParFiniteElementSpace. */
132  /** This method changes the ParFiniteElementSpace associated with the
133  ParGridFunction and sets the pointer @a v as external data in the
134  ParGridFunction. */
135  void MakeRef(ParFiniteElementSpace *f, double *v);
136 
137  /** @brief Make the ParGridFunction reference external data on a new
138  FiniteElementSpace. */
139  /** This method changes the FiniteElementSpace associated with the
140  ParGridFunction and sets the data of the Vector @a v (plus the @a
141  v_offset) as external data in the ParGridFunction. The new space @a f is
142  expected to be a ParFiniteElementSpace.
143  @note This version of the method will also perform bounds checks when
144  the build option MFEM_DEBUG is enabled. */
145  virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset);
146 
147  /** @brief Make the ParGridFunction reference external data on a new
148  ParFiniteElementSpace. */
149  /** This method changes the ParFiniteElementSpace associated with the
150  ParGridFunction and sets the data of the Vector @a v (plus the
151  @a v_offset) as external data in the ParGridFunction.
152  @note This version of the method will also perform bounds checks when
153  the build option MFEM_DEBUG is enabled. */
154  void MakeRef(ParFiniteElementSpace *f, Vector &v, int v_offset);
155 
156  /** Set the grid function on (all) dofs from a given vector on the
157  true dofs, i.e. P tv. */
158  void Distribute(const Vector *tv);
159  void Distribute(const Vector &tv) { Distribute(&tv); }
160  void AddDistribute(double a, const Vector *tv);
161  void AddDistribute(double a, const Vector &tv) { AddDistribute(a, &tv); }
162 
163  /// Set the GridFunction from the given true-dof vector.
164  virtual void SetFromTrueDofs(const Vector &tv) { Distribute(tv); }
165 
166  /// Short semantic for Distribute()
168  { Distribute(&tv); return (*this); }
169 
171 
172  /// Returns the true dofs in a new HypreParVector
173  HypreParVector *GetTrueDofs() const;
174 
175  /// Returns the vector averaged on the true dofs.
176  void ParallelAverage(Vector &tv) const;
177 
178  /// Returns the vector averaged on the true dofs.
179  void ParallelAverage(HypreParVector &tv) const;
180 
181  /// Returns a new vector averaged on the true dofs.
183 
184  /// Returns the vector restricted to the true dofs.
185  void ParallelProject(Vector &tv) const;
186 
187  /// Returns the vector restricted to the true dofs.
188  void ParallelProject(HypreParVector &tv) const;
189 
190  /// Returns a new vector restricted to the true dofs.
192 
193  /// Returns the vector assembled on the true dofs.
194  void ParallelAssemble(Vector &tv) const;
195 
196  /// Returns the vector assembled on the true dofs.
197  void ParallelAssemble(HypreParVector &tv) const;
198 
199  /// Returns a new vector assembled on the true dofs.
201 
202  void ExchangeFaceNbrData();
204  const Vector &FaceNbrData() const { return face_nbr_data; }
205 
206  // Redefine to handle the case when i is a face-neighbor element
207  virtual double GetValue(int i, const IntegrationPoint &ip,
208  int vdim = 1) const;
210  { return GetValue(T.ElementNo, T.GetIntPoint()); }
211 
212  // Redefine to handle the case when T describes a face-neighbor element
213  virtual double GetValue(ElementTransformation &T, const IntegrationPoint &ip,
214  int comp = 0, Vector *tr = NULL) const;
215 
216  virtual void GetVectorValue(int i, const IntegrationPoint &ip,
217  Vector &val) const;
218 
219  // Redefine to handle the case when T describes a face-neighbor element
220  virtual void GetVectorValue(ElementTransformation &T,
221  const IntegrationPoint &ip,
222  Vector &val, Vector *tr = NULL) const;
223 
225  virtual void ProjectCoefficient(Coefficient &coeff);
226 
228  /** @brief Project a discontinuous vector coefficient as a grid function on
229  a continuous finite element space. The values in shared dofs are
230  determined from the element with maximal attribute. */
231  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff);
232 
233  virtual void ProjectDiscCoefficient(Coefficient &coeff, AvgType type);
234 
235  virtual void ProjectDiscCoefficient(VectorCoefficient &vcoeff, AvgType type);
236 
238 
239  // Only the values in the master are guaranteed to be correct!
241  Array<int> &attr)
242  { ProjectBdrCoefficient(NULL, &vcoeff, attr); }
243 
244  // Only the values in the master are guaranteed to be correct!
245  virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array<int> &attr)
246  { ProjectBdrCoefficient(coeff, NULL, attr); }
247 
248  // Only the values in the master are guaranteed to be correct!
249  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff,
250  Array<int> &bdr_attr);
251 
252  virtual double ComputeL1Error(Coefficient *exsol[],
253  const IntegrationRule *irs[] = NULL) const
254  {
256  *exsol, NULL, 1, NULL, irs), pfes->GetComm());
257  }
258 
259  virtual double ComputeL1Error(Coefficient &exsol,
260  const IntegrationRule *irs[] = NULL) const
261  { return ComputeLpError(1.0, exsol, NULL, irs); }
262 
263  virtual double ComputeL1Error(VectorCoefficient &exsol,
264  const IntegrationRule *irs[] = NULL) const
265  { return ComputeLpError(1.0, exsol, NULL, NULL, irs); }
266 
267  virtual double ComputeL2Error(Coefficient *exsol[],
268  const IntegrationRule *irs[] = NULL) const
269  {
270  return GlobalLpNorm(2.0, GridFunction::ComputeL2Error(exsol, irs),
271  pfes->GetComm());
272  }
273 
274  virtual double ComputeL2Error(Coefficient &exsol,
275  const IntegrationRule *irs[] = NULL) const
276  { return ComputeLpError(2.0, exsol, NULL, irs); }
277 
278  virtual double ComputeL2Error(VectorCoefficient &exsol,
279  const IntegrationRule *irs[] = NULL,
280  Array<int> *elems = NULL) const
281  {
282  return GlobalLpNorm(2.0, GridFunction::ComputeL2Error(exsol, irs, elems),
283  pfes->GetComm());
284  }
285 
286  /// Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements
287  virtual double ComputeGradError(VectorCoefficient *exgrad,
288  const IntegrationRule *irs[] = NULL) const
289  {
290  return GlobalLpNorm(2.0, GridFunction::ComputeGradError(exgrad,irs),
291  pfes->GetComm());
292  }
293 
294  /// Returns ||curl u_ex - curl u_h||_L2 for ND elements
295  virtual double ComputeCurlError(VectorCoefficient *excurl,
296  const IntegrationRule *irs[] = NULL) const
297  {
298  return GlobalLpNorm(2.0, GridFunction::ComputeCurlError(excurl,irs),
299  pfes->GetComm());
300  }
301 
302  /// Returns ||div u_ex - div u_h||_L2 for RT elements
303  virtual double ComputeDivError(Coefficient *exdiv,
304  const IntegrationRule *irs[] = NULL) const
305  {
306  return GlobalLpNorm(2.0, GridFunction::ComputeDivError(exdiv,irs),
307  pfes->GetComm());
308  }
309 
310  /// Returns the Face Jumps error for L2 elements
311  virtual double ComputeDGFaceJumpError(Coefficient *exsol,
312  Coefficient *ell_coeff,
313  double Nu,
314  const IntegrationRule *irs[]=NULL)
315  const;
316 
317  /// Returns either the H1-seminorm or the DG Face Jumps error or both
318  /// depending on norm_type = 1, 2, 3
319  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
320  Coefficient *ell_coef, double Nu,
321  int norm_type) const
322  {
323  return GlobalLpNorm(2.0,
324  GridFunction::ComputeH1Error(exsol,exgrad,ell_coef,
325  Nu, norm_type),
326  pfes->GetComm());
327  }
328 
329  /// Returns the error measured in H1-norm for H1 elements or in "broken"
330  /// H1-norm for L2 elements
331  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
332  const IntegrationRule *irs[] = NULL) const
333  {
334  return GlobalLpNorm(2.0, GridFunction::ComputeH1Error(exsol,exgrad,irs),
335  pfes->GetComm());
336  }
337 
338  /// Returns the error measured H(div)-norm for RT elements
339  virtual double ComputeHDivError(VectorCoefficient *exsol,
340  Coefficient *exdiv,
341  const IntegrationRule *irs[] = NULL) const
342  {
343  return GlobalLpNorm(2.0, GridFunction::ComputeHDivError(exsol,exdiv,irs),
344  pfes->GetComm());
345  }
346 
347  /// Returns the error measured H(curl)-norm for ND elements
348  virtual double ComputeHCurlError(VectorCoefficient *exsol,
349  VectorCoefficient *excurl,
350  const IntegrationRule *irs[] = NULL) const
351  {
352  return GlobalLpNorm(2.0,
353  GridFunction::ComputeHCurlError(exsol,excurl,irs),
354  pfes->GetComm());
355  }
356 
357  virtual double ComputeMaxError(Coefficient *exsol[],
358  const IntegrationRule *irs[] = NULL) const
359  {
360  return GlobalLpNorm(infinity(),
361  GridFunction::ComputeMaxError(exsol, irs),
362  pfes->GetComm());
363  }
364 
365  virtual double ComputeMaxError(Coefficient &exsol,
366  const IntegrationRule *irs[] = NULL) const
367  {
368  return ComputeLpError(infinity(), exsol, NULL, irs);
369  }
370 
371  virtual double ComputeMaxError(VectorCoefficient &exsol,
372  const IntegrationRule *irs[] = NULL) const
373  {
374  return ComputeLpError(infinity(), exsol, NULL, NULL, irs);
375  }
376 
377  virtual double ComputeLpError(const double p, Coefficient &exsol,
378  Coefficient *weight = NULL,
379  const IntegrationRule *irs[] = NULL) const
380  {
382  p, exsol, weight, irs), pfes->GetComm());
383  }
384 
385  /** When given a vector weight, compute the pointwise (scalar) error as the
386  dot product of the vector error with the vector weight. Otherwise, the
387  scalar error is the l_2 norm of the vector error. */
388  virtual double ComputeLpError(const double p, VectorCoefficient &exsol,
389  Coefficient *weight = NULL,
390  VectorCoefficient *v_weight = NULL,
391  const IntegrationRule *irs[] = NULL) const
392  {
394  p, exsol, weight, v_weight, irs), pfes->GetComm());
395  }
396 
397  virtual void ComputeFlux(BilinearFormIntegrator &blfi,
398  GridFunction &flux,
399  bool wcoef = true, int subdomain = -1);
400 
401  /** Save the local portion of the ParGridFunction. This differs from the
402  serial GridFunction::Save in that it takes into account the signs of
403  the local dofs. */
404  virtual void Save(std::ostream &out) const;
405 
406 #ifdef MFEM_USE_ADIOS2
407  /** Save the local portion of the ParGridFunction. This differs from the
408  serial GridFunction::Save in that it takes into account the signs of
409  the local dofs. */
410  virtual void Save(
411  adios2stream &out, const std::string &variable_name,
413 #endif
414 
415  /// Merge the local grid functions
416  void SaveAsOne(std::ostream &out = mfem::out);
417 
418  virtual ~ParGridFunction() { }
419 };
420 
421 
422 /** Performs a global L2 projection (through a HypreBoomerAMG solve) of flux
423  from supplied discontinuous space into supplied smooth (continuous, or at
424  least conforming) space, and computes the Lp norms of the differences
425  between them on each element. This is one approach to handling conforming
426  and non-conforming elements in parallel. Returns the total error estimate. */
427 double L2ZZErrorEstimator(BilinearFormIntegrator &flux_integrator,
428  const ParGridFunction &x,
429  ParFiniteElementSpace &smooth_flux_fes,
430  ParFiniteElementSpace &flux_fes,
431  Vector &errors, int norm_p = 2, double solver_tol = 1e-12,
432  int solver_max_it = 200);
433 
434 }
435 
436 #endif // MFEM_USE_MPI
437 
438 #endif
virtual double ComputeLpError(const double p, VectorCoefficient &exsol, Coefficient *weight=NULL, VectorCoefficient *v_weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:388
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:30
Memory< double > data
Definition: vector.hpp:55
void AddDistribute(double a, const Vector &tv)
Definition: pgridfunc.hpp:161
ParGridFunction(const ParGridFunction &orig)
Copy constructor. The internal vector face_nbr_data is not copied.
Definition: pgridfunc.hpp:53
HypreParVector * ParallelAssemble() const
Returns a new vector assembled on the true dofs.
Definition: pgridfunc.cpp:198
Base class for vector Coefficients that optionally depend on time and space.
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:431
HypreParVector * ParallelProject() const
Returns a new vector restricted to the true dofs.
Definition: pgridfunc.cpp:181
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
Definition: gridfunc.cpp:2701
virtual ~ParGridFunction()
Definition: pgridfunc.hpp:418
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
Definition: gridfunc.cpp:2387
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:167
void SaveAsOne(std::ostream &out=mfem::out)
Merge the local grid functions.
Definition: pgridfunc.cpp:855
ParGridFunction(ParFiniteElementSpace *pf)
Definition: pgridfunc.hpp:56
virtual double GetValue(int i, const IntegrationPoint &ip, int vdim=1) const
Definition: pgridfunc.cpp:264
virtual void Save(std::ostream &out) const
Definition: pgridfunc.cpp:819
virtual void ProjectDiscCoefficient(VectorCoefficient &coeff)
Project a discontinuous vector coefficient as a grid function on a continuous finite element space...
Definition: pgridfunc.cpp:495
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 double ComputeHDivError(VectorCoefficient *exsol, Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns the error measured in H(div)-norm for RT elements.
Definition: gridfunc.cpp:2909
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: pgridfunc.cpp:474
Vector send_data
Vector used as an MPI buffer to send face-neighbor data in ExchangeFaceNbrData() to neighboring proce...
Definition: pgridfunc.hpp:44
const IntegrationPoint & GetIntPoint()
Get a const reference to the currently set integration point. This will return NULL if no integration...
Definition: eltrans.hpp:90
virtual double ComputeDGFaceJumpError(Coefficient *exsol, Coefficient *ell_coeff, double Nu, const IntegrationRule *irs[]=NULL) const
Returns the Face Jumps error for L2 elements.
Definition: pgridfunc.cpp:658
void Distribute(const Vector &tv)
Definition: pgridfunc.hpp:159
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:240
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
Definition: pgridfunc.hpp:295
virtual double ComputeHDivError(VectorCoefficient *exsol, Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns the error measured H(div)-norm for RT elements.
Definition: pgridfunc.hpp:339
double GetValue(ElementTransformation &T)
Definition: pgridfunc.hpp:209
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:101
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
Definition: pgridfunc.hpp:164
const Vector & FaceNbrData() const
Definition: pgridfunc.hpp:204
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: pgridfunc.cpp:81
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
Definition: gridfunc.cpp:2887
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:263
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:321
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:377
GridFunction & operator=(const GridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: gridfunc.hpp:109
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3089
ParGridFunction & operator=(const Vector &v)
Copy the data from a Vector to the ParGridFunction data.
Definition: pgridfunc.hpp:105
Vector & FaceNbrData()
Definition: pgridfunc.hpp:203
virtual void GetVectorValue(int i, const IntegrationPoint &ip, Vector &val) const
Definition: pgridfunc.cpp:321
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:70
virtual double ComputeL1Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:252
virtual double ComputeGradError(VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements.
Definition: gridfunc.cpp:2661
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:989
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:31
ParGridFunction & operator=(const ParGridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: pgridfunc.hpp:97
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:274
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Definition: pgridfunc.cpp:582
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:87
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
virtual void ComputeFlux(BilinearFormIntegrator &blfi, GridFunction &flux, bool wcoef=true, int subdomain=-1)
Definition: pgridfunc.cpp:1024
virtual double ComputeMaxError(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:357
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:486
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:1061
void Distribute(const Vector *tv)
Definition: pgridfunc.cpp:132
virtual double ComputeGradError(VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements.
Definition: pgridfunc.hpp:287
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Definition: gridfunc.cpp:2743
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:365
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:371
virtual double ComputeHCurlError(VectorCoefficient *exsol, VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns the error measured in H(curl)-norm for ND elements.
Definition: gridfunc.cpp:2918
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:245
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:267
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:331
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:620
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:45
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: gridfunc.cpp:2252
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:259
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:2983
virtual double ComputeL2Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL, Array< int > *elems=NULL) const
Definition: pgridfunc.hpp:278
Vector data type.
Definition: vector.hpp:51
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
Definition: pgridfunc.hpp:319
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Definition: pgridfunc.hpp:303
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:65
virtual double ComputeHCurlError(VectorCoefficient *exsol, VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns the error measured H(curl)-norm for ND elements.
Definition: pgridfunc.hpp:348
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:399
double f(const Vector &p)
ParFiniteElementSpace * ParFESpace() const
Definition: pgridfunc.hpp:108