12 #ifndef MFEM_PGRIDFUNC
13 #define MFEM_PGRIDFUNC
15 #include "../config/config.hpp"
19 #include "../general/globals.hpp"
29 double GlobalLpNorm(
const double p,
double loc_norm, MPI_Comm comm);
82 const int *partitioning = NULL);
214 int comp = 0,
Vector *tr = NULL)
const;
394 p, exsol, weight, v_weight, irs),
pfes->
GetComm());
399 bool wcoef =
true,
int subdomain = -1);
404 virtual void Save(std::ostream &
out)
const;
406 #ifdef MFEM_USE_ADIOS2
428 const ParGridFunction &x,
429 ParFiniteElementSpace &smooth_flux_fes,
430 ParFiniteElementSpace &flux_fes,
431 Vector &errors,
int norm_p = 2,
double solver_tol = 1e-12,
432 int solver_max_it = 200);
436 #endif // MFEM_USE_MPI
virtual double ComputeLpError(const double p, VectorCoefficient &exsol, Coefficient *weight=NULL, VectorCoefficient *v_weight=NULL, const IntegrationRule *irs[]=NULL) const
Class for an integration rule - an Array of IntegrationPoint.
Class for grid function - Vector with associated FE space.
void ExchangeFaceNbrData()
void AddDistribute(double a, const Vector &tv)
ParGridFunction(const ParGridFunction &orig)
Copy constructor. The internal vector face_nbr_data is not copied.
HypreParVector * ParallelAssemble() const
Returns a new vector assembled on the true dofs.
Base class for vector Coefficients that optionally depend on time and space.
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
HypreParVector * ParallelProject() const
Returns a new vector restricted to the true dofs.
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
virtual ~ParGridFunction()
void ProjectDiscCoefficient(VectorCoefficient &coeff, Array< int > &dof_attr)
ParFiniteElementSpace * pfes
Points to the same object as fes.
ParGridFunction & operator=(const HypreParVector &tv)
Short semantic for Distribute()
void SaveAsOne(std::ostream &out=mfem::out)
Merge the local grid functions.
ParGridFunction(ParFiniteElementSpace *pf)
virtual double GetValue(int i, const IntegrationPoint &ip, int vdim=1) const
virtual void Save(std::ostream &out) const
virtual void ProjectDiscCoefficient(VectorCoefficient &coeff)
Project a discontinuous vector coefficient as a grid function on a continuous finite element space...
virtual void SetSpace(FiniteElementSpace *f)
Associate a new FiniteElementSpace with the ParGridFunction.
Abstract parallel finite element space.
virtual double ComputeHDivError(VectorCoefficient *exsol, Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns the error measured in H(div)-norm for RT elements.
virtual void ProjectCoefficient(Coefficient &coeff)
Vector send_data
Vector used as an MPI buffer to send face-neighbor data in ExchangeFaceNbrData() to neighboring proce...
virtual double ComputeDGFaceJumpError(Coefficient *exsol, Coefficient *ell_coeff, double Nu, const IntegrationRule *irs[]=NULL) const
Returns the Face Jumps error for L2 elements.
void Distribute(const Vector &tv)
virtual void ProjectBdrCoefficient(VectorCoefficient &vcoeff, Array< int > &attr)
Project a VectorCoefficient on the GridFunction, modifying only DOFs on the boundary associated with ...
virtual double ComputeCurlError(VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns ||curl u_ex - curl u_h||_L2 for ND elements.
virtual double ComputeHDivError(VectorCoefficient *exsol, Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns the error measured H(div)-norm for RT elements.
double GetValue(ElementTransformation &T)
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the ParGridFunction reference external data on a new FiniteElementSpace.
ParGridFunction & operator=(double value)
Assign constant values to the ParGridFunction data.
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
const Vector & FaceNbrData() const
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
virtual double ComputeL1Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
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...
virtual double ComputeLpError(const double p, Coefficient &exsol, Coefficient *weight=NULL, const IntegrationRule *irs[]=NULL) const
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
ParGridFunction & operator=(const Vector &v)
Copy the data from a Vector to the ParGridFunction data.
virtual void GetVectorValue(int i, const IntegrationPoint &ip, Vector &val) const
Wrapper for hypre's parallel vector class.
virtual double ComputeL1Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
virtual double ComputeGradError(VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements.
double GlobalLpNorm(const double p, double loc_norm, MPI_Comm comm)
Compute a global Lp norm from the local Lp norms computed by each processor.
double p(const Vector &x, double t)
ParGridFunction & operator=(const ParGridFunction &rhs)
Copy assignment. Only the data of the base class Vector is copied.
virtual double ComputeL2Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
virtual void ComputeFlux(BilinearFormIntegrator &blfi, GridFunction &flux, bool wcoef=true, int subdomain=-1)
virtual double ComputeMaxError(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
double L2ZZErrorEstimator(BilinearFormIntegrator &flux_integrator, const ParGridFunction &x, ParFiniteElementSpace &smooth_flux_fes, ParFiniteElementSpace &flux_fes, Vector &errors, int norm_p, double solver_tol, int solver_max_it)
void Distribute(const Vector *tv)
virtual double ComputeGradError(VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
Returns ||grad u_ex - grad u_h||_L2 for H1 or L2 elements.
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
virtual double ComputeMaxError(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
virtual double ComputeMaxError(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL) const
virtual double ComputeHCurlError(VectorCoefficient *exsol, VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns the error measured in H(curl)-norm for ND elements.
virtual void ProjectBdrCoefficient(Coefficient *coeff[], Array< int > &attr)
Project a set of Coefficients on the components of the GridFunction, modifying only DOFs on the bound...
Vector face_nbr_data
Vector used to store data from face-neighbor processors, initialized by ExchangeFaceNbrData().
Class for integration point with weight.
virtual double ComputeL2Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL) const
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, const IntegrationRule *irs[]=NULL) const
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
HypreParVector * ParallelAverage() const
Returns a new vector averaged on the true dofs.
HypreParVector * GetTrueDofs() const
Returns the true dofs in a new HypreParVector.
void AddDistribute(double a, const Vector *tv)
double infinity()
Define a shortcut for std::numeric_limits<double>::infinity()
virtual void ProjectCoefficient(Coefficient &coeff)
virtual double ComputeL1Error(Coefficient &exsol, const IntegrationRule *irs[]=NULL) const
virtual double ComputeW11Error(Coefficient *exsol, VectorCoefficient *exgrad, int norm_type, Array< int > *elems=NULL, const IntegrationRule *irs[]=NULL) const
virtual double ComputeL2Error(VectorCoefficient &exsol, const IntegrationRule *irs[]=NULL, Array< int > *elems=NULL) const
virtual double ComputeH1Error(Coefficient *exsol, VectorCoefficient *exgrad, Coefficient *ell_coef, double Nu, int norm_type) const
virtual double ComputeDivError(Coefficient *exdiv, const IntegrationRule *irs[]=NULL) const
Returns ||div u_ex - div u_h||_L2 for RT elements.
Class for parallel grid function.
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Class for parallel meshes.
ParGridFunction(ParFiniteElementSpace *pf, double *data)
Construct a ParGridFunction using previously allocated array data.
virtual double ComputeHCurlError(VectorCoefficient *exsol, VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const
Returns the error measured H(curl)-norm for ND elements.
void ProjectBdrCoefficient(Coefficient &coeff, Array< int > &attr)
Project a Coefficient on the GridFunction, modifying only DOFs on the boundary associated with the bo...
double f(const Vector &p)
ParFiniteElementSpace * ParFESpace() const