MFEM  v4.3.0
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-2021, 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  /// Construct a GridFunction on the given Mesh, using the data from @a input.
99  /** The content of @a input should be in the format created by the method
100  Save(). The reconstructed FiniteElementSpace and FiniteElementCollection
101  are owned by the GridFunction. */
102  GridFunction(Mesh *m, std::istream &input);
103 
104  GridFunction(Mesh *m, GridFunction *gf_array[], int num_pieces);
105 
106  /// Copy assignment. Only the data of the base class Vector is copied.
107  /** It is assumed that this object and @a rhs use FiniteElementSpace%s that
108  have the same size.
109 
110  @note Defining this method overwrites the implicitly defined copy
111  assignment operator. */
113  { return operator=((const Vector &)rhs); }
114 
115  /// Make the GridFunction the owner of #fec and #fes.
116  /** If the new FiniteElementCollection, @a fec_, is NULL, ownership of #fec
117  and #fes is taken away. */
118  void MakeOwner(FiniteElementCollection *fec_) { fec = fec_; }
119 
121 
122  int VectorDim() const;
123 
124  /// Read only access to the (optional) internal true-dof Vector.
125  /** Note that the returned Vector may be empty, if not previously allocated
126  or set. */
127  const Vector &GetTrueVector() const { return t_vec; }
128  /// Read and write access to the (optional) internal true-dof Vector.
129  /** Note that the returned Vector may be empty, if not previously allocated
130  or set. */
131  Vector &GetTrueVector() { return t_vec; }
132 
133  /// Extract the true-dofs from the GridFunction.
134  void GetTrueDofs(Vector &tv) const;
135 
136  /// Shortcut for calling GetTrueDofs() with GetTrueVector() as argument.
138 
139  /// Set the GridFunction from the given true-dof vector.
140  virtual void SetFromTrueDofs(const Vector &tv);
141 
142  /// Shortcut for calling SetFromTrueDofs() with GetTrueVector() as argument.
144 
145  /// Returns the values in the vertices of i'th element for dimension vdim.
146  void GetNodalValues(int i, Array<double> &nval, int vdim = 1) const;
147 
148  /** @name Element index Get Value Methods
149 
150  These methods take an element index and return the interpolated value of
151  the field at a given reference point within the element.
152 
153  @warning These methods retrieve and use the ElementTransformation object
154  from the mfem::Mesh. This can alter the state of the element
155  transformation object and can also lead to unexpected results when the
156  ElementTransformation object is already in use such as when these methods
157  are called from within an integration loop. Consider using
158  GetValue(ElementTransformation &T, ...) instead.
159  */
160  ///@{
161  /** Return a scalar value from within the given element. */
162  virtual double GetValue(int i, const IntegrationPoint &ip,
163  int vdim = 1) const;
164 
165  /** Return a vector value from within the given element. */
166  virtual void GetVectorValue(int i, const IntegrationPoint &ip,
167  Vector &val) const;
168  ///@}
169 
170  /** @name Element Index Get Values Methods
171 
172  These are convenience methods for repeatedly calling GetValue for
173  multiple points within a given element. The GetValues methods are
174  optimized and should perform better than repeatedly calling GetValue. The
175  GetVectorValues method simply calls GetVectorValue repeatedly.
176 
177  @warning These methods retrieve and use the ElementTransformation object
178  from the mfem::Mesh. This can alter the state of the element
179  transformation object and can also lead to unexpected results when the
180  ElementTransformation object is already in use such as when these methods
181  are called from within an integration loop. Consider using
182  GetValues(ElementTransformation &T, ...) instead.
183  */
184  ///@{
185  /** Compute a collection of scalar values from within the element indicated
186  by the index i. */
187  void GetValues(int i, const IntegrationRule &ir, Vector &vals,
188  int vdim = 1) const;
189 
190  /** Compute a collection of vector values from within the element indicated
191  by the index i. */
192  void GetValues(int i, const IntegrationRule &ir, Vector &vals,
193  DenseMatrix &tr, int vdim = 1) const;
194 
195  void GetVectorValues(int i, const IntegrationRule &ir,
196  DenseMatrix &vals, DenseMatrix &tr) const;
197  ///@}
198 
199  /** @name ElementTransformation Get Value Methods
200 
201  These member functions are designed for use within
202  GridFunctionCoefficient objects. These can be used with
203  ElementTransformation objects coming from either
204  Mesh::GetElementTransformation() or Mesh::GetBdrElementTransformation().
205 
206  @note These methods do not reset the ElementTransformation object so they
207  should be safe to use within integration loops or other contexts where
208  the ElementTransformation is already in use.
209  */
210  ///@{
211  /** Return a scalar value from within the element indicated by the
212  ElementTransformation Object. */
213  virtual double GetValue(ElementTransformation &T, const IntegrationPoint &ip,
214  int comp = 0, Vector *tr = NULL) const;
215 
216  /** Return a vector value from within the element indicated by the
217  ElementTransformation Object. */
218  virtual void GetVectorValue(ElementTransformation &T,
219  const IntegrationPoint &ip,
220  Vector &val, Vector *tr = NULL) const;
221  ///@}
222 
223  /** @name ElementTransformation Get Values Methods
224 
225  These are convenience methods for repeatedly calling GetValue for
226  multiple points within a given element. They work by calling either the
227  ElementTransformation or FaceElementTransformations versions described
228  above. Consequently, these methods should not be expected to run faster
229  than calling the above methods in an external loop.
230 
231  @note These methods do not reset the ElementTransformation object so they
232  should be safe to use within integration loops or other contexts where
233  the ElementTransformation is already in use.
234 
235  @note These methods can also be used with FaceElementTransformations
236  objects.
237  */
238  ///@{
239  /** Compute a collection of scalar values from within the element indicated
240  by the ElementTransformation object. */
242  Vector &vals, int comp = 0, DenseMatrix *tr = NULL) const;
243 
244  /** Compute a collection of vector values from within the element indicated
245  by the ElementTransformation object. */
247  DenseMatrix &vals, DenseMatrix *tr = NULL) const;
248  ///@}
249 
250  /** @name Face Index Get Values Methods
251 
252  These methods are designed to work with Discontinuous Galerkin basis
253  functions. They compute field values on the interface between elements,
254  or on boundary elements, by interpolating the field in a neighboring
255  element. The \a side argument indices which neighboring element should be
256  used: 0, 1, or 2 (automatically chosen).
257 
258  @warning These methods retrieve and use the FaceElementTransformations
259  object from the mfem::Mesh. This can alter the state of the face element
260  transformations object and can also lead to unexpected results when the
261  FaceElementTransformations object is already in use such as when these
262  methods are called from within an integration loop. Consider using
263  GetValues(ElementTransformation &T, ...) instead.
264  */
265  ///@{
266  /** Compute a collection of scalar values from within the face
267  indicated by the index i. */
268  int GetFaceValues(int i, int side, const IntegrationRule &ir, Vector &vals,
269  DenseMatrix &tr, int vdim = 1) const;
270 
271  /** Compute a collection of vector values from within the face
272  indicated by the index i. */
273  int GetFaceVectorValues(int i, int side, const IntegrationRule &ir,
274  DenseMatrix &vals, DenseMatrix &tr) const;
275  ///@}
276 
277  void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps,
278  int vdim = 1) const;
279 
280  void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps,
281  DenseMatrix &tr, int vdim = 1) const;
282 
283  void GetHessians(int i, const IntegrationRule &ir, DenseMatrix &hess,
284  int vdim = 1) const;
285 
286  void GetHessians(int i, const IntegrationRule &ir, DenseMatrix &hess,
287  DenseMatrix &tr, int vdim = 1) const;
288 
289  void GetValuesFrom(const GridFunction &orig_func);
290 
291  void GetBdrValuesFrom(const GridFunction &orig_func);
292 
293  void GetVectorFieldValues(int i, const IntegrationRule &ir,
294  DenseMatrix &vals,
295  DenseMatrix &tr, int comp = 0) const;
296 
297  /// For a vector grid function, makes sure that the ordering is byNODES.
298  void ReorderByNodes();
299 
300  /// Return the values as a vector on mesh vertices for dimension vdim.
301  void GetNodalValues(Vector &nval, int vdim = 1) const;
302 
303  void GetVectorFieldNodalValues(Vector &val, int comp) const;
304 
305  void ProjectVectorFieldOn(GridFunction &vec_field, int comp = 0);
306 
307  void GetDerivative(int comp, int der_comp, GridFunction &der);
308 
309  double GetDivergence(ElementTransformation &tr) const;
310 
311  void GetCurl(ElementTransformation &tr, Vector &curl) const;
312 
313  void GetGradient(ElementTransformation &tr, Vector &grad) const;
314 
316  DenseMatrix &grad) const;
317 
318  void GetGradients(const int elem, const IntegrationRule &ir,
319  DenseMatrix &grad) const
320  { GetGradients(*fes->GetElementTransformation(elem), ir, grad); }
321 
322  void GetVectorGradient(ElementTransformation &tr, DenseMatrix &grad) const;
323 
324  /** Compute \f$ (\int_{\Omega} (*this) \psi_i)/(\int_{\Omega} \psi_i) \f$,
325  where \f$ \psi_i \f$ are the basis functions for the FE space of avgs.
326  Both FE spaces should be scalar and on the same mesh. */
327  void GetElementAverages(GridFunction &avgs) const;
328 
329  /** Sets the output vector @a dof_vals to the values of the degrees of
330  freedom of element @a el. */
331  virtual void GetElementDofValues(int el, Vector &dof_vals) const;
332 
333  /** Impose the given bounds on the function's DOFs while preserving its local
334  * integral (described in terms of the given weights) on the i'th element
335  * through SLBPQ optimization.
336  * Intended to be used for discontinuous FE functions. */
337  void ImposeBounds(int i, const Vector &weights,
338  const Vector &lo_, const Vector &hi_);
339  void ImposeBounds(int i, const Vector &weights,
340  double min_ = 0.0, double max_ = infinity());
341 
342  /** On a non-conforming mesh, make sure the function lies in the conforming
343  space by multiplying with R and then with P, the conforming restriction
344  and prolongation matrices of the space, respectively. */
345  void RestrictConforming();
346 
347  /** @brief Project the @a src GridFunction to @a this GridFunction, both of
348  which must be on the same mesh. */
349  /** The current implementation assumes that all elements use the same
350  projection matrix. */
351  void ProjectGridFunction(const GridFunction &src);
352 
353  /** @brief Project @a coeff Coefficient to @a this GridFunction. The
354  projection computation depends on the choice of the FiniteElementSpace
355  #fes. Note that this is usually interpolation at the degrees of freedom
356  in each element (not L2 projection). */
357  virtual void ProjectCoefficient(Coefficient &coeff);
358 
359  /** @brief Project @a coeff Coefficient to @a this GridFunction, using one
360  element for each degree of freedom in @a dofs and nodal interpolation on
361  that element. */
362  void ProjectCoefficient(Coefficient &coeff, Array<int> &dofs, int vd = 0);
363 
364  /** @brief Project @a vcoeff VectorCoefficient to @a this GridFunction. The
365  projection computation depends on the choice of the FiniteElementSpace
366  #fes. Note that this is usually interpolation at the degrees of freedom
367  in each element (not L2 projection).*/
369 
370  /** @brief Project @a vcoeff VectorCoefficient to @a this GridFunction, using
371  one element for each degree of freedom in @a dofs and nodal interpolation
372  on that element. */
373  void ProjectCoefficient(VectorCoefficient &vcoeff, Array<int> &dofs);
374 
375  /** @brief Project @a vcoeff VectorCoefficient to @a this GridFunction, only
376  projecting onto elements with the given @a attribute */
377  void ProjectCoefficient(VectorCoefficient &vcoeff, int attribute);
378 
379  /** @brief Analogous to the version with argument @a vcoeff VectorCoefficient
380  but using an array of scalar coefficients for each component. */
381  void ProjectCoefficient(Coefficient *coeff[]);
382 
383  /** @brief Project a discontinuous vector coefficient as a grid function on
384  a continuous finite element space. The values in shared dofs are
385  determined from the element with maximal attribute. */
386  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff);
387 
389  /** @brief Projects a discontinuous coefficient so that the values in shared
390  vdofs are computed by taking an average of the possible values. */
391  virtual void ProjectDiscCoefficient(Coefficient &coeff, AvgType type);
392  /** @brief Projects a discontinuous _vector_ coefficient so that the values
393  in shared vdofs are computed by taking an average of the possible values.
394  */
395  virtual void ProjectDiscCoefficient(VectorCoefficient &coeff, AvgType type);
396 
397 protected:
398  /** @brief Accumulates (depending on @a type) the values of @a coeff at all
399  shared vdofs and counts in how many zones each vdof appears. */
400  void AccumulateAndCountZones(Coefficient &coeff, AvgType type,
401  Array<int> &zones_per_vdof);
402 
403  /** @brief Accumulates (depending on @a type) the values of @a vcoeff at all
404  shared vdofs and counts in how many zones each vdof appears. */
406  Array<int> &zones_per_vdof);
407 
409  VectorCoefficient *vcoeff, Array<int> &attr,
410  Array<int> &values_counter);
411 
413  Array<int> &bdr_attr,
414  Array<int> &values_counter);
415 
416  // Complete the computation of averages; called e.g. after
417  // AccumulateAndCountZones().
418  void ComputeMeans(AvgType type, Array<int> &zones_per_vdof);
419 
420 public:
421  /** @brief Project a Coefficient on the GridFunction, modifying only DOFs on
422  the boundary associated with the boundary attributes marked in the
423  @a attr array. */
425  {
426  Coefficient *coeff_p = &coeff;
427  ProjectBdrCoefficient(&coeff_p, attr);
428  }
429 
430  /** @brief Project a VectorCoefficient on the GridFunction, modifying only
431  DOFs on the boundary associated with the boundary attributes marked in
432  the @a attr array. */
433  virtual void ProjectBdrCoefficient(VectorCoefficient &vcoeff,
434  Array<int> &attr);
435 
436  /** @brief Project a set of Coefficient%s on the components of the
437  GridFunction, modifying only DOFs on the boundary associated with the
438  boundary attributed marked in the @a attr array. */
439  /** If a Coefficient pointer in the array @a coeff is NULL, that component
440  will not be touched. */
441  virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array<int> &attr);
442 
443  /** Project the normal component of the given VectorCoefficient on
444  the boundary. Only boundary attributes that are marked in
445  'bdr_attr' are projected. Assumes RT-type VectorFE GridFunction. */
447  Array<int> &bdr_attr);
448 
449  /** @brief Project the tangential components of the given VectorCoefficient
450  on the boundary. Only boundary attributes that are marked in @a bdr_attr
451  are projected. Assumes ND-type VectorFE GridFunction. */
452  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff,
453  Array<int> &bdr_attr);
454 
455 
456  virtual double ComputeL2Error(Coefficient &exsol,
457  const IntegrationRule *irs[] = NULL) const
458  { return ComputeLpError(2.0, exsol, NULL, irs); }
459 
460  virtual double ComputeL2Error(Coefficient *exsol[],
461  const IntegrationRule *irs[] = NULL) const;
462 
463  virtual double ComputeL2Error(VectorCoefficient &exsol,
464  const IntegrationRule *irs[] = NULL,
465  Array<int> *elems = NULL) const;
466 
467  /// Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements
468  virtual double ComputeGradError(VectorCoefficient *exgrad,
469  const IntegrationRule *irs[] = NULL) const;
470 
471  /// Returns ||curl u_ex - curl u_h||_L2 for ND elements
472  virtual double ComputeCurlError(VectorCoefficient *excurl,
473  const IntegrationRule *irs[] = NULL) const;
474 
475  /// Returns ||div u_ex - div u_h||_L2 for RT elements
476  virtual double ComputeDivError(Coefficient *exdiv,
477  const IntegrationRule *irs[] = NULL) const;
478 
479  /// Returns the Face Jumps error for L2 elements. The error can be weighted
480  /// by a constant nu, by nu/h, or nu*p^2/h, depending on the value of
481  /// @a jump_scaling.
482  virtual double ComputeDGFaceJumpError(Coefficient *exsol,
483  Coefficient *ell_coeff,
484  class JumpScaling jump_scaling,
485  const IntegrationRule *irs[] = NULL)
486  const;
487 
488  /// Returns the Face Jumps error for L2 elements, with 1/h scaling.
489  MFEM_DEPRECATED
490  double ComputeDGFaceJumpError(Coefficient *exsol,
491  Coefficient *ell_coeff,
492  double Nu,
493  const IntegrationRule *irs[] = NULL) const;
494 
495  /** This method is kept for backward compatibility.
496 
497  Returns either the H1-seminorm, or the DG face jumps error, or both
498  depending on norm_type = 1, 2, 3. Additional arguments for the DG face
499  jumps norm: ell_coeff: mesh-depended coefficient (weight) Nu: scalar
500  constant weight */
501  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
502  Coefficient *ell_coef, double Nu,
503  int norm_type) const;
504 
505  /// Returns the error measured in H1-norm for H1 elements or in "broken"
506  /// H1-norm for L2 elements
507  virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad,
508  const IntegrationRule *irs[] = NULL) const;
509 
510  /// Returns the error measured in H(div)-norm for RT elements
511  virtual double ComputeHDivError(VectorCoefficient *exsol,
512  Coefficient *exdiv,
513  const IntegrationRule *irs[] = NULL) const;
514 
515  /// Returns the error measured in H(curl)-norm for ND elements
516  virtual double ComputeHCurlError(VectorCoefficient *exsol,
517  VectorCoefficient *excurl,
518  const IntegrationRule *irs[] = NULL) const;
519 
520  virtual double ComputeMaxError(Coefficient &exsol,
521  const IntegrationRule *irs[] = NULL) const
522  {
523  return ComputeLpError(infinity(), exsol, NULL, irs);
524  }
525 
526  virtual double ComputeMaxError(Coefficient *exsol[],
527  const IntegrationRule *irs[] = NULL) const;
528 
529  virtual double ComputeMaxError(VectorCoefficient &exsol,
530  const IntegrationRule *irs[] = NULL) const
531  {
532  return ComputeLpError(infinity(), exsol, NULL, NULL, irs);
533  }
534 
535  virtual double ComputeL1Error(Coefficient &exsol,
536  const IntegrationRule *irs[] = NULL) const
537  { return ComputeLpError(1.0, exsol, NULL, irs); }
538 
539  virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad,
540  int norm_type, Array<int> *elems = NULL,
541  const IntegrationRule *irs[] = NULL) const;
542 
543  virtual double ComputeL1Error(VectorCoefficient &exsol,
544  const IntegrationRule *irs[] = NULL) const
545  { return ComputeLpError(1.0, exsol, NULL, NULL, irs); }
546 
547  virtual double ComputeLpError(const double p, Coefficient &exsol,
548  Coefficient *weight = NULL,
549  const IntegrationRule *irs[] = NULL) const;
550 
551  /** Compute the Lp error in each element of the mesh and store the results in
552  the Vector @a error. The result should be of length number of elements,
553  for example an L2 GridFunction of order zero using map type VALUE. */
554  virtual void ComputeElementLpErrors(const double p, Coefficient &exsol,
555  Vector &error,
556  Coefficient *weight = NULL,
557  const IntegrationRule *irs[] = NULL
558  ) const;
559 
560  virtual void ComputeElementL1Errors(Coefficient &exsol,
561  Vector &error,
562  const IntegrationRule *irs[] = NULL
563  ) const
564  { ComputeElementLpErrors(1.0, exsol, error, NULL, irs); }
565 
566  virtual void ComputeElementL2Errors(Coefficient &exsol,
567  Vector &error,
568  const IntegrationRule *irs[] = NULL
569  ) const
570  { ComputeElementLpErrors(2.0, exsol, error, NULL, irs); }
571 
573  Vector &error,
574  const IntegrationRule *irs[] = NULL
575  ) const
576  { ComputeElementLpErrors(infinity(), exsol, error, NULL, irs); }
577 
578  /** When given a vector weight, compute the pointwise (scalar) error as the
579  dot product of the vector error with the vector weight. Otherwise, the
580  scalar error is the l_2 norm of the vector error. */
581  virtual double ComputeLpError(const double p, VectorCoefficient &exsol,
582  Coefficient *weight = NULL,
583  VectorCoefficient *v_weight = NULL,
584  const IntegrationRule *irs[] = NULL) const;
585 
586  /** Compute the Lp error in each element of the mesh and store the results in
587  the Vector @ error. The result should be of length number of elements,
588  for example an L2 GridFunction of order zero using map type VALUE. */
589  virtual void ComputeElementLpErrors(const double p, VectorCoefficient &exsol,
590  Vector &error,
591  Coefficient *weight = NULL,
592  VectorCoefficient *v_weight = NULL,
593  const IntegrationRule *irs[] = NULL
594  ) const;
595 
597  Vector &error,
598  const IntegrationRule *irs[] = NULL
599  ) const
600  { ComputeElementLpErrors(1.0, exsol, error, NULL, NULL, irs); }
601 
603  Vector &error,
604  const IntegrationRule *irs[] = NULL
605  ) const
606  { ComputeElementLpErrors(2.0, exsol, error, NULL, NULL, irs); }
607 
609  Vector &error,
610  const IntegrationRule *irs[] = NULL
611  ) const
612  { ComputeElementLpErrors(infinity(), exsol, error, NULL, NULL, irs); }
613 
614  virtual void ComputeFlux(BilinearFormIntegrator &blfi,
615  GridFunction &flux,
616  bool wcoef = true, int subdomain = -1);
617 
618  /// Redefine '=' for GridFunction = constant.
619  GridFunction &operator=(double value);
620 
621  /// Copy the data from @a v.
622  /** The size of @a v must be equal to the size of the associated
623  FiniteElementSpace #fes. */
624  GridFunction &operator=(const Vector &v);
625 
626  /// Transform by the Space UpdateMatrix (e.g., on Mesh change).
627  virtual void Update();
628 
630  const FiniteElementSpace *FESpace() const { return fes; }
631 
632  /// Associate a new FiniteElementSpace with the GridFunction.
633  /** The GridFunction is resized using the SetSize() method. */
634  virtual void SetSpace(FiniteElementSpace *f);
635 
636  using Vector::MakeRef;
637 
638  /** @brief Make the GridFunction reference external data on a new
639  FiniteElementSpace. */
640  /** This method changes the FiniteElementSpace associated with the
641  GridFunction and sets the pointer @a v as external data in the
642  GridFunction. */
643  virtual void MakeRef(FiniteElementSpace *f, double *v);
644 
645  /** @brief Make the GridFunction reference external data on a new
646  FiniteElementSpace. */
647  /** This method changes the FiniteElementSpace associated with the
648  GridFunction and sets the data of the Vector @a v (plus the @a v_offset)
649  as external data in the GridFunction.
650  @note This version of the method will also perform bounds checks when
651  the build option MFEM_DEBUG is enabled. */
652  virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset);
653 
654  /** @brief Associate a new FiniteElementSpace and new true-dof data with the
655  GridFunction. */
656  /** - If the prolongation matrix of @a f is trivial (i.e. its method
657  FiniteElementSpace::GetProlongationMatrix() returns NULL), then the
658  method MakeRef() is called with the same arguments.
659  - Otherwise, the method SetSpace() is called with argument @a f.
660  - The internal true-dof vector is set to reference @a tv. */
661  void MakeTRef(FiniteElementSpace *f, double *tv);
662 
663  /** @brief Associate a new FiniteElementSpace and new true-dof data with the
664  GridFunction. */
665  /** - If the prolongation matrix of @a f is trivial (i.e. its method
666  FiniteElementSpace::GetProlongationMatrix() returns NULL), this method
667  calls MakeRef() with the same arguments.
668  - Otherwise, this method calls SetSpace() with argument @a f.
669  - The internal true-dof vector is set to reference the sub-vector of
670  @a tv starting at the offset @a tv_offset. */
671  void MakeTRef(FiniteElementSpace *f, Vector &tv, int tv_offset);
672 
673  /// Save the GridFunction to an output stream.
674  virtual void Save(std::ostream &out) const;
675 
676  /// Save the GridFunction to a file. The given @a precision will be used for
677  /// ASCII output.
678  virtual void Save(const char *fname, int precision=16) const;
679 
680 #ifdef MFEM_USE_ADIOS2
681  /// Save the GridFunction to a binary output stream using adios2 bp format.
682  virtual void Save(adios2stream &out, const std::string& variable_name,
685 #endif
686 
687  /** @brief Write the GridFunction in VTK format. Note that Mesh::PrintVTK
688  must be called first. The parameter ref > 0 must match the one used in
689  Mesh::PrintVTK. */
690  void SaveVTK(std::ostream &out, const std::string &field_name, int ref);
691 
692  /** @brief Write the GridFunction in STL format. Note that the mesh dimension
693  must be 2 and that quad elements will be broken into two triangles.*/
694  void SaveSTL(std::ostream &out, int TimesToRefine = 1);
695 
696  /// Destroys grid function.
697  virtual ~GridFunction() { Destroy(); }
698 };
699 
700 
701 /** Overload operator<< for std::ostream and GridFunction; valid also for the
702  derived class ParGridFunction */
703 std::ostream &operator<<(std::ostream &out, const GridFunction &sol);
704 
705 /// Class used to specify how the jump terms in
706 /// GridFunction::ComputeDGFaceJumpError are scaled.
708 {
709 public:
711  {
715  };
716 private:
717  double nu;
718  JumpScalingType type;
719 public:
720  JumpScaling(double nu_=1.0, JumpScalingType type_=CONSTANT)
721  : nu(nu_), type(type_) { }
722  double Eval(double h, int p) const
723  {
724  double val = nu;
725  if (type != CONSTANT) { val /= h; }
726  if (type == P_SQUARED_OVER_H) { val *= p*p; }
727  return val;
728  }
729 };
730 
731 
732 /** @brief Class representing a function through its values (scalar or vector)
733  at quadrature points. */
735 {
736 protected:
737  QuadratureSpace *qspace; ///< Associated QuadratureSpace
738  int vdim; ///< Vector dimension
739  bool own_qspace; ///< QuadratureSpace ownership flag
740 
741 public:
742  /// Create an empty QuadratureFunction.
743  /** The object can be initialized later using the SetSpace() methods. */
745  : qspace(NULL), vdim(0), own_qspace(false) { }
746 
747  /** @brief Copy constructor. The QuadratureSpace ownership flag, #own_qspace,
748  in the new object is set to false. */
750  : Vector(orig),
751  qspace(orig.qspace), vdim(orig.vdim), own_qspace(false) { }
752 
753  /// Create a QuadratureFunction based on the given QuadratureSpace.
754  /** The QuadratureFunction does not assume ownership of the QuadratureSpace.
755  @note The Vector data is not initialized. */
756  QuadratureFunction(QuadratureSpace *qspace_, int vdim_ = 1)
757  : Vector(vdim_*qspace_->GetSize()),
758  qspace(qspace_), vdim(vdim_), own_qspace(false) { }
759 
760  /** @brief Create a QuadratureFunction based on the given QuadratureSpace,
761  using the external data, @a qf_data. */
762  /** The QuadratureFunction does not assume ownership of neither the
763  QuadratureSpace nor the external data. */
764  QuadratureFunction(QuadratureSpace *qspace_, double *qf_data, int vdim_ = 1)
765  : Vector(qf_data, vdim_*qspace_->GetSize()),
766  qspace(qspace_), vdim(vdim_), own_qspace(false) { }
767 
768  /// Read a QuadratureFunction from the stream @a in.
769  /** The QuadratureFunction assumes ownership of the read QuadratureSpace. */
770  QuadratureFunction(Mesh *mesh, std::istream &in);
771 
772  virtual ~QuadratureFunction() { if (own_qspace) { delete qspace; } }
773 
774  /// Get the associated QuadratureSpace.
775  QuadratureSpace *GetSpace() const { return qspace; }
776 
777  /// Change the QuadratureSpace and optionally the vector dimension.
778  /** If the new QuadratureSpace is different from the current one, the
779  QuadratureFunction will not assume ownership of the new space; otherwise,
780  the ownership flag remains the same.
781 
782  If the new vector dimension @a vdim_ < 0, the vector dimension remains
783  the same.
784 
785  The data size is updated by calling Vector::SetSize(). */
786  inline void SetSpace(QuadratureSpace *qspace_, int vdim_ = -1);
787 
788  /** @brief Change the QuadratureSpace, the data array, and optionally the
789  vector dimension. */
790  /** If the new QuadratureSpace is different from the current one, the
791  QuadratureFunction will not assume ownership of the new space; otherwise,
792  the ownership flag remains the same.
793 
794  If the new vector dimension @a vdim_ < 0, the vector dimension remains
795  the same.
796 
797  The data array is replaced by calling Vector::NewDataAndSize(). */
798  inline void SetSpace(QuadratureSpace *qspace_, double *qf_data,
799  int vdim_ = -1);
800 
801  /// Get the vector dimension.
802  int GetVDim() const { return vdim; }
803 
804  /// Set the vector dimension, updating the size by calling Vector::SetSize().
805  void SetVDim(int vdim_)
806  { vdim = vdim_; SetSize(vdim*qspace->GetSize()); }
807 
808  /// Get the QuadratureSpace ownership flag.
809  bool OwnsSpace() { return own_qspace; }
810 
811  /// Set the QuadratureSpace ownership flag.
812  void SetOwnsSpace(bool own) { own_qspace = own; }
813 
814  /// Redefine '=' for QuadratureFunction = constant.
815  QuadratureFunction &operator=(double value);
816 
817  /// Copy the data from @a v.
818  /** The size of @a v must be equal to the size of the associated
819  QuadratureSpace #qspace times the QuadratureFunction dimension
820  i.e. QuadratureFunction::Size(). */
822 
823  /// Copy assignment. Only the data of the base class Vector is copied.
824  /** The QuadratureFunctions @a v and @a *this must have QuadratureSpaces with
825  the same size.
826 
827  @note Defining this method overwrites the implicitly defined copy
828  assignment operator. */
830 
831  /// Get the IntegrationRule associated with mesh element @a idx.
832  const IntegrationRule &GetElementIntRule(int idx) const
833  { return qspace->GetElementIntRule(idx); }
834 
835  /// Return all values associated with mesh element @a idx in a Vector.
836  /** The result is stored in the Vector @a values as a reference to the
837  global values.
838 
839  Inside the Vector @a values, the index `i+vdim*j` corresponds to the
840  `i`-th vector component at the `j`-th quadrature point.
841  */
842  inline void GetElementValues(int idx, Vector &values);
843 
844  /// Return all values associated with mesh element @a idx in a Vector.
845  /** The result is stored in the Vector @a values as a copy of the
846  global values.
847 
848  Inside the Vector @a values, the index `i+vdim*j` corresponds to the
849  `i`-th vector component at the `j`-th quadrature point.
850  */
851  inline void GetElementValues(int idx, Vector &values) const;
852 
853  /// Return the quadrature function values at an integration point.
854  /** The result is stored in the Vector @a values as a reference to the
855  global values. */
856  inline void GetElementValues(int idx, const int ip_num, Vector &values);
857 
858  /// Return the quadrature function values at an integration point.
859  /** The result is stored in the Vector @a values as a copy to the
860  global values. */
861  inline void GetElementValues(int idx, const int ip_num, Vector &values) const;
862 
863  /// Return all values associated with mesh element @a idx in a DenseMatrix.
864  /** The result is stored in the DenseMatrix @a values as a reference to the
865  global values.
866 
867  Inside the DenseMatrix @a values, the `(i,j)` entry corresponds to the
868  `i`-th vector component at the `j`-th quadrature point.
869  */
870  inline void GetElementValues(int idx, DenseMatrix &values);
871 
872  /// Return all values associated with mesh element @a idx in a const DenseMatrix.
873  /** The result is stored in the DenseMatrix @a values as a copy of the
874  global values.
875 
876  Inside the DenseMatrix @a values, the `(i,j)` entry corresponds to the
877  `i`-th vector component at the `j`-th quadrature point.
878  */
879  inline void GetElementValues(int idx, DenseMatrix &values) const;
880 
881  /// Write the QuadratureFunction to the stream @a out.
882  void Save(std::ostream &out) const;
883 };
884 
885 /// Overload operator<< for std::ostream and QuadratureFunction.
886 std::ostream &operator<<(std::ostream &out, const QuadratureFunction &qf);
887 
888 
889 double ZZErrorEstimator(BilinearFormIntegrator &blfi,
890  GridFunction &u,
891  GridFunction &flux,
892  Vector &error_estimates,
893  Array<int> *aniso_flags = NULL,
894  int with_subdomains = 1,
895  bool with_coeff = false);
896 
897 /// Compute the Lp distance between two grid functions on the given element.
898 double ComputeElementLpDistance(double p, int i,
899  GridFunction& gf1, GridFunction& gf2);
900 
901 
902 /// Class used for extruding scalar GridFunctions
904 {
905 private:
906  int n;
907  Mesh *mesh_in;
908  Coefficient &sol_in;
909 public:
911  : n(n_), mesh_in(m), sol_in(s) { }
912  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip);
913  virtual ~ExtrudeCoefficient() { }
914 };
915 
916 /// Extrude a scalar 1D GridFunction, after extruding the mesh with Extrude1D.
917 GridFunction *Extrude1DGridFunction(Mesh *mesh, Mesh *mesh2d,
918  GridFunction *sol, const int ny);
919 
920 
921 // Inline methods
922 
923 inline void QuadratureFunction::SetSpace(QuadratureSpace *qspace_, int vdim_)
924 {
925  if (qspace_ != qspace)
926  {
927  if (own_qspace) { delete qspace; }
928  qspace = qspace_;
929  own_qspace = false;
930  }
931  vdim = (vdim_ < 0) ? vdim : vdim_;
933 }
934 
936  double *qf_data, int vdim_)
937 {
938  if (qspace_ != qspace)
939  {
940  if (own_qspace) { delete qspace; }
941  qspace = qspace_;
942  own_qspace = false;
943  }
944  vdim = (vdim_ < 0) ? vdim : vdim_;
945  NewDataAndSize(qf_data, vdim*qspace->GetSize());
946 }
947 
948 inline void QuadratureFunction::GetElementValues(int idx, Vector &values)
949 {
950  const int s_offset = qspace->element_offsets[idx];
951  const int sl_size = qspace->element_offsets[idx+1] - s_offset;
952  values.NewDataAndSize(data + vdim*s_offset, vdim*sl_size);
953 }
954 
955 inline void QuadratureFunction::GetElementValues(int idx, Vector &values) const
956 {
957  const int s_offset = qspace->element_offsets[idx];
958  const int sl_size = qspace->element_offsets[idx+1] - s_offset;
959  values.SetSize(vdim*sl_size);
960  const double *q = data + vdim*s_offset;
961  for (int i = 0; i<values.Size(); i++)
962  {
963  values(i) = *(q++);
964  }
965 }
966 
967 inline void QuadratureFunction::GetElementValues(int idx, const int ip_num,
968  Vector &values)
969 {
970  const int s_offset = qspace->element_offsets[idx] * vdim + ip_num * vdim;
971  values.NewDataAndSize(data + s_offset, vdim);
972 }
973 
974 inline void QuadratureFunction::GetElementValues(int idx, const int ip_num,
975  Vector &values) const
976 {
977  const int s_offset = qspace->element_offsets[idx] * vdim + ip_num * vdim;
978  values.SetSize(vdim);
979  const double *q = data + s_offset;
980  for (int i = 0; i < values.Size(); i++)
981  {
982  values(i) = *(q++);
983  }
984 }
985 
987 {
988  const int s_offset = qspace->element_offsets[idx];
989  const int sl_size = qspace->element_offsets[idx+1] - s_offset;
990  values.Reset(data + vdim*s_offset, vdim, sl_size);
991 }
992 
994  DenseMatrix &values) const
995 {
996  const int s_offset = qspace->element_offsets[idx];
997  const int sl_size = qspace->element_offsets[idx+1] - s_offset;
998  values.SetSize(vdim, sl_size);
999  const double *q = data + vdim*s_offset;
1000  for (int j = 0; j<sl_size; j++)
1001  {
1002  for (int i = 0; i<vdim; i++)
1003  {
1004  values(i,j) = *(q++);
1005  }
1006  }
1007 }
1008 
1009 } // namespace mfem
1010 
1011 #endif
void SaveSTLTri(std::ostream &out, double p1[], double p2[], double p3[])
Definition: gridfunc.cpp:3604
int GetVDim() const
Get the vector dimension.
Definition: gridfunc.hpp:802
virtual ~GridFunction()
Destroys grid function.
Definition: gridfunc.hpp:697
QuadratureFunction(QuadratureSpace *qspace_, int vdim_=1)
Create a QuadratureFunction based on the given QuadratureSpace.
Definition: gridfunc.hpp:756
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
void GetElementAverages(GridFunction &avgs) const
Definition: gridfunc.cpp:1723
void GetVectorValues(int i, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr) const
Definition: gridfunc.cpp:655
void NewDataAndSize(double *d, int s)
Set the Vector data and size, deleting the old data, if owned.
Definition: vector.hpp:153
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
const IntegrationRule & GetElementIntRule(int idx) const
Get the IntegrationRule associated with mesh element idx.
Definition: gridfunc.hpp:832
Memory< double > data
Definition: vector.hpp:64
Class used for extruding scalar GridFunctions.
Definition: gridfunc.hpp:903
void SumFluxAndCount(BilinearFormIntegrator &blfi, GridFunction &flux, Array< int > &counts, bool wcoef, int subdomain)
Definition: gridfunc.cpp:249
void SetFromTrueVector()
Shortcut for calling SetFromTrueDofs() with GetTrueVector() as argument.
Definition: gridfunc.hpp:143
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:456
void ProjectDeltaCoefficient(DeltaCoefficient &delta_coeff, double &integral)
Definition: gridfunc.cpp:2215
const IntegrationRule & GetElementIntRule(int idx) const
Get the IntegrationRule associated with mesh element idx.
Definition: fespace.hpp:950
virtual void ComputeElementLpErrors(const double p, Coefficient &exsol, Vector &error, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3221
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
Definition: gridfunc.cpp:2749
virtual void GetVectorValue(int i, const IntegrationPoint &ip, Vector &val) const
Definition: gridfunc.cpp:425
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:513
void RestrictConforming()
Definition: gridfunc.cpp:1862
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:3993
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
Definition: gridfunc.cpp:2433
bool own_qspace
QuadratureSpace ownership flag.
Definition: gridfunc.hpp:739
int GetSize() const
Return the total number of quadrature points.
Definition: fespace.hpp:938
int VectorDim() const
Definition: gridfunc.cpp:311
virtual ~QuadratureFunction()
Definition: gridfunc.hpp:772
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:190
void AccumulateAndCountBdrTangentValues(VectorCoefficient &vcoeff, Array< int > &bdr_attr, Array< int > &values_counter)
Definition: gridfunc.cpp:2139
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:2976
int vdim
Vector dimension.
Definition: gridfunc.hpp:738
void GetVectorGradientHat(ElementTransformation &T, DenseMatrix &gh) const
Definition: gridfunc.cpp:1361
Vector & GetTrueVector()
Read and write access to the (optional) internal true-dof Vector.
Definition: gridfunc.hpp:131
void MakeOwner(FiniteElementCollection *fec_)
Make the GridFunction the owner of fec and fes.
Definition: gridfunc.hpp:118
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:529
ExtrudeCoefficient(Mesh *m, Coefficient &s, int n_)
Definition: gridfunc.hpp:910
void SetSpace(QuadratureSpace *qspace_, int vdim_=-1)
Change the QuadratureSpace and optionally the vector dimension.
Definition: gridfunc.hpp:923
void ComputeMeans(AvgType type, Array< int > &zones_per_vdof)
Definition: gridfunc.cpp:2190
bool OwnsSpace()
Get the QuadratureSpace ownership flag.
Definition: gridfunc.hpp:809
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:361
virtual void ComputeElementL2Errors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:566
QuadratureFunction(QuadratureSpace *qspace_, double *qf_data, int vdim_=1)
Create a QuadratureFunction based on the given QuadratureSpace, using the external data...
Definition: gridfunc.hpp:764
void GetVectorFieldNodalValues(Vector &val, int comp) const
Definition: gridfunc.cpp:1230
void ImposeBounds(int i, const Vector &weights, const Vector &lo_, const Vector &hi_)
Definition: gridfunc.cpp:1804
virtual void Save(std::ostream &out) const
Save the GridFunction to an output stream.
Definition: gridfunc.cpp:3484
const FiniteElementSpace * FESpace() const
Definition: gridfunc.hpp:630
virtual void GetElementDofValues(int el, Vector &dof_vals) const
Definition: gridfunc.cpp:1753
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
Definition: gridfunc.cpp:3983
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)
Definition: gridfunc.cpp:1301
void ProjectBdrCoefficientNormal(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Definition: gridfunc.cpp:2527
void SetTrueVector()
Shortcut for calling GetTrueDofs() with GetTrueVector() as argument.
Definition: gridfunc.hpp:137
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:549
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:535
void GetGradients(const int elem, const IntegrationRule &ir, DenseMatrix &grad) const
Definition: gridfunc.hpp:318
void Reset(double *d, int h, int w)
Change the data array and the size of the DenseMatrix.
Definition: densemat.hpp:90
void GetVectorFieldValues(int i, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr, int comp=0) const
Definition: gridfunc.cpp:1161
void SetOwnsSpace(bool own)
Set the QuadratureSpace ownership flag.
Definition: gridfunc.hpp:812
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
Definition: gridfunc.cpp:2950
GridFunction(const GridFunction &orig)
Copy constructor. The internal true-dof vector t_vec is not copied.
Definition: gridfunc.hpp:80
JumpScaling(double nu_=1.0, JumpScalingType type_=CONSTANT)
Definition: gridfunc.hpp:720
void MakeTRef(FiniteElementSpace *f, double *tv)
Associate a new FiniteElementSpace and new true-dof data with the GridFunction.
Definition: gridfunc.cpp:219
QuadratureSpace * GetSpace() const
Get the associated QuadratureSpace.
Definition: gridfunc.hpp:775
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:543
QuadratureFunction()
Create an empty QuadratureFunction.
Definition: gridfunc.hpp:744
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:1901
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
Definition: gridfunc.cpp:201
void SetVDim(int vdim_)
Set the vector dimension, updating the size by calling Vector::SetSize().
Definition: gridfunc.hpp:805
void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps, int vdim=1) const
Definition: gridfunc.cpp:510
void GetTrueDofs(Vector &tv) const
Extract the true-dofs from the GridFunction.
Definition: gridfunc.cpp:332
void GetValues(int i, const IntegrationRule &ir, Vector &vals, int vdim=1) const
Definition: gridfunc.cpp:466
FiniteElementCollection * OwnFEC()
Definition: gridfunc.hpp:120
virtual ~ExtrudeCoefficient()
Definition: gridfunc.hpp:913
FiniteElementSpace * FESpace()
Definition: gridfunc.hpp:629
GridFunction & operator=(const GridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
Definition: gridfunc.hpp:112
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3156
virtual void ComputeElementL1Errors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:560
QuadratureFunction & operator=(double value)
Redefine &#39;=&#39; for QuadratureFunction = constant.
Definition: gridfunc.cpp:3801
const Vector & GetTrueVector() const
Read only access to the (optional) internal true-dof Vector.
Definition: gridfunc.hpp:127
ElementTransformation * GetElementTransformation(int i) const
Returns ElementTransformation for the i-th element.
Definition: fespace.hpp:600
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:2709
virtual void ComputeElementL2Errors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:602
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:34
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:87
virtual void ComputeElementL1Errors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:596
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
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:1760
QuadratureSpace * qspace
Associated QuadratureSpace.
Definition: gridfunc.hpp:737
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:520
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:1661
double GetDivergence(ElementTransformation &tr) const
Definition: gridfunc.cpp:1380
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: gridfunc.cpp:164
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:3926
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Definition: gridfunc.cpp:2791
int GetFaceVectorValues(int i, int side, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr) const
Definition: gridfunc.cpp:1044
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:2985
double ZZErrorEstimator(BilinearFormIntegrator &blfi, GridFunction &u, GridFunction &flux, Vector &error_estimates, Array< int > *aniso_flags, int with_subdomains, bool with_coeff)
Definition: gridfunc.cpp:3835
virtual void ComputeElementMaxErrors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:572
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:3624
int GetFaceValues(int i, int side, const IntegrationRule &ir, Vector &vals, DenseMatrix &tr, int vdim=1) const
Definition: gridfunc.cpp:605
void GetElementValues(int idx, Vector &values)
Return all values associated with mesh element idx in a Vector.
Definition: gridfunc.hpp:948
void GetGradients(ElementTransformation &tr, const IntegrationRule &ir, DenseMatrix &grad) const
Definition: gridfunc.cpp:1636
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:2598
void LegacyNCReorder()
Loading helper.
Definition: gridfunc.cpp:3724
QuadratureFunction(const QuadratureFunction &orig)
Copy constructor. The QuadratureSpace ownership flag, own_qspace, in the new object is set to false...
Definition: gridfunc.hpp:749
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:2278
virtual void SetSpace(FiniteElementSpace *f)
Associate a new FiniteElementSpace with the GridFunction.
Definition: gridfunc.cpp:193
void ReorderByNodes()
For a vector grid function, makes sure that the ordering is byNODES.
Definition: gridfunc.cpp:1203
std::ostream & operator<<(std::ostream &os, SparseMatrix const &mat)
Definition: sparsemat.hpp:655
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.cpp:3050
void GetCurl(ElementTransformation &tr, Vector &curl) const
Definition: gridfunc.cpp:1466
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:111
virtual void ComputeFlux(BilinearFormIntegrator &blfi, GridFunction &flux, bool wcoef=true, int subdomain=-1)
Definition: gridfunc.cpp:296
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:2827
Class representing the storage layout of a QuadratureFunction.
Definition: fespace.hpp:909
double u(const Vector &xvec)
Definition: lor_mms.hpp:24
void SetSize(int s)
Change the size of the DenseMatrix to s x s.
Definition: densemat.hpp:105
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:872
double Eval(double h, int p) const
Definition: gridfunc.hpp:722
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
Definition: gridfunc.cpp:347
Class representing a function through its values (scalar or vector) at quadrature points...
Definition: gridfunc.hpp:734
void GetValuesFrom(const GridFunction &orig_func)
Definition: gridfunc.cpp:1087
void GetBdrValuesFrom(const GridFunction &orig_func)
Definition: gridfunc.cpp:1124
virtual void ComputeElementMaxErrors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Definition: gridfunc.hpp:608
void AccumulateAndCountBdrValues(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr, Array< int > &values_counter)
Definition: gridfunc.cpp:1997
void Save(std::ostream &out) const
Write the QuadratureFunction to the stream out.
Definition: gridfunc.cpp:3819
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:3525
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:424
virtual double GetValue(int i, const IntegrationPoint &ip, int vdim=1) const
Definition: gridfunc.cpp:402
void GetGradient(ElementTransformation &tr, Vector &grad) const
Definition: gridfunc.cpp:1568
void ProjectVectorFieldOn(GridFunction &vec_field, int comp=0)
Definition: gridfunc.cpp:1261