15 #include "../config/config.hpp"
19 #ifdef MFEM_USE_ADIOS2
20 #include "../general/adios2stream.hpp"
49 void SaveSTLTri(std::ostream &
out,
double p1[],
double p2[],
double p3[]);
92 :
Vector(data, f->GetVSize())
213 int comp = 0,
Vector *tr = NULL)
const;
335 double _min = 0.0,
double _max =
infinity());
469 int norm_type)
const;
582 bool wcoef =
true,
int subdomain = -1);
640 virtual void Save(std::ostream &
out)
const;
642 #ifdef MFEM_USE_ADIOS2
652 void SaveVTK(std::ostream &
out,
const std::string &field_name,
int ref);
656 void SaveSTL(std::ostream &
out,
int TimesToRefine = 1);
665 std::ostream &
operator<<(std::ostream &
out,
const GridFunction &sol);
693 :
Vector(vdim_*qspace_->GetSize()),
701 :
Vector(qf_data, vdim_*qspace_->GetSize()),
817 void Save(std::ostream &
out)
const;
821 std::ostream &
operator<<(std::ostream &
out,
const QuadratureFunction &qf);
827 Vector &error_estimates,
828 Array<int> *aniso_flags = NULL,
829 int with_subdomains = 1,
830 bool with_coeff =
false);
834 GridFunction& gf1, GridFunction& gf2);
846 : n(_n), mesh_in(m), sol_in(s) { }
853 GridFunction *sol,
const int ny);
871 double *qf_data,
int vdim_)
895 const double *q =
data +
vdim*s_offset;
896 for (
int i = 0; i<values.
Size(); i++)
914 const double *q =
data + s_offset;
915 for (
int i = 0; i < values.
Size(); i++)
934 const double *q =
data +
vdim*s_offset;
935 for (
int j = 0; j<sl_size; j++)
937 for (
int i = 0; i<
vdim; i++)
939 values(i,j) = *(q++);
void SaveSTLTri(std::ostream &out, double p1[], double p2[], double p3[])
int GetVDim() const
Get the vector dimension.
virtual ~GridFunction()
Destroys grid function.
QuadratureFunction(QuadratureSpace *qspace_, int vdim_=1)
Create a QuadratureFunction based on the given QuadratureSpace.
Class for an integration rule - an Array of IntegrationPoint.
void GetElementAverages(GridFunction &avgs) const
void GetVectorValues(int i, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr) const
void NewDataAndSize(double *d, int s)
Set the Vector data and size, deleting the old data, if owned.
Class for grid function - Vector with associated FE space.
const IntegrationRule & GetElementIntRule(int idx) const
Get the IntegrationRule associated with mesh element idx.
Class used for extruding scalar GridFunctions.
virtual double ComputeDGFaceJumpError(Coefficient *exsol, Coefficient *ell_coeff, double Nu, const IntegrationRule *irs[]=NULL) const
Returns the Face Jumps error for L2 elements.
void SumFluxAndCount(BilinearFormIntegrator &blfi, GridFunction &flux, Array< int > &counts, bool wcoef, int subdomain)
void SetFromTrueVector()
Shortcut for calling SetFromTrueDofs() with GetTrueVector() as argument.
Base class for vector Coefficients that optionally depend on time and space.
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
void ProjectDeltaCoefficient(DeltaCoefficient &delta_coeff, double &integral)
const IntegrationRule & GetElementIntRule(int idx) const
Get the IntegrationRule associated with mesh element idx.
virtual void ComputeElementLpErrors(const double p, Coefficient &exsol, Vector &error, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
virtual void GetVectorValue(int i, const IntegrationPoint &ip, Vector &val) const
void SetSize(int s)
Resize the vector to size s.
void RestrictConforming()
GridFunction * Extrude1DGridFunction(Mesh *mesh, Mesh *mesh2d, GridFunction *sol, const int ny)
Extrude a scalar 1D GridFunction, after extruding the mesh with Extrude1D.
void MakeOwner(FiniteElementCollection *_fec)
Make the GridFunction the owner of fec and fes.
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
bool own_qspace
QuadratureSpace ownership flag.
int GetSize() const
Return the total number of quadrature points.
virtual ~QuadratureFunction()
Data type dense matrix using column-major storage.
Delta function coefficient optionally multiplied by a weight coefficient and a scaled time dependent ...
int Size() const
Returns the size of the vector.
void ImposeBounds(int i, const Vector &weights, const Vector &_lo, const Vector &_hi)
void AccumulateAndCountBdrTangentValues(VectorCoefficient &vcoeff, Array< int > &bdr_attr, Array< int > &values_counter)
virtual double ComputeHDivError(VectorCoefficient *exsol, Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns the error measured in H(div)-norm for RT elements.
int vdim
Vector dimension.
void GetVectorGradientHat(ElementTransformation &T, DenseMatrix &gh) const
Vector & GetTrueVector()
Read and write access to the (optional) internal true-dof Vector.
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
void SetSpace(QuadratureSpace *qspace_, int vdim_=-1)
Change the QuadratureSpace and optionally the vector dimension.
void ComputeMeans(AvgType type, Array< int > &zones_per_vdof)
bool OwnsSpace()
Get the QuadratureSpace ownership flag.
void GetNodalValues(int i, Array< double > &nval, int vdim=1) const
Returns the values in the vertices of i'th element for dimension vdim.
virtual void ComputeElementL2Errors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
QuadratureFunction(QuadratureSpace *qspace_, double *qf_data, int vdim_=1)
Create a QuadratureFunction based on the given QuadratureSpace, using the external data...
void GetVectorFieldNodalValues(Vector &val, int comp) const
virtual void Save(std::ostream &out) const
Save the GridFunction to an output stream.
const FiniteElementSpace * FESpace() const
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
GridFunction(FiniteElementSpace *f, double *data)
Construct a GridFunction using previously allocated array data.
void GetDerivative(int comp, int der_comp, GridFunction &der)
void ProjectBdrCoefficientNormal(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
void SetTrueVector()
Shortcut for calling GetTrueDofs() with GetTrueVector() as argument.
void GetHessians(int i, const IntegrationRule &ir, DenseMatrix &hess, int vdim=1) const
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
bool UseDevice() const
Return the device flag of the Memory object used by the Vector.
void GetGradients(const int elem, const IntegrationRule &ir, DenseMatrix &grad) const
void Reset(double *d, int h, int w)
Change the data array and the size of the DenseMatrix.
void GetVectorFieldValues(int i, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr, int comp=0) const
void SetOwnsSpace(bool own)
Set the QuadratureSpace ownership flag.
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
GridFunction(const GridFunction &orig)
Copy constructor. The internal true-dof vector t_vec is not copied.
void MakeTRef(FiniteElementSpace *f, double *tv)
Associate a new FiniteElementSpace and new true-dof data with the GridFunction.
QuadratureSpace * GetSpace() const
Get the associated QuadratureSpace.
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
QuadratureFunction()
Create an empty QuadratureFunction.
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 ...
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
void SetVDim(int vdim_)
Set the vector dimension, updating the size by calling Vector::SetSize().
void GetLaplacians(int i, const IntegrationRule &ir, Vector &laps, int vdim=1) const
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...
void GetValues(int i, const IntegrationRule &ir, Vector &vals, int vdim=1) const
FiniteElementCollection * OwnFEC()
virtual ~ExtrudeCoefficient()
FiniteElementSpace * FESpace()
GridFunction & operator=(const GridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
virtual void ComputeElementL1Errors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
QuadratureFunction & operator=(double value)
Redefine '=' for QuadratureFunction = constant.
const Vector & GetTrueVector() const
Read only access to the (optional) internal true-dof Vector.
ElementTransformation * GetElementTransformation(int i) const
Returns ElementTransformation for the i-th element.
virtual double ComputeGradError(VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements.
virtual void ComputeElementL2Errors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
double p(const Vector &x, double t)
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
virtual void ComputeElementL1Errors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
FiniteElementSpace * fes
FE space on which the grid function lives. Owned if fec is not NULL.
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
void ProjectGridFunction(const GridFunction &src)
Project the src GridFunction to this GridFunction, both of which must be on the same mesh...
QuadratureSpace * qspace
Associated QuadratureSpace.
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
FiniteElementCollection * fec
Used when the grid function is read from a file. It can also be set explicitly, see MakeOwner()...
void GetVectorGradient(ElementTransformation &tr, DenseMatrix &grad) const
double GetDivergence(ElementTransformation &tr) const
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
double ComputeElementLpDistance(double p, int i, GridFunction &gf1, GridFunction &gf2)
Compute the Lp distance between two grid functions on the given element.
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
int GetFaceVectorValues(int i, int side, const IntegrationRule &ir, DenseMatrix &vals, DenseMatrix &tr) const
virtual double ComputeHCurlError(VectorCoefficient *exsol, VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns the error measured in H(curl)-norm for ND elements.
double ZZErrorEstimator(BilinearFormIntegrator &blfi, GridFunction &u, GridFunction &flux, Vector &error_estimates, Array< int > *aniso_flags, int with_subdomains, bool with_coeff)
virtual void ComputeElementMaxErrors(Coefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
Class for integration point with weight.
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...
int GetFaceValues(int i, int side, const IntegrationRule &ir, Vector &vals, DenseMatrix &tr, int vdim=1) const
void GetElementValues(int idx, Vector &values)
Return all values associated with mesh element idx in a Vector.
void GetGradients(ElementTransformation &tr, const IntegrationRule &ir, DenseMatrix &grad) const
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
QuadratureFunction(const QuadratureFunction &orig)
Copy constructor. The QuadratureSpace ownership flag, own_qspace, in the new object is set to false...
double infinity()
Define a shortcut for std::numeric_limits<double>::infinity()
virtual void ProjectCoefficient(Coefficient &coeff)
virtual void SetSpace(FiniteElementSpace *f)
Associate a new FiniteElementSpace with the GridFunction.
void ReorderByNodes()
For a vector grid function, makes sure that the ordering is byNODES.
std::ostream & operator<<(std::ostream &os, SparseMatrix const &mat)
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
ExtrudeCoefficient(Mesh *m, Coefficient &s, int _n)
void GetCurl(ElementTransformation &tr, Vector &curl) const
void MakeRef(Vector &base, int offset, int size)
Reset the Vector to be a reference to a sub-vector of base.
virtual void ComputeFlux(BilinearFormIntegrator &blfi, GridFunction &flux, bool wcoef=true, int subdomain=-1)
GridFunction(FiniteElementSpace *f)
Construct a GridFunction associated with the FiniteElementSpace *f.
Class representing the storage layout of a QuadratureFunction.
void SetSize(int s)
Change the size of the DenseMatrix to s x s.
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
long GetSequence() const
Return update counter (see Mesh::sequence)
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
Class representing a function through its values (scalar or vector) at quadrature points...
void GetValuesFrom(const GridFunction &orig_func)
void GetBdrValuesFrom(const GridFunction &orig_func)
virtual void ComputeElementMaxErrors(VectorCoefficient &exsol, Vector &error, const IntegrationRule *irs[]=NULL) const
void AccumulateAndCountBdrValues(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr, Array< int > &values_counter)
void Save(std::ostream &out) const
Write the QuadratureFunction to the stream out.
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 > 0 must match the one used in Mesh::PrintVTK.
void ProjectBdrCoefficient(Coefficient &coeff, Array< int > &attr)
Project a Coefficient on the GridFunction, modifying only DOFs on the boundary associated with the bo...
virtual double GetValue(int i, const IntegrationPoint &ip, int vdim=1) const
void GetGradient(ElementTransformation &tr, Vector &grad) const
double f(const Vector &p)
void ProjectVectorFieldOn(GridFunction &vec_field, int comp=0)