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