MFEM  v4.4.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-2022, 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  /** @brief Construct a ParGridFunction using previously allocated Vector
69  @a base starting at the given offset, @a base_offset. */
70  ParGridFunction(ParFiniteElementSpace *pf, Vector &base, int base_offset = 0)
71  : GridFunction(pf, base, base_offset), pfes(pf) { }
72 
73  /// Construct a ParGridFunction using a GridFunction as external data.
74  /** The parallel space @a *pf and the space used by @a *gf should match. The
75  data from @a *gf is used as the local data of the ParGridFunction on each
76  processor. The ParGridFunction does not assume ownership of the data. */
78 
79  /** @brief Creates grid function on (all) dofs from a given vector on the
80  true dofs, i.e. P tv. */
82 
83  /** @brief Construct a local ParGridFunction from the given *global*
84  GridFunction. If @a partitioning is NULL (default), the data from @a gf
85  is NOT copied. */
86  ParGridFunction(ParMesh *pmesh, const GridFunction *gf,
87  const int *partitioning = NULL);
88 
89  /** @brief Construct a ParGridFunction on a given ParMesh, @a pmesh, reading
90  from an std::istream.
91 
92  In the process, a ParFiniteElementSpace and a FiniteElementCollection are
93  constructed. The new ParGridFunction assumes ownership of both. */
94  ParGridFunction(ParMesh *pmesh, std::istream &input);
95 
96  /// Copy assignment. Only the data of the base class Vector is copied.
97  /** It is assumed that this object and @a rhs use ParFiniteElementSpace%s
98  that have the same size.
99 
100  @note Defining this method overwrites the implicitly defined copy
101  assignment operator. */
103  { return operator=((const Vector &)rhs); }
104 
105  /// Assign constant values to the ParGridFunction data.
106  ParGridFunction &operator=(double value)
107  { GridFunction::operator=(value); return *this; }
108 
109  /// Copy the data from a Vector to the ParGridFunction data.
111  { GridFunction::operator=(v); return *this; }
112 
113  ParFiniteElementSpace *ParFESpace() const { return pfes; }
114 
115  virtual void Update();
116 
117  /// Associate a new FiniteElementSpace with the ParGridFunction.
118  /** The ParGridFunction is resized using the SetSize() method. The new space
119  @a f is expected to be a ParFiniteElementSpace. */
120  virtual void SetSpace(FiniteElementSpace *f);
121 
122  /// Associate a new parallel space with the ParGridFunction.
124 
125  using GridFunction::MakeRef;
126 
127  /** @brief Make the ParGridFunction reference external data on a new
128  FiniteElementSpace. */
129  /** This method changes the FiniteElementSpace associated with the
130  ParGridFunction and sets the pointer @a v as external data in the
131  ParGridFunction. The new space @a f is expected to be a
132  ParFiniteElementSpace. */
133  virtual void MakeRef(FiniteElementSpace *f, double *v);
134 
135  /** @brief Make the ParGridFunction reference external data on a new
136  ParFiniteElementSpace. */
137  /** This method changes the ParFiniteElementSpace associated with the
138  ParGridFunction and sets the pointer @a v as external data in the
139  ParGridFunction. */
140  void MakeRef(ParFiniteElementSpace *f, double *v);
141 
142  /** @brief Make the ParGridFunction reference external data on a new
143  FiniteElementSpace. */
144  /** This method changes the FiniteElementSpace associated with the
145  ParGridFunction and sets the data of the Vector @a v (plus the @a
146  v_offset) as external data in the ParGridFunction. The new space @a f is
147  expected to be a ParFiniteElementSpace.
148  @note This version of the method will also perform bounds checks when
149  the build option MFEM_DEBUG is enabled. */
150  virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset);
151 
152  /** @brief Make the ParGridFunction reference external data on a new
153  ParFiniteElementSpace. */
154  /** This method changes the ParFiniteElementSpace associated with the
155  ParGridFunction and sets the data of the Vector @a v (plus the
156  @a v_offset) as external data in the ParGridFunction.
157  @note This version of the method will also perform bounds checks when
158  the build option MFEM_DEBUG is enabled. */
159  void MakeRef(ParFiniteElementSpace *f, Vector &v, int v_offset);
160 
161  /** Set the grid function on (all) dofs from a given vector on the
162  true dofs, i.e. P tv. */
163  void Distribute(const Vector *tv);
164  void Distribute(const Vector &tv) { Distribute(&tv); }
165  void AddDistribute(double a, const Vector *tv);
166  void AddDistribute(double a, const Vector &tv) { AddDistribute(a, &tv); }
167 
168  /// Set the GridFunction from the given true-dof vector.
169  virtual void SetFromTrueDofs(const Vector &tv) { Distribute(tv); }
170 
171  /// Short semantic for Distribute()
173  { Distribute(&tv); return (*this); }
174 
176 
177  /// Returns the true dofs in a new HypreParVector
178  HypreParVector *GetTrueDofs() const;
179 
180  /// Returns the vector averaged on the true dofs.
181  void ParallelAverage(Vector &tv) const;
182 
183  /// Returns the vector averaged on the true dofs.
184  void ParallelAverage(HypreParVector &tv) const;
185 
186  /// Returns a new vector averaged on the true dofs.
188 
189  /// Returns the vector restricted to the true dofs.
190  void ParallelProject(Vector &tv) const;
191 
192  /// Returns the vector restricted to the true dofs.
193  void ParallelProject(HypreParVector &tv) const;
194 
195  /// Returns a new vector restricted to the true dofs.
197 
198  /// Returns the vector assembled on the true dofs.
199  void ParallelAssemble(Vector &tv) const;
200 
201  /// Returns the vector assembled on the true dofs.
202  void ParallelAssemble(HypreParVector &tv) const;
203 
204  /// Returns a new vector assembled on the true dofs.
206 
207  void ExchangeFaceNbrData();
209  const Vector &FaceNbrData() const { return face_nbr_data; }
210 
211  // Redefine to handle the case when i is a face-neighbor element
212  virtual double GetValue(int i, const IntegrationPoint &ip,
213  int vdim = 1) const;
215  { return GetValue(T.ElementNo, T.GetIntPoint()); }
216 
217  // Redefine to handle the case when T describes a face-neighbor element
218  virtual double GetValue(ElementTransformation &T, const IntegrationPoint &ip,
219  int comp = 0, Vector *tr = NULL) const;
220 
221  virtual void GetVectorValue(int i, const IntegrationPoint &ip,
222  Vector &val) const;
223 
224  // Redefine to handle the case when T describes a face-neighbor element
225  virtual void GetVectorValue(ElementTransformation &T,
226  const IntegrationPoint &ip,
227  Vector &val, Vector *tr = NULL) const;
228 
229  /// Parallel version of GridFunction::GetDerivative(); see its documentation.
230  void GetDerivative(int comp, int der_comp, ParGridFunction &der);
231 
232  /** Sets the output vector @a dof_vals to the values of the degrees of
233  freedom of element @a el. If @a el is greater than or equal to the number
234  of local elements, it will be interpreted as a shifted index of a face
235  neighbor element. */
236  virtual void GetElementDofValues(int el, Vector &dof_vals) const;
237 
239  virtual void ProjectCoefficient(Coefficient &coeff);
240 
242  /** @brief Project a discontinuous vector coefficient as a grid function on
243  a continuous finite element space. The values in shared dofs are
244  determined from the element with maximal attribute. */
245  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff);
246 
247  virtual void ProjectDiscCoefficient(Coefficient &coeff, AvgType type);
248 
249  virtual void ProjectDiscCoefficient(VectorCoefficient &vcoeff, AvgType type);
250 
252 
253  // Only the values in the master are guaranteed to be correct!
255  Array<int> &attr)
256  { ProjectBdrCoefficient(NULL, &vcoeff, attr); }
257 
258  // Only the values in the master are guaranteed to be correct!
259  virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array<int> &attr)
260  { ProjectBdrCoefficient(coeff, NULL, attr); }
261 
262  // Only the values in the master are guaranteed to be correct!
263  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff,
264  Array<int> &bdr_attr);
265 
266  virtual double ComputeL1Error(Coefficient *exsol[],
267  const IntegrationRule *irs[] = NULL) const
268  {
270  *exsol, NULL, 1, NULL, irs), pfes->GetComm());
271  }
272 
273  virtual double ComputeL1Error(Coefficient &exsol,
274  const IntegrationRule *irs[] = NULL) const
275  { return ComputeLpError(1.0, exsol, NULL, irs); }
276 
277  virtual double ComputeL1Error(VectorCoefficient &exsol,
278  const IntegrationRule *irs[] = NULL) const
279  { return ComputeLpError(1.0, exsol, NULL, NULL, irs); }
280 
281  virtual double ComputeL2Error(Coefficient *exsol[],
282  const IntegrationRule *irs[] = NULL) const
283  {
284  return GlobalLpNorm(2.0, GridFunction::ComputeL2Error(exsol, irs),
285  pfes->GetComm());
286  }
287 
288  virtual double ComputeL2Error(Coefficient &exsol,
289  const IntegrationRule *irs[] = NULL) const
290  { return ComputeLpError(2.0, exsol, NULL, irs); }
291 
292  virtual double ComputeL2Error(VectorCoefficient &exsol,
293  const IntegrationRule *irs[] = NULL,
294  Array<int> *elems = NULL) const
295  {
296  return GlobalLpNorm(2.0, GridFunction::ComputeL2Error(exsol, irs, elems),
297  pfes->GetComm());
298  }
299 
300  /// Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements
301  virtual double ComputeGradError(VectorCoefficient *exgrad,
302  const IntegrationRule *irs[] = NULL) const
303  {
304  return GlobalLpNorm(2.0, GridFunction::ComputeGradError(exgrad,irs),
305  pfes->GetComm());
306  }
307 
308  /// Returns ||curl u_ex - curl u_h||_L2 for ND elements
309  virtual double ComputeCurlError(VectorCoefficient *excurl,
310  const IntegrationRule *irs[] = NULL) const
311  {
312  return GlobalLpNorm(2.0, GridFunction::ComputeCurlError(excurl,irs),
313  pfes->GetComm());
314  }
315 
316  /// Returns ||div u_ex - div u_h||_L2 for RT elements
317  virtual double ComputeDivError(Coefficient *exdiv,
318  const IntegrationRule *irs[] = NULL) const
319  {
320  return GlobalLpNorm(2.0, GridFunction::ComputeDivError(exdiv,irs),
321  pfes->GetComm());
322  }
323 
324  /// Returns the Face Jumps error for L2 elements
325  virtual double ComputeDGFaceJumpError(Coefficient *exsol,
326  Coefficient *ell_coeff,
327  JumpScaling jump_scaling,
328  const IntegrationRule *irs[]=NULL)
329  const;
330 
331  /// Returns either the H1-seminorm or the DG Face Jumps error or both
332  /// depending on norm_type = 1, 2, 3
333  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
334  Coefficient *ell_coef, double Nu,
335  int norm_type) const
336  {
337  return GlobalLpNorm(2.0,
338  GridFunction::ComputeH1Error(exsol,exgrad,ell_coef,
339  Nu, norm_type),
340  pfes->GetComm());
341  }
342 
343  /// Returns the error measured in H1-norm for H1 elements or in "broken"
344  /// H1-norm for L2 elements
345  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
346  const IntegrationRule *irs[] = NULL) const
347  {
348  return GlobalLpNorm(2.0, GridFunction::ComputeH1Error(exsol,exgrad,irs),
349  pfes->GetComm());
350  }
351 
352  /// Returns the error measured H(div)-norm for RT elements
353  virtual double ComputeHDivError(VectorCoefficient *exsol,
354  Coefficient *exdiv,
355  const IntegrationRule *irs[] = NULL) const
356  {
357  return GlobalLpNorm(2.0, GridFunction::ComputeHDivError(exsol,exdiv,irs),
358  pfes->GetComm());
359  }
360 
361  /// Returns the error measured H(curl)-norm for ND elements
362  virtual double ComputeHCurlError(VectorCoefficient *exsol,
363  VectorCoefficient *excurl,
364  const IntegrationRule *irs[] = NULL) const
365  {
366  return GlobalLpNorm(2.0,
367  GridFunction::ComputeHCurlError(exsol,excurl,irs),
368  pfes->GetComm());
369  }
370 
371  virtual double ComputeMaxError(Coefficient *exsol[],
372  const IntegrationRule *irs[] = NULL) const
373  {
374  return GlobalLpNorm(infinity(),
375  GridFunction::ComputeMaxError(exsol, irs),
376  pfes->GetComm());
377  }
378 
379  virtual double ComputeMaxError(Coefficient &exsol,
380  const IntegrationRule *irs[] = NULL) const
381  {
382  return ComputeLpError(infinity(), exsol, NULL, irs);
383  }
384 
385  virtual double ComputeMaxError(VectorCoefficient &exsol,
386  const IntegrationRule *irs[] = NULL) const
387  {
388  return ComputeLpError(infinity(), exsol, NULL, NULL, irs);
389  }
390 
391  virtual double ComputeLpError(const double p, Coefficient &exsol,
392  Coefficient *weight = NULL,
393  const IntegrationRule *irs[] = NULL) const
394  {
396  p, exsol, weight, irs), pfes->GetComm());
397  }
398 
399  /** When given a vector weight, compute the pointwise (scalar) error as the
400  dot product of the vector error with the vector weight. Otherwise, the
401  scalar error is the l_2 norm of the vector error. */
402  virtual double ComputeLpError(const double p, VectorCoefficient &exsol,
403  Coefficient *weight = NULL,
404  VectorCoefficient *v_weight = NULL,
405  const IntegrationRule *irs[] = NULL) const
406  {
408  p, exsol, weight, v_weight, irs), pfes->GetComm());
409  }
410 
411  virtual void ComputeFlux(BilinearFormIntegrator &blfi,
412  GridFunction &flux,
413  bool wcoef = true, int subdomain = -1);
414 
415  /** Save the local portion of the ParGridFunction. This differs from the
416  serial GridFunction::Save in that it takes into account the signs of
417  the local dofs. */
418  virtual void Save(std::ostream &out) const;
419 
420  /// Save the ParGridFunction to a single file (written using MPI rank 0). The
421  /// given @a precision will be used for ASCII output.
422  void SaveAsOne(const char *fname, int precision=16) const;
423 
424  /// Save the ParGridFunction to files (one for each MPI rank). The files will
425  /// be given suffixes according to the MPI rank. The given @a precision will
426  /// be used for ASCII output.
427  virtual void Save(const char *fname, int precision=16) const;
428 
429 #ifdef MFEM_USE_ADIOS2
430  /** Save the local portion of the ParGridFunction. This differs from the
431  serial GridFunction::Save in that it takes into account the signs of
432  the local dofs. */
433  virtual void Save(
434  adios2stream &out, const std::string &variable_name,
436 #endif
437 
438  /// Merge the local grid functions
439  void SaveAsOne(std::ostream &out = mfem::out) const;
440 
441  virtual ~ParGridFunction() { }
442 };
443 
444 
445 /** Performs a global L2 projection (through a HypreBoomerAMG solve) of flux
446  from supplied discontinuous space into supplied smooth (continuous, or at
447  least conforming) space, and computes the Lp norms of the differences
448  between them on each element. This is one approach to handling conforming
449  and non-conforming elements in parallel. Returns the total error estimate. */
450 double L2ZZErrorEstimator(BilinearFormIntegrator &flux_integrator,
451  const ParGridFunction &x,
452  ParFiniteElementSpace &smooth_flux_fes,
453  ParFiniteElementSpace &flux_fes,
454  Vector &errors, int norm_p = 2, double solver_tol = 1e-12,
455  int solver_max_it = 200);
456 
457 }
458 
459 #endif // MFEM_USE_MPI
460 
461 #endif
virtual double ComputeLpError(const double p, VectorCoefficient &exsol, Coefficient *weight=NULL, VectorCoefficient *v_weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:402
ParGridFunction(ParFiniteElementSpace *pf, Vector &base, int base_offset=0)
Construct a ParGridFunction using previously allocated Vector base starting at the given offset...
Definition: pgridfunc.hpp:70
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
void GetDerivative(int comp, int der_comp, ParGridFunction &der)
Parallel version of GridFunction::GetDerivative(); see its documentation.
Definition: pgridfunc.cpp:485
Memory< double > data
Definition: vector.hpp:64
void AddDistribute(double a, const Vector &tv)
Definition: pgridfunc.hpp:166
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:479
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:2885
virtual void GetElementDofValues(int el, Vector &dof_vals) const
Definition: pgridfunc.cpp:506
virtual ~ParGridFunction()
Definition: pgridfunc.hpp:441
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
Definition: gridfunc.cpp:2569
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:172
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:873
virtual void ProjectDiscCoefficient(VectorCoefficient &coeff)
Project a discontinuous vector coefficient as a grid function on a continuous finite element space...
Definition: pgridfunc.cpp:546
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:3111
virtual void ProjectCoefficient(Coefficient &coeff)
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
Definition: pgridfunc.cpp:525
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:98
void Distribute(const Vector &tv)
Definition: pgridfunc.hpp:164
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:254
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
Definition: pgridfunc.hpp:309
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:353
double GetValue(ElementTransformation &T)
Definition: pgridfunc.hpp:214
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:273
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
ParGridFunction & operator=(double value)
Assign constant values to the ParGridFunction data.
Definition: pgridfunc.hpp:106
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
Definition: pgridfunc.hpp:169
const Vector & FaceNbrData() const
Definition: pgridfunc.hpp:209
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:3085
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:277
virtual double ComputeDGFaceJumpError(Coefficient *exsol, Coefficient *ell_coeff, JumpScaling jump_scaling, const IntegrationRule *irs[]=NULL) const
Returns the Face Jumps error for L2 elements.
Definition: pgridfunc.cpp:709
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
Definition: gridfunc.cpp:201
void GetTrueDofs(Vector &tv) const
Extract the true-dofs from the GridFunction.
Definition: gridfunc.cpp:364
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:391
GridFunction & operator=(const GridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: gridfunc.hpp:118
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3291
ParGridFunction & operator=(const Vector &v)
Copy the data from a Vector to the ParGridFunction data.
Definition: pgridfunc.hpp:110
Vector & FaceNbrData()
Definition: pgridfunc.hpp:208
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:148
virtual double ComputeL1Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:266
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:2845
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:1065
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:34
ParGridFunction & operator=(const ParGridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: pgridfunc.hpp:102
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:288
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Definition: pgridfunc.cpp:633
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:88
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:1100
virtual double ComputeMaxError(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:371
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:543
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:1137
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:301
void SaveAsOne(const char *fname, int precision=16) const
Definition: pgridfunc.cpp:899
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Definition: gridfunc.cpp:2926
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:379
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:385
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:3120
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:259
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:281
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:345
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:671
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:46
virtual void ProjectCoefficient(Coefficient &coeff)
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
Definition: gridfunc.cpp:2395
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: pgridfunc.hpp:273
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3185
virtual double ComputeL2Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL, Array< int > *elems=NULL) const
Definition: pgridfunc.hpp:292
Vector data type.
Definition: vector.hpp:60
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
Definition: pgridfunc.hpp:333
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Definition: pgridfunc.hpp:317
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:362
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:447
ParFiniteElementSpace * ParFESpace() const
Definition: pgridfunc.hpp:113