MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
gridfunc.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_GRIDFUNC
13 #define MFEM_GRIDFUNC
14 
15 #include "../config/config.hpp"
16 #include "fespace.hpp"
17 #include "coefficient.hpp"
18 #include "bilininteg.hpp"
19 #ifdef MFEM_USE_ADIOS2
20 #include "../general/adios2stream.hpp"
21 #endif
22 #include <limits>
23 #include <ostream>
24 #include <string>
25 
26 namespace mfem
27 {
28 
29 /// Class for grid function - Vector with associated FE space.
30 class GridFunction : public Vector
31 {
32 protected:
33  /// FE space on which the grid function lives. Owned if #fec is not NULL.
35 
36  /** @brief Used when the grid function is read from a file. It can also be
37  set explicitly, see MakeOwner().
38 
39  If not NULL, this pointer is owned by the GridFunction. */
41 
42  long fes_sequence; // see FiniteElementSpace::sequence, Mesh::sequence
43 
44  /** Optional, internal true-dof vector: if the FiniteElementSpace #fes has a
45  non-trivial (i.e. not NULL) prolongation operator, this Vector may hold
46  associated true-dof values - either owned or external. */
48 
49  void SaveSTLTri(std::ostream &out, double p1[], double p2[], double p3[]);
50 
52 
53  // Project the delta coefficient without scaling and return the (local)
54  // integral of the projection.
55  void ProjectDeltaCoefficient(DeltaCoefficient &delta_coeff,
56  double &integral);
57 
58  // Sum fluxes to vertices and count element contributions
60  GridFunction &flux,
61  Array<int>& counts,
62  bool wcoef,
63  int subdomain);
64 
65  /** Project a discontinuous vector coefficient in a continuous space and
66  return in dof_attr the maximal attribute of the elements containing each
67  degree of freedom. */
68  void ProjectDiscCoefficient(VectorCoefficient &coeff, Array<int> &dof_attr);
69 
70  /// Loading helper.
71  void LegacyNCReorder();
72 
73  void Destroy();
74 
75 public:
76 
77  GridFunction() { fes = NULL; fec = NULL; fes_sequence = 0; UseDevice(true); }
78 
79  /// Copy constructor. The internal true-dof vector #t_vec is not copied.
81  : Vector(orig), fes(orig.fes), fec(NULL), fes_sequence(orig.fes_sequence)
82  { UseDevice(true); }
83 
84  /// Construct a GridFunction associated with the FiniteElementSpace @a *f.
86  { fes = f; fec = NULL; fes_sequence = f->GetSequence(); UseDevice(true); }
87 
88  /// Construct a GridFunction using previously allocated array @a data.
89  /** The GridFunction does not assume ownership of @a data which is assumed to
90  be of size at least `f->GetVSize()`. Similar to the Vector constructor
91  for externally allocated array, the pointer @a data can be NULL. The data
92  array can be replaced later using the method SetData().
93  */
95  : Vector(data, f->GetVSize())
96  { fes = f; fec = NULL; fes_sequence = f->GetSequence(); UseDevice(true); }
97 
98  /** @brief Construct a GridFunction using previously allocated Vector @a base
99  starting at the given offset, @a base_offset. */
100  GridFunction(FiniteElementSpace *f, Vector &base, int base_offset = 0)
101  : Vector(base, base_offset, f->GetVSize())
102  { fes = f; fec = NULL; fes_sequence = f->GetSequence(); UseDevice(true); }
103 
104  /// Construct a GridFunction on the given Mesh, using the data from @a input.
105  /** The content of @a input should be in the format created by the method
106  Save(). The reconstructed FiniteElementSpace and FiniteElementCollection
107  are owned by the GridFunction. */
108  GridFunction(Mesh *m, std::istream &input);
109 
110  GridFunction(Mesh *m, GridFunction *gf_array[], int num_pieces);
111 
112  /// Copy assignment. Only the data of the base class Vector is copied.
113  /** It is assumed that this object and @a rhs use FiniteElementSpace%s that
114  have the same size.
115 
116  @note Defining this method overwrites the implicitly defined copy
117  assignment operator. */
119  { return operator=((const Vector &)rhs); }
120 
121  /// Make the GridFunction the owner of #fec and #fes.
122  /** If the new FiniteElementCollection, @a fec_, is NULL, ownership of #fec
123  and #fes is taken away. */
124  void MakeOwner(FiniteElementCollection *fec_) { fec = fec_; }
125 
127 
128  int VectorDim() const;
129  int CurlDim() const;
130 
131  /// Read only access to the (optional) internal true-dof Vector.
132  const Vector &GetTrueVector() const
133  {
134  MFEM_VERIFY(t_vec.Size() > 0, "SetTrueVector() before GetTrueVector()");
135  return t_vec;
136  }
137  /// Read and write access to the (optional) internal true-dof Vector.
138  /** Note that @a t_vec is set if it is not allocated or set already.*/
140  { if (t_vec.Size() == 0) { SetTrueVector(); } return t_vec; }
141 
142  /// Extract the true-dofs from the GridFunction.
143  void GetTrueDofs(Vector &tv) const;
144 
145  /// Shortcut for calling GetTrueDofs() with GetTrueVector() as argument.
147 
148  /// Set the GridFunction from the given true-dof vector.
149  virtual void SetFromTrueDofs(const Vector &tv);
150 
151  /// Shortcut for calling SetFromTrueDofs() with GetTrueVector() as argument.
153 
154  /// Returns the values in the vertices of i'th element for dimension vdim.
155  void GetNodalValues(int i, Array<double> &nval, int vdim = 1) const;
156 
157  /** @name Element index Get Value Methods
158 
159  These methods take an element index and return the interpolated value of
160  the field at a given reference point within the element.
161 
162  @warning These methods retrieve and use the ElementTransformation object
163  from the mfem::Mesh. This can alter the state of the element
164  transformation object and can also lead to unexpected results when the
165  ElementTransformation object is already in use such as when these methods
166  are called from within an integration loop. Consider using
167  GetValue(ElementTransformation &T, ...) instead.
168  */
169  ///@{
170  /** Return a scalar value from within the given element. */
171  virtual double GetValue(int i, const IntegrationPoint &ip,
172  int vdim = 1) const;
173 
174  /** Return a vector value from within the given element. */
175  virtual void GetVectorValue(int i, const IntegrationPoint &ip,
176  Vector &val) const;
177  ///@}
178 
179  /** @name Element Index Get Values Methods
180 
181  These are convenience methods for repeatedly calling GetValue for
182  multiple points within a given element. The GetValues methods are
183  optimized and should perform better than repeatedly calling GetValue. The
184  GetVectorValues method simply calls GetVectorValue repeatedly.
185 
186  @warning These methods retrieve and use the ElementTransformation object
187  from the mfem::Mesh. This can alter the state of the element
188  transformation object and can also lead to unexpected results when the
189  ElementTransformation object is already in use such as when these methods
190  are called from within an integration loop. Consider using
191  GetValues(ElementTransformation &T, ...) instead.
192  */
193  ///@{
194  /** Compute a collection of scalar values from within the element indicated
195  by the index i. */
196  void GetValues(int i, const IntegrationRule &ir, Vector &vals,
197  int vdim = 1) const;
198 
199  /** Compute a collection of vector values from within the element indicated
200  by the index i. */
201  void GetValues(int i, const IntegrationRule &ir, Vector &vals,
202  DenseMatrix &tr, int vdim = 1) const;
203 
204  void GetVectorValues(int i, const IntegrationRule &ir,
205  DenseMatrix &vals, DenseMatrix &tr) const;
206  ///@}
207 
208  /** @name ElementTransformation Get Value Methods
209 
210  These member functions are designed for use within
211  GridFunctionCoefficient objects. These can be used with
212  ElementTransformation objects coming from either
213  Mesh::GetElementTransformation() or Mesh::GetBdrElementTransformation().
214 
215  @note These methods do not reset the ElementTransformation object so they
216  should be safe to use within integration loops or other contexts where
217  the ElementTransformation is already in use.
218  */
219  ///@{
220  /** Return a scalar value from within the element indicated by the
221  ElementTransformation Object. */
222  virtual double GetValue(ElementTransformation &T, const IntegrationPoint &ip,
223  int comp = 0, Vector *tr = NULL) const;
224 
225  /** Return a vector value from within the element indicated by the
226  ElementTransformation Object. */
227  virtual void GetVectorValue(ElementTransformation &T,
228  const IntegrationPoint &ip,
229  Vector &val, Vector *tr = NULL) const;
230  ///@}
231 
232  /** @name ElementTransformation Get Values Methods
233 
234  These are convenience methods for repeatedly calling GetValue for
235  multiple points within a given element. They work by calling either the
236  ElementTransformation or FaceElementTransformations versions described
237  above. Consequently, these methods should not be expected to run faster
238  than calling the above methods in an external loop.
239 
240  @note These methods do not reset the ElementTransformation object so they
241  should be safe to use within integration loops or other contexts where
242  the ElementTransformation is already in use.
243 
244  @note These methods can also be used with FaceElementTransformations
245  objects.
246  */
247  ///@{
248  /** Compute a collection of scalar values from within the element indicated
249  by the ElementTransformation object. */
251  Vector &vals, int comp = 0, DenseMatrix *tr = NULL) const;
252 
253  /** Compute a collection of vector values from within the element indicated
254  by the ElementTransformation object. */
256  DenseMatrix &vals, DenseMatrix *tr = NULL) const;
257  ///@}
258 
259  /** @name Face Index Get Values Methods
260 
261  These methods are designed to work with Discontinuous Galerkin basis
262  functions. They compute field values on the interface between elements,
263  or on boundary elements, by interpolating the field in a neighboring
264  element. The \a side argument indices which neighboring element should be
265  used: 0, 1, or 2 (automatically chosen).
266 
267  @warning These methods retrieve and use the FaceElementTransformations
268  object from the mfem::Mesh. This can alter the state of the face element
269  transformations object and can also lead to unexpected results when the
270  FaceElementTransformations object is already in use such as when these
271  methods are called from within an integration loop. Consider using
272  GetValues(ElementTransformation &T, ...) instead.
273  */
274  ///@{
275  /** Compute a collection of scalar values from within the face
276  indicated by the index i. */
277  int GetFaceValues(int i, int side, const IntegrationRule &ir, Vector &vals,
278  DenseMatrix &tr, int vdim = 1) const;
279 
280  /** Compute a collection of vector values from within the face
281  indicated by the index i. */
282  int GetFaceVectorValues(int i, int side, const IntegrationRule &ir,
283  DenseMatrix &vals, DenseMatrix &tr) const;
284  ///@}
285 
286  void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps,
287  int vdim = 1) const;
288 
289  void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps,
290  DenseMatrix &tr, int vdim = 1) const;
291 
292  void GetHessians(int i, const IntegrationRule &ir, DenseMatrix &hess,
293  int vdim = 1) const;
294 
295  void GetHessians(int i, const IntegrationRule &ir, DenseMatrix &hess,
296  DenseMatrix &tr, int vdim = 1) const;
297 
298  void GetValuesFrom(const GridFunction &orig_func);
299 
300  void GetBdrValuesFrom(const GridFunction &orig_func);
301 
302  void GetVectorFieldValues(int i, const IntegrationRule &ir,
303  DenseMatrix &vals,
304  DenseMatrix &tr, int comp = 0) const;
305 
306  /// For a vector grid function, makes sure that the ordering is byNODES.
307  void ReorderByNodes();
308 
309  /// Return the values as a vector on mesh vertices for dimension vdim.
310  void GetNodalValues(Vector &nval, int vdim = 1) const;
311 
312  void GetVectorFieldNodalValues(Vector &val, int comp) const;
313 
314  void ProjectVectorFieldOn(GridFunction &vec_field, int comp = 0);
315 
316  /** @brief Compute a certain derivative of a function's component.
317  Derivatives of the function are computed at the DOF locations of @a der,
318  and averaged over overlapping DOFs. Thus this function projects the
319  derivative to the FiniteElementSpace of @a der.
320  @param[in] comp Index of the function's component to be differentiated.
321  The index is 1-based, i.e., use 1 for scalar functions.
322  @param[in] der_comp Use 0/1/2 for derivatives in x/y/z directions.
323  @param[out] der The resulting derivative (scalar function). The
324  FiniteElementSpace of this function must be set
325  before the call. */
326  void GetDerivative(int comp, int der_comp, GridFunction &der);
327 
328  double GetDivergence(ElementTransformation &tr) const;
329 
330  void GetCurl(ElementTransformation &tr, Vector &curl) const;
331 
332  void GetGradient(ElementTransformation &tr, Vector &grad) const;
333 
335  DenseMatrix &grad) const;
336 
337  void GetGradients(const int elem, const IntegrationRule &ir,
338  DenseMatrix &grad) const
339  { GetGradients(*fes->GetElementTransformation(elem), ir, grad); }
340 
341  void GetVectorGradient(ElementTransformation &tr, DenseMatrix &grad) const;
342 
343  /** Compute \f$ (\int_{\Omega} (*this) \psi_i)/(\int_{\Omega} \psi_i) \f$,
344  where \f$ \psi_i \f$ are the basis functions for the FE space of avgs.
345  Both FE spaces should be scalar and on the same mesh. */
346  void GetElementAverages(GridFunction &avgs) const;
347 
348  /** Sets the output vector @a dof_vals to the values of the degrees of
349  freedom of element @a el. */
350  virtual void GetElementDofValues(int el, Vector &dof_vals) const;
351 
352  /** Impose the given bounds on the function's DOFs while preserving its local
353  * integral (described in terms of the given weights) on the i'th element
354  * through SLBPQ optimization.
355  * Intended to be used for discontinuous FE functions. */
356  void ImposeBounds(int i, const Vector &weights,
357  const Vector &lo_, const Vector &hi_);
358  void ImposeBounds(int i, const Vector &weights,
359  double min_ = 0.0, double max_ = infinity());
360 
361  /** On a non-conforming mesh, make sure the function lies in the conforming
362  space by multiplying with R and then with P, the conforming restriction
363  and prolongation matrices of the space, respectively. */
364  void RestrictConforming();
365 
366  /** @brief Project the @a src GridFunction to @a this GridFunction, both of
367  which must be on the same mesh. */
368  /** The current implementation assumes that all elements use the same
369  projection matrix. */
370  void ProjectGridFunction(const GridFunction &src);
371 
372  /** @brief Project @a coeff Coefficient to @a this GridFunction. The
373  projection computation depends on the choice of the FiniteElementSpace
374  #fes. Note that this is usually interpolation at the degrees of freedom
375  in each element (not L2 projection). */
376  virtual void ProjectCoefficient(Coefficient &coeff);
377 
378  /** @brief Project @a coeff Coefficient to @a this GridFunction, using one
379  element for each degree of freedom in @a dofs and nodal interpolation on
380  that element. */
381  void ProjectCoefficient(Coefficient &coeff, Array<int> &dofs, int vd = 0);
382 
383  /** @brief Project @a vcoeff VectorCoefficient to @a this GridFunction. The
384  projection computation depends on the choice of the FiniteElementSpace
385  #fes. Note that this is usually interpolation at the degrees of freedom
386  in each element (not L2 projection).*/
388 
389  /** @brief Project @a vcoeff VectorCoefficient to @a this GridFunction, using
390  one element for each degree of freedom in @a dofs and nodal interpolation
391  on that element. */
392  void ProjectCoefficient(VectorCoefficient &vcoeff, Array<int> &dofs);
393 
394  /** @brief Project @a vcoeff VectorCoefficient to @a this GridFunction, only
395  projecting onto elements with the given @a attribute */
396  void ProjectCoefficient(VectorCoefficient &vcoeff, int attribute);
397 
398  /** @brief Analogous to the version with argument @a vcoeff VectorCoefficient
399  but using an array of scalar coefficients for each component. */
400  void ProjectCoefficient(Coefficient *coeff[]);
401 
402  /** @brief Project a discontinuous vector coefficient as a grid function on
403  a continuous finite element space. The values in shared dofs are
404  determined from the element with maximal attribute. */
405  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff);
406 
408  /** @brief Projects a discontinuous coefficient so that the values in shared
409  vdofs are computed by taking an average of the possible values. */
410  virtual void ProjectDiscCoefficient(Coefficient &coeff, AvgType type);
411  /** @brief Projects a discontinuous _vector_ coefficient so that the values
412  in shared vdofs are computed by taking an average of the possible values.
413  */
414  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff, AvgType type);
415 
416 protected:
417  /** @brief Accumulates (depending on @a type) the values of @a coeff at all
418  shared vdofs and counts in how many zones each vdof appears. */
419  void AccumulateAndCountZones(Coefficient &coeff, AvgType type,
420  Array<int> &zones_per_vdof);
421 
422  /** @brief Accumulates (depending on @a type) the values of @a vcoeff at all
423  shared vdofs and counts in how many zones each vdof appears. */
425  Array<int> &zones_per_vdof);
426 
427  /** @brief Used for the serial and parallel implementations of the
428  GetDerivative() method; see its documentation. */
429  void AccumulateAndCountDerivativeValues(int comp, int der_comp,
430  GridFunction &der,
431  Array<int> &zones_per_dof);
432 
434  VectorCoefficient *vcoeff, Array<int> &attr,
435  Array<int> &values_counter);
436 
438  Array<int> &bdr_attr,
439  Array<int> &values_counter);
440 
441  // Complete the computation of averages; called e.g. after
442  // AccumulateAndCountZones().
443  void ComputeMeans(AvgType type, Array<int> &zones_per_vdof);
444 
445 public:
446  /** @brief Project a Coefficient on the GridFunction, modifying only DOFs on
447  the boundary associated with the boundary attributes marked in the
448  @a attr array. */
450  {
451  Coefficient *coeff_p = &coeff;
452  ProjectBdrCoefficient(&coeff_p, attr);
453  }
454 
455  /** @brief Project a VectorCoefficient on the GridFunction, modifying only
456  DOFs on the boundary associated with the boundary attributes marked in
457  the @a attr array. */
458  virtual void ProjectBdrCoefficient(VectorCoefficient &vcoeff,
459  Array<int> &attr);
460 
461  /** @brief Project a set of Coefficient%s on the components of the
462  GridFunction, modifying only DOFs on the boundary associated with the
463  boundary attributed marked in the @a attr array. */
464  /** If a Coefficient pointer in the array @a coeff is NULL, that component
465  will not be touched. */
466  virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array<int> &attr);
467 
468  /** Project the normal component of the given VectorCoefficient on
469  the boundary. Only boundary attributes that are marked in
470  'bdr_attr' are projected. Assumes RT-type VectorFE GridFunction. */
472  Array<int> &bdr_attr);
473 
474  /** @brief Project the tangential components of the given VectorCoefficient
475  on the boundary. Only boundary attributes that are marked in @a bdr_attr
476  are projected. Assumes ND-type VectorFE GridFunction. */
477  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff,
478  Array<int> &bdr_attr);
479 
480 
481  virtual double ComputeL2Error(Coefficient *exsol[],
482  const IntegrationRule *irs[] = NULL,
483  const Array<int> *elems = NULL) const;
484 
485  /// Returns ||grad u_ex - grad u_h||_L2 in element ielem for H1 or L2 elements
486  virtual double ComputeElementGradError(int ielem, VectorCoefficient *exgrad,
487  const IntegrationRule *irs[] = NULL) const;
488 
489  /// Returns ||u_ex - u_h||_L2 for H1 or L2 elements
490  /* The @a elems input variable expects a list of markers:
491  an elem marker equal to 1 will compute the L2 error on that element
492  an elem marker equal to 0 will not compute the L2 error on that element */
493  virtual double ComputeL2Error(Coefficient &exsol,
494  const IntegrationRule *irs[] = NULL,
495  const Array<int> *elems = NULL) const
496  { return GridFunction::ComputeLpError(2.0, exsol, NULL, irs, elems); }
497 
498  virtual double ComputeL2Error(VectorCoefficient &exsol,
499  const IntegrationRule *irs[] = NULL,
500  const Array<int> *elems = NULL) const;
501 
502  /// Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements
503  virtual double ComputeGradError(VectorCoefficient *exgrad,
504  const IntegrationRule *irs[] = NULL) const;
505 
506  /// Returns ||curl u_ex - curl u_h||_L2 for ND elements
507  virtual double ComputeCurlError(VectorCoefficient *excurl,
508  const IntegrationRule *irs[] = NULL) const;
509 
510  /// Returns ||div u_ex - div u_h||_L2 for RT elements
511  virtual double ComputeDivError(Coefficient *exdiv,
512  const IntegrationRule *irs[] = NULL) const;
513 
514  /// Returns the Face Jumps error for L2 elements. The error can be weighted
515  /// by a constant nu, by nu/h, or nu*p^2/h, depending on the value of
516  /// @a jump_scaling.
517  virtual double ComputeDGFaceJumpError(Coefficient *exsol,
518  Coefficient *ell_coeff,
519  class JumpScaling jump_scaling,
520  const IntegrationRule *irs[] = NULL)
521  const;
522 
523  /// Returns the Face Jumps error for L2 elements, with 1/h scaling.
524  MFEM_DEPRECATED
525  double ComputeDGFaceJumpError(Coefficient *exsol,
526  Coefficient *ell_coeff,
527  double Nu,
528  const IntegrationRule *irs[] = NULL) const;
529 
530  /** This method is kept for backward compatibility.
531 
532  Returns either the H1-seminorm, or the DG face jumps error, or both
533  depending on norm_type = 1, 2, 3. Additional arguments for the DG face
534  jumps norm: ell_coeff: mesh-depended coefficient (weight) Nu: scalar
535  constant weight */
536  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
537  Coefficient *ell_coef, double Nu,
538  int norm_type) const;
539 
540  /// Returns the error measured in H1-norm for H1 elements or in "broken"
541  /// H1-norm for L2 elements
542  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
543  const IntegrationRule *irs[] = NULL) const;
544 
545  /// Returns the error measured in H(div)-norm for RT elements
546  virtual double ComputeHDivError(VectorCoefficient *exsol,
547  Coefficient *exdiv,
548  const IntegrationRule *irs[] = NULL) const;
549 
550  /// Returns the error measured in H(curl)-norm for ND elements
551  virtual double ComputeHCurlError(VectorCoefficient *exsol,
552  VectorCoefficient *excurl,
553  const IntegrationRule *irs[] = NULL) const;
554 
555  virtual double ComputeMaxError(Coefficient &exsol,
556  const IntegrationRule *irs[] = NULL) const
557  {
558  return ComputeLpError(infinity(), exsol, NULL, irs);
559  }
560 
561  virtual double ComputeMaxError(Coefficient *exsol[],
562  const IntegrationRule *irs[] = NULL) const;
563 
564  virtual double ComputeMaxError(VectorCoefficient &exsol,
565  const IntegrationRule *irs[] = NULL) const
566  {
567  return ComputeLpError(infinity(), exsol, NULL, NULL, irs);
568  }
569 
570  virtual double ComputeL1Error(Coefficient &exsol,
571  const IntegrationRule *irs[] = NULL) const
572  { return ComputeLpError(1.0, exsol, NULL, irs); }
573 
574  virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad,
575  int norm_type, const Array<int> *elems = NULL,
576  const IntegrationRule *irs[] = NULL) const;
577 
578  virtual double ComputeL1Error(VectorCoefficient &exsol,
579  const IntegrationRule *irs[] = NULL) const
580  { return ComputeLpError(1.0, exsol, NULL, NULL, irs); }
581 
582  /* The @a elems input variable expects a list of markers:
583  an elem marker equal to 1 will compute the L2 error on that element
584  an elem marker equal to 0 will not compute the L2 error on that element */
585  virtual double ComputeLpError(const double p, Coefficient &exsol,
586  Coefficient *weight = NULL,
587  const IntegrationRule *irs[] = NULL,
588  const Array<int> *elems = NULL) const;
589 
590  /** Compute the Lp error in each element of the mesh and store the results in
591  the Vector @a error. The result should be of length number of elements,
592  for example an L2 GridFunction of order zero using map type VALUE. */
593  virtual void ComputeElementLpErrors(const double p, Coefficient &exsol,
594  Vector &error,
595  Coefficient *weight = NULL,
596  const IntegrationRule *irs[] = NULL
597  ) const;
598 
599  virtual void ComputeElementL1Errors(Coefficient &exsol,
600  Vector &error,
601  const IntegrationRule *irs[] = NULL
602  ) const
603  { ComputeElementLpErrors(1.0, exsol, error, NULL, irs); }
604 
605  virtual void ComputeElementL2Errors(Coefficient &exsol,
606  Vector &error,
607  const IntegrationRule *irs[] = NULL
608  ) const
609  { ComputeElementLpErrors(2.0, exsol, error, NULL, irs); }
610 
612  Vector &error,
613  const IntegrationRule *irs[] = NULL
614  ) const
615  { ComputeElementLpErrors(infinity(), exsol, error, NULL, irs); }
616 
617  /** When given a vector weight, compute the pointwise (scalar) error as the
618  dot product of the vector error with the vector weight. Otherwise, the
619  scalar error is the l_2 norm of the vector error. */
620  virtual double ComputeLpError(const double p, VectorCoefficient &exsol,
621  Coefficient *weight = NULL,
622  VectorCoefficient *v_weight = NULL,
623  const IntegrationRule *irs[] = NULL) const;
624 
625  /** Compute the Lp error in each element of the mesh and store the results in
626  the Vector @ error. The result should be of length number of elements,
627  for example an L2 GridFunction of order zero using map type VALUE. */
628  virtual void ComputeElementLpErrors(const double p, VectorCoefficient &exsol,
629  Vector &error,
630  Coefficient *weight = NULL,
631  VectorCoefficient *v_weight = NULL,
632  const IntegrationRule *irs[] = NULL
633  ) const;
634 
636  Vector &error,
637  const IntegrationRule *irs[] = NULL
638  ) const
639  { ComputeElementLpErrors(1.0, exsol, error, NULL, NULL, irs); }
640 
642  Vector &error,
643  const IntegrationRule *irs[] = NULL
644  ) const
645  { ComputeElementLpErrors(2.0, exsol, error, NULL, NULL, irs); }
646 
648  Vector &error,
649  const IntegrationRule *irs[] = NULL
650  ) const
651  { ComputeElementLpErrors(infinity(), exsol, error, NULL, NULL, irs); }
652 
653  virtual void ComputeFlux(BilinearFormIntegrator &blfi,
654  GridFunction &flux,
655  bool wcoef = true, int subdomain = -1);
656 
657  /// Redefine '=' for GridFunction = constant.
658  GridFunction &operator=(double value);
659 
660  /// Copy the data from @a v.
661  /** The size of @a v must be equal to the size of the associated
662  FiniteElementSpace #fes. */
663  GridFunction &operator=(const Vector &v);
664 
665  /// Transform by the Space UpdateMatrix (e.g., on Mesh change).
666  virtual void Update();
667 
669  const FiniteElementSpace *FESpace() const { return fes; }
670 
671  /// Associate a new FiniteElementSpace with the GridFunction.
672  /** The GridFunction is resized using the SetSize() method. */
673  virtual void SetSpace(FiniteElementSpace *f);
674 
675  using Vector::MakeRef;
676 
677  /** @brief Make the GridFunction reference external data on a new
678  FiniteElementSpace. */
679  /** This method changes the FiniteElementSpace associated with the
680  GridFunction and sets the pointer @a v as external data in the
681  GridFunction. */
682  virtual void MakeRef(FiniteElementSpace *f, double *v);
683 
684  /** @brief Make the GridFunction reference external data on a new
685  FiniteElementSpace. */
686  /** This method changes the FiniteElementSpace associated with the
687  GridFunction and sets the data of the Vector @a v (plus the @a v_offset)
688  as external data in the GridFunction.
689  @note This version of the method will also perform bounds checks when
690  the build option MFEM_DEBUG is enabled. */
691  virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset);
692 
693  /** @brief Associate a new FiniteElementSpace and new true-dof data with the
694  GridFunction. */
695  /** - If the prolongation matrix of @a f is trivial (i.e. its method
696  FiniteElementSpace::GetProlongationMatrix() returns NULL), then the
697  method MakeRef() is called with the same arguments.
698  - Otherwise, the method SetSpace() is called with argument @a f.
699  - The internal true-dof vector is set to reference @a tv. */
700  void MakeTRef(FiniteElementSpace *f, double *tv);
701 
702  /** @brief Associate a new FiniteElementSpace and new true-dof data with the
703  GridFunction. */
704  /** - If the prolongation matrix of @a f is trivial (i.e. its method
705  FiniteElementSpace::GetProlongationMatrix() returns NULL), this method
706  calls MakeRef() with the same arguments.
707  - Otherwise, this method calls SetSpace() with argument @a f.
708  - The internal true-dof vector is set to reference the sub-vector of
709  @a tv starting at the offset @a tv_offset. */
710  void MakeTRef(FiniteElementSpace *f, Vector &tv, int tv_offset);
711 
712  /// Save the GridFunction to an output stream.
713  virtual void Save(std::ostream &out) const;
714 
715  /// Save the GridFunction to a file. The given @a precision will be used for
716  /// ASCII output.
717  virtual void Save(const char *fname, int precision=16) const;
718 
719 #ifdef MFEM_USE_ADIOS2
720  /// Save the GridFunction to a binary output stream using adios2 bp format.
721  virtual void Save(adios2stream &out, const std::string& variable_name,
724 #endif
725 
726  /** @brief Write the GridFunction in VTK format. Note that Mesh::PrintVTK
727  must be called first. The parameter ref > 0 must match the one used in
728  Mesh::PrintVTK. */
729  void SaveVTK(std::ostream &out, const std::string &field_name, int ref);
730 
731  /** @brief Write the GridFunction in STL format. Note that the mesh dimension
732  must be 2 and that quad elements will be broken into two triangles.*/
733  void SaveSTL(std::ostream &out, int TimesToRefine = 1);
734 
735  /// Destroys grid function.
736  virtual ~GridFunction() { Destroy(); }
737 };
738 
739 
740 /** Overload operator<< for std::ostream and GridFunction; valid also for the
741  derived class ParGridFunction */
742 std::ostream &operator<<(std::ostream &out, const GridFunction &sol);
743 
744 /// Class used to specify how the jump terms in
745 /// GridFunction::ComputeDGFaceJumpError are scaled.
747 {
748 public:
750  {
754  };
755 private:
756  double nu;
757  JumpScalingType type;
758 public:
759  JumpScaling(double nu_=1.0, JumpScalingType type_=CONSTANT)
760  : nu(nu_), type(type_) { }
761  double Eval(double h, int p) const
762  {
763  double val = nu;
764  if (type != CONSTANT) { val /= h; }
765  if (type == P_SQUARED_OVER_H) { val *= p*p; }
766  return val;
767  }
768 };
769 
770 /// Overload operator<< for std::ostream and QuadratureFunction.
771 std::ostream &operator<<(std::ostream &out, const QuadratureFunction &qf);
772 
773 
774 double ZZErrorEstimator(BilinearFormIntegrator &blfi,
775  GridFunction &u,
776  GridFunction &flux,
777  Vector &error_estimates,
778  Array<int> *aniso_flags = NULL,
779  int with_subdomains = 1,
780  bool with_coeff = false);
781 
782 /// Defines the global tensor product polynomial space used by NewZZErorrEstimator
783 /**
784  * See BoundingBox(...) for a description of @a angle and @a midpoint
785  */
786 void TensorProductLegendre(int dim, // input
787  int order, // input
788  const Vector &x_in, // input
789  const Vector &xmax, // input
790  const Vector &xmin, // input
791  Vector &poly, // output
792  double angle=0.0, // input (optional)
793  const Vector *midpoint=NULL); // input (optional)
794 
795 /// Defines the bounding box for the face patches used by NewZZErorrEstimator
796 /**
797  * By default, BoundingBox(...) computes the parameters of a minimal bounding box
798  * for the given @a face_patch that is aligned with the physical (i.e. global)
799  * Cartesian axes. This means that the size of the bounding box will depend on the
800  * orientation of the patch. It is better to construct an orientation-independent box.
801  * This is implemented for 2D patches. The parameters @a angle and @a midpoint encode
802  * the necessary additional geometric information.
803  *
804  * @a iface : Index of the face that the patch corresponds to.
805  * This is used to compute @a angle and @a midpoint.
806  *
807  * @a angle : The angle the patch face makes with the x-axis.
808  * @a midpoint : The midpoint of the face.
809  */
810 void BoundingBox(const Array<int> &face_patch, // input
811  FiniteElementSpace *ufes, // input
812  int order, // input
813  Vector &xmin, // output
814  Vector &xmax, // output
815  double &angle, // output
816  Vector &midpoint, // output
817  int iface=-1); // input (optional)
818 
819 /// A ``true'' ZZ error estimator that uses face-based patches for flux reconstruction.
820 /**
821  * Only two-element face patches are ever used:
822  * - For conforming faces, the face patch consists of its two neighboring elements.
823  * - In the non-conforming setting, only the face patches associated to fine-scale
824  * element faces are used. These face patches always consist of two elements
825  * delivered by mesh::GetFaceElements(Face, *Elem1, *Elem2).
826  */
827 double LSZZErrorEstimator(BilinearFormIntegrator &blfi, // input
828  GridFunction &u, // input
829  Vector &error_estimates, // output
830  bool subdomain_reconstruction = true, // input (optional)
831  bool with_coeff = false, // input (optional)
832  double tichonov_coeff = 0.0); // input (optional)
833 
834 /// Compute the Lp distance between two grid functions on the given element.
835 double ComputeElementLpDistance(double p, int i,
836  GridFunction& gf1, GridFunction& gf2);
837 
838 
839 /// Class used for extruding scalar GridFunctions
841 {
842 private:
843  int n;
844  Mesh *mesh_in;
845  Coefficient &sol_in;
846 public:
848  : n(n_), mesh_in(m), sol_in(s) { }
849  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip);
850  virtual ~ExtrudeCoefficient() { }
851 };
852 
853 /// Extrude a scalar 1D GridFunction, after extruding the mesh with Extrude1D.
854 GridFunction *Extrude1DGridFunction(Mesh *mesh, Mesh *mesh2d,
855  GridFunction *sol, const int ny);
856 
857 } // namespace mfem
858 
859 #endif
void SaveSTLTri(std::ostream &out, double p1[], double p2[], double p3[])
Definition: gridfunc.cpp:3796
virtual ~GridFunction()
Destroys grid function.
Definition: gridfunc.hpp:736
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
void GetElementAverages(GridFunction &avgs) const
Definition: gridfunc.cpp:1804
void GetVectorValues(int i, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr) const
Definition: gridfunc.cpp:706
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
Memory< double > data
Definition: vector.hpp:64
Class used for extruding scalar GridFunctions.
Definition: gridfunc.hpp:840
void SumFluxAndCount(BilinearFormIntegrator &blfi, GridFunction &flux, Array< int > &counts, bool wcoef, int subdomain)
Definition: gridfunc.cpp:252
void SetFromTrueVector()
Shortcut for calling SetFromTrueDofs() with GetTrueVector() as argument.
Definition: gridfunc.hpp:152
Base class for vector Coefficients that optionally depend on time and space.
void ProjectDeltaCoefficient(DeltaCoefficient &delta_coeff, double &integral)
Definition: gridfunc.cpp:2336
int CurlDim() const
Definition: gridfunc.cpp:347
virtual void ComputeElementLpErrors(const double p, Coefficient &exsol, Vector &error, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3413
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
Definition: gridfunc.cpp:2940
virtual void GetVectorValue(int i, const IntegrationPoint &ip, Vector &val) const
Definition: gridfunc.cpp:468
void RestrictConforming()
Definition: gridfunc.cpp:1982
GridFunction * Extrude1DGridFunction(Mesh *mesh, Mesh *mesh2d, GridFunction *sol, const int ny)
Extrude a scalar 1D GridFunction, after extruding the mesh with Extrude1D.
Definition: gridfunc.cpp:4497
void BoundingBox(const Array< int > &patch, FiniteElementSpace *ufes, int order, Vector &xmin, Vector &xmax, double &angle, Vector &midpoint, int iface)
Defines the bounding box for the face patches used by NewZZErorrEstimator.
Definition: gridfunc.cpp:4160
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
Definition: gridfunc.cpp:2573
int VectorDim() const
Definition: gridfunc.cpp:324
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
Delta function coefficient optionally multiplied by a weight coefficient and a scaled time dependent ...
int Size() const
Returns the size of the vector.
Definition: vector.hpp:200
void AccumulateAndCountBdrTangentValues(VectorCoefficient &vcoeff, Array< int > &bdr_attr, Array< int > &values_counter)
Definition: gridfunc.cpp:2259
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:3166
void GetVectorGradientHat(ElementTransformation &T, DenseMatrix &gh) const
Definition: gridfunc.cpp:1436
Vector & GetTrueVector()
Read and write access to the (optional) internal true-dof Vector.
Definition: gridfunc.hpp:139
void MakeOwner(FiniteElementCollection *fec_)
Make the GridFunction the owner of fec and fes.
Definition: gridfunc.hpp:124
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:564
ExtrudeCoefficient(Mesh *m, Coefficient &s, int n_)
Definition: gridfunc.hpp:847
void ComputeMeans(AvgType type, Array< int > &zones_per_vdof)
Definition: gridfunc.cpp:2311
void GetNodalValues(int i, Array< double > &nval, int vdim=1) const
Returns the values in the vertices of i&#39;th element for dimension vdim.
Definition: gridfunc.cpp:396
virtual void ComputeElementL2Errors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:605
void GetVectorFieldNodalValues(Vector &val, int comp) const
Definition: gridfunc.cpp:1300
void ImposeBounds(int i, const Vector &weights, const Vector &lo_, const Vector &hi_)
Definition: gridfunc.cpp:1907
virtual void Save(std::ostream &out) const
Save the GridFunction to an output stream.
Definition: gridfunc.cpp:3676
const FiniteElementSpace * FESpace() const
Definition: gridfunc.hpp:669
virtual void GetElementDofValues(int el, Vector &dof_vals) const
Definition: gridfunc.cpp:1844
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
Definition: gridfunc.cpp:4487
GridFunction(FiniteElementSpace *f, double *data)
Construct a GridFunction using previously allocated array data.
Definition: gridfunc.hpp:94
void GetDerivative(int comp, int der_comp, GridFunction &der)
Compute a certain derivative of a function&#39;s component. Derivatives of the function are computed at t...
Definition: gridfunc.cpp:1424
void ProjectBdrCoefficientNormal(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Definition: gridfunc.cpp:2678
void SetTrueVector()
Shortcut for calling GetTrueDofs() with GetTrueVector() as argument.
Definition: gridfunc.hpp:146
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
void GetHessians(int i, const IntegrationRule &ir, DenseMatrix &hess, int vdim=1) const
Definition: gridfunc.cpp:600
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:570
void GetGradients(const int elem, const IntegrationRule &ir, DenseMatrix &grad) const
Definition: gridfunc.hpp:337
void GetVectorFieldValues(int i, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr, int comp=0) const
Definition: gridfunc.cpp:1235
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
Definition: gridfunc.cpp:3140
GridFunction(const GridFunction &orig)
Copy constructor. The internal true-dof vector t_vec is not copied.
Definition: gridfunc.hpp:80
virtual double ComputeElementGradError(int ielem, VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Returns ||grad u_ex - grad u_h||_L2 in element ielem for H1 or L2 elements.
Definition: gridfunc.cpp:2862
JumpScaling(double nu_=1.0, JumpScalingType type_=CONSTANT)
Definition: gridfunc.hpp:759
void MakeTRef(FiniteElementSpace *f, double *tv)
Associate a new FiniteElementSpace and new true-dof data with the GridFunction.
Definition: gridfunc.cpp:222
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:578
void AccumulateAndCountDerivativeValues(int comp, int der_comp, GridFunction &der, Array< int > &zones_per_dof)
Used for the serial and parallel implementations of the GetDerivative() method; see its documentation...
Definition: gridfunc.cpp:1371
void AccumulateAndCountZones(Coefficient &coeff, AvgType type, Array< int > &zones_per_vdof)
Accumulates (depending on type) the values of coeff at all shared vdofs and counts in how many zones ...
Definition: gridfunc.cpp:2021
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
Definition: gridfunc.cpp:204
void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps, int vdim=1) const
Definition: gridfunc.cpp:561
void GetTrueDofs(Vector &tv) const
Extract the true-dofs from the GridFunction.
Definition: gridfunc.cpp:367
void GetValues(int i, const IntegrationRule &ir, Vector &vals, int vdim=1) const
Definition: gridfunc.cpp:513
FiniteElementCollection * OwnFEC()
Definition: gridfunc.hpp:126
virtual ~ExtrudeCoefficient()
Definition: gridfunc.hpp:850
FiniteElementSpace * FESpace()
Definition: gridfunc.hpp:668
GridFunction & operator=(const GridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: gridfunc.hpp:118
virtual void ComputeElementL1Errors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:599
const Vector & GetTrueVector() const
Read only access to the (optional) internal true-dof Vector.
Definition: gridfunc.hpp:132
ElementTransformation * GetElementTransformation(int i) const
Returns ElementTransformation for the i-th element.
Definition: fespace.hpp:647
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:2900
virtual void ComputeElementL2Errors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:641
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:35
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:96
virtual void ComputeElementL1Errors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:635
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:41
FiniteElementSpace * fes
FE space on which the grid function lives. Owned if fec is not NULL.
Definition: gridfunc.hpp:34
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Definition: fe_coll.hpp:26
void ProjectGridFunction(const GridFunction &src)
Project the src GridFunction to this GridFunction, both of which must be on the same mesh...
Definition: gridfunc.cpp:1855
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:555
FiniteElementCollection * fec
Used when the grid function is read from a file. It can also be set explicitly, see MakeOwner()...
Definition: gridfunc.hpp:40
void GetVectorGradient(ElementTransformation &tr, DenseMatrix &grad) const
Definition: gridfunc.cpp:1742
virtual double ComputeL2Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL, const Array< int > *elems=NULL) const
Definition: gridfunc.cpp:2766
double GetDivergence(ElementTransformation &tr) const
Definition: gridfunc.cpp:1459
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: gridfunc.cpp:165
double ComputeElementLpDistance(double p, int i, GridFunction &gf1, GridFunction &gf2)
Compute the Lp distance between two grid functions on the given element.
Definition: gridfunc.cpp:4430
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Definition: gridfunc.cpp:2981
int GetFaceVectorValues(int i, int side, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr) const
Definition: gridfunc.cpp:1106
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:3175
double ZZErrorEstimator(BilinearFormIntegrator &blfi, GridFunction &u, GridFunction &flux, Vector &error_estimates, Array< int > *aniso_flags, int with_subdomains, bool with_coeff)
Definition: gridfunc.cpp:3978
virtual void ComputeElementMaxErrors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:611
Class for integration point with weight.
Definition: intrules.hpp:25
void SaveSTL(std::ostream &out, int TimesToRefine=1)
Write the GridFunction in STL format. Note that the mesh dimension must be 2 and that quad elements w...
Definition: gridfunc.cpp:3816
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL, const Array< int > *elems=NULL) const
Definition: gridfunc.cpp:3346
double LSZZErrorEstimator(BilinearFormIntegrator &blfi, GridFunction &u, Vector &error_estimates, bool subdomain_reconstruction, bool with_coeff, double tichonov_coeff)
A ``true&#39;&#39; ZZ error estimator that uses face-based patches for flux reconstruction.
Definition: gridfunc.cpp:4227
int dim
Definition: ex24.cpp:53
int GetFaceValues(int i, int side, const IntegrationRule &ir, Vector &vals, DenseMatrix &tr, int vdim=1) const
Definition: gridfunc.cpp:656
void GetGradients(ElementTransformation &tr, const IntegrationRule &ir, DenseMatrix &grad) const
Definition: gridfunc.cpp:1713
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
Definition: gridfunc.cpp:2749
void LegacyNCReorder()
Loading helper.
Definition: gridfunc.cpp:3916
void TensorProductLegendre(int dim, int order, const Vector &x_in, const Vector &xmax, const Vector &xmin, Vector &poly, double angle, const Vector *midpoint)
Defines the global tensor product polynomial space used by NewZZErorrEstimator.
Definition: gridfunc.cpp:4068
double infinity()
Define a shortcut for std::numeric_limits&lt;double&gt;::infinity()
Definition: vector.hpp:46
GridFunction(FiniteElementSpace *f, Vector &base, int base_offset=0)
Construct a GridFunction using previously allocated Vector base starting at the given offset...
Definition: gridfunc.hpp:100
virtual void ProjectCoefficient(Coefficient &coeff)
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
Definition: gridfunc.cpp:2399
virtual void SetSpace(FiniteElementSpace *f)
Associate a new FiniteElementSpace with the GridFunction.
Definition: gridfunc.cpp:196
void ReorderByNodes()
For a vector grid function, makes sure that the ordering is byNODES.
Definition: gridfunc.cpp:1273
std::ostream & operator<<(std::ostream &os, SparseMatrix const &mat)
Definition: sparsemat.hpp:710
void GetCurl(ElementTransformation &tr, Vector &curl) const
Definition: gridfunc.cpp:1549
Vector data type.
Definition: vector.hpp:60
void MakeRef(Vector &base, int offset, int size)
Reset the Vector to be a reference to a sub-vector of base.
Definition: vector.hpp:577
virtual bool UseDevice() const
Return the device flag of the Memory object used by the Vector.
Definition: vector.hpp:121
virtual void ComputeFlux(BilinearFormIntegrator &blfi, GridFunction &flux, bool wcoef=true, int subdomain=-1)
Definition: gridfunc.cpp:309
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, const Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3240
GridFunction(FiniteElementSpace *f)
Construct a GridFunction associated with the FiniteElementSpace *f.
Definition: gridfunc.hpp:85
RefCoord s[3]
virtual double ComputeDGFaceJumpError(Coefficient *exsol, Coefficient *ell_coeff, class JumpScaling jump_scaling, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3017
double u(const Vector &xvec)
Definition: lor_mms.hpp:24
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL, const Array< int > *elems=NULL) const
Returns ||u_ex - u_h||_L2 for H1 or L2 elements.
Definition: gridfunc.hpp:493
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
long GetSequence() const
Definition: fespace.hpp:922
double Eval(double h, int p) const
Definition: gridfunc.hpp:761
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
Definition: gridfunc.cpp:382
void GetValuesFrom(const GridFunction &orig_func)
Definition: gridfunc.cpp:1149
void GetBdrValuesFrom(const GridFunction &orig_func)
Definition: gridfunc.cpp:1196
virtual void ComputeElementMaxErrors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:647
void AccumulateAndCountBdrValues(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr, Array< int > &values_counter)
Definition: gridfunc.cpp:2117
void SaveVTK(std::ostream &out, const std::string &field_name, int ref)
Write the GridFunction in VTK format. Note that Mesh::PrintVTK must be called first. The parameter ref &gt; 0 must match the one used in Mesh::PrintVTK.
Definition: gridfunc.cpp:3717
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:449
virtual double GetValue(int i, const IntegrationPoint &ip, int vdim=1) const
Definition: gridfunc.cpp:441
void GetGradient(ElementTransformation &tr, Vector &grad) const
Definition: gridfunc.cpp:1645
void ProjectVectorFieldOn(GridFunction &vec_field, int comp=0)
Definition: gridfunc.cpp:1331