12#ifndef MFEM_TRANSFER_HPP
13#define MFEM_TRANSFER_HPP
157 bool own_mass_integ_ =
true);
425 std::pair<std::unique_ptr<SparseMatrix>,
454 std::unique_ptr<SparseMatrix>
AllocR();
462 std::unique_ptr<Operator>
R;
464 std::unique_ptr<Operator>
M_LH;
518 bool force_l2_space_ =
false,
608 const Operator* elem_restrict_lex_l;
609 const Operator* elem_restrict_lex_h;
Conjugate gradient method.
Data type dense matrix using column-major storage.
static MemoryType GetHostMemoryType()
Get the current Host MemoryType. This is the MemoryType used by most MFEM classes when allocating mem...
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Abstract class for all finite elements.
Base class for transfer algorithms that construct transfer Operators between two finite element (FE) ...
virtual ~GridTransfer()
Virtual destructor.
FiniteElementSpace & dom_fes
Domain FE space.
OperatorHandle fw_t_oper
Forward true-dof operator.
FiniteElementSpace & ran_fes
Range FE space.
virtual const Operator & TrueBackwardOperator()
Return an Operator that transfers true-dof Vectors from the range FE space back to true-dof Vectors i...
OperatorHandle bw_t_oper
Backward true-dof operator.
virtual bool SupportsBackwardsOperator() const
void SetMemType(MemoryType d_mt_)
virtual const Operator & ForwardOperator()=0
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
void SetOperatorType(Operator::Type type)
Set the desired Operator::Type for the construction of all operators defined by the underlying transf...
virtual const Operator & BackwardOperator()=0
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Operator::Type oper_type
Desired Operator::Type for the construction of all operators defined by the underlying transfer algor...
virtual const Operator & TrueForwardOperator()
Return an Operator that transfers true-dof Vectors from the domain FE space to true-dof Vectors in th...
const Operator & MakeTrueOperator(FiniteElementSpace &fes_in, FiniteElementSpace &fes_out, const Operator &oper, OperatorHandle &t_oper)
GridTransfer(FiniteElementSpace &dom_fes_, FiniteElementSpace &ran_fes_)
Transfer data between a coarse mesh and an embedded refined mesh using interpolation.
virtual ~InterpolationGridTransfer()
InterpolationGridTransfer(FiniteElementSpace &coarse_fes, FiniteElementSpace &fine_fes)
BilinearFormIntegrator * mass_integ
Ownership depends on own_mass_integ.
OperatorHandle F
Forward, coarse-to-fine, operator.
bool own_mass_integ
Ownership flag for mass_integ.
const Operator & BackwardOperator() override
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
OperatorHandle B
Backward, fine-to-coarse, operator.
void SetMassIntegrator(BilinearFormIntegrator *mass_integ_, bool own_mass_integ_=true)
Assign a mass integrator to be used in the construction of the backward, fine-to-coarse,...
const Operator & ForwardOperator() override
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
const FiniteElementSpace * fes_lor
H1SpaceLumpedMassOperator(const FiniteElementSpace *fes_ho_, const FiniteElementSpace *fes_lor_, Vector &ML_inv_)
const FiniteElementSpace * fes_ho
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
const FiniteElementSpace * fes_lor
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
const FiniteElementSpace * fes_ho
H1SpaceMixedMassOperator(const FiniteElementSpace *fes_ho_, const FiniteElementSpace *fes_lor_, Table *ho2lor_, Vector *M_LH_ea_)
Class below must be public as we now have device code.
void Mult(const Vector &x, Vector &y) const override
void EAL2ProjectionH1Space()
void SetFromTDofsTranspose(const FiniteElementSpace &fes, const Vector &X, Vector &x) const
Sets dual field coefficients given a vector of dual field coefficients on the tdofs and a finite elem...
std::unique_ptr< FiniteElementSpace > fes_ho_scalar
std::unique_ptr< ParFiniteElementSpace > pfes_ho_scalar
std::unique_ptr< Operator > ML_inv_vea
void TDofsListByVDim(const FiniteElementSpace &fes, int vdim, Array< int > &vdofs_list) const
Fills the vdofs_list array with a list of vdofs for a given vdim and a given finite element space.
std::unique_ptr< Solver > precon
void SetupPCG()
Sets up the PCG solver (sets parameters, operator, and preconditioner)
void Prolongate(const Vector &x, Vector &y) const override
void GetTDofs(const FiniteElementSpace &fes, const Vector &x, Vector &X) const
Recovers vector of tdofs given a vector of dofs and a finite element space.
void SetFromTDofs(const FiniteElementSpace &fes, const Vector &X, Vector &x) const
Sets dof values given a vector of tdofs and a finite element space.
std::unique_ptr< Operator > RTxM_LH
void GetTDofsTranspose(const FiniteElementSpace &fes, const Vector &x, Vector &X) const
Recovers a vector of dual field coefficients on the tdofs given a vector of dual coefficients and a f...
L2ProjectionH1Space(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, const bool use_ea_, MemoryType d_mt_=Device::GetHostMemoryType())
std::unique_ptr< SparseMatrix > AllocR()
Computes sparsity pattern and initializes R matrix. Based on BilinearForm::AllocMat(),...
std::unique_ptr< Operator > M_LH
void MultTranspose(const Vector &x, Vector &y) const override
void LumpedMassInverse(Vector &ML_inv) const
Returns the inverse of an on-rank lumped mass matrix.
std::unique_ptr< FiniteElementSpace > fes_lor_scalar
void SetAbsTol(real_t p_atol_) override
Sets absolute tolerance in preconditioned conjugate gradient solver.
void ProlongateTranspose(const Vector &x, Vector &y) const override
std::pair< std::unique_ptr< SparseMatrix >, std::unique_ptr< SparseMatrix > > ComputeSparseRAndM_LH()
Computes on-rank R and M_LH matrices. If true, computes mixed mass and/or inverse lumped mass matrix ...
std::unique_ptr< Operator > R
void SetRelTol(real_t p_rtol_) override
Sets relative tolerance in preconditioned conjugate gradient solver.
std::unique_ptr< ParFiniteElementSpace > pfes_lor_scalar
void ProlongateTranspose(const Vector &x, Vector &y) const override
L2ProjectionL2Space(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, const bool use_ea_, MemoryType d_mt_=Device::GetHostMemoryType())
void EAMultTranspose(const Vector &x, Vector &y) const
void EAProlongate(const Vector &x, Vector &y) const
void MultTranspose(const Vector &x, Vector &y) const override
void EAMult(const Vector &x, Vector &y) const
Perform mult on the device (same as above)
void EAL2ProjectionL2Space()
void SetRelTol(real_t p_rtol_) override
No-op.
void EAProlongateTranspose(const Vector &x, Vector &y) const
void SetAbsTol(real_t p_atol_) override
No-op.
void Prolongate(const Vector &x, Vector &y) const override
void Mult(const Vector &x, Vector &y) const override
void BuildHo2Lor(int nel_ho, int nel_lor, const CoarseFineTransformations &cf_tr)
const FiniteElementSpace & fes_ho
void MixedMassEA(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, Vector &M_LH, MemoryType d_mt_=Device::GetHostMemoryType())
virtual void ProlongateTranspose(const Vector &x, Vector &y) const =0
virtual void SetRelTol(real_t p_rtol_)=0
Sets relative tolerance in preconditioned conjugate gradient solver.
const FiniteElementSpace & fes_lor
L2Projection(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, MemoryType d_mt_=Device::GetHostMemoryType())
void ElemMixedMass(Geometry::Type geom, const FiniteElement &fe_ho, const FiniteElement &fe_lor, ElementTransformation *tr_ho, ElementTransformation *tr_lor, IntegrationPointTransformation &ip_tr, DenseMatrix &M_mixed_el) const
virtual void SetAbsTol(real_t p_atol_)=0
Sets absolute tolerance in preconditioned conjugate gradient solver.
virtual void Prolongate(const Vector &x, Vector &y) const =0
L2Prolongation(const L2Projection &l2proj_)
void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
virtual ~L2Prolongation()
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Transfer data in L2 and H1 finite element spaces between a coarse mesh and an embedded refined mesh u...
L2Projection * F
Forward, coarse-to-fine, operator.
const Operator & BackwardOperator() override
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
bool SupportsBackwardsOperator() const override
L2Prolongation * B
Backward, fine-to-coarse, operator.
const Operator & ForwardOperator() override
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
L2ProjectionGridTransfer(FiniteElementSpace &coarse_fes_, FiniteElementSpace &fine_fes_, bool force_l2_space_=false, MemoryType d_mt_=Device::GetHostMemoryType())
virtual ~L2ProjectionGridTransfer()
Pointer to an Operator of a specified type.
int Height() const
Get the height (size of output) of the Operator. Synonym with NumRows().
int Width() const
Get the width (size of input) of the Operator. Synonym with NumCols().
Type
Enumeration defining IDs for some classes derived from Operator.
Matrix-free transfer operator between finite element spaces on the same mesh.
PRefinementTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator from lFESpace to hFESpace which have different FE collections.
virtual ~PRefinementTransferOperator()
Destructor.
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
Abstract parallel finite element space.
Matrix-free transfer operator between finite element spaces on the same mesh exploiting the tensor pr...
virtual ~TensorProductPRefinementTransferOperator()
Destructor.
void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
TensorProductPRefinementTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator from lFESpace to hFESpace which have different FE collections.
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
Matrix-free transfer operator between finite element spaces.
virtual ~TransferOperator()
Destructor.
void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
TransferOperator(const FiniteElementSpace &lFESpace, const FiniteElementSpace &hFESpace)
Constructs a transfer operator from lFESpace to hFESpace.
Matrix-free transfer operator between finite element spaces working on true degrees of freedom.
~TrueTransferOperator()
Destructor.
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
TrueTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator working on true degrees of freedom from lFESpace to hFESpace.
void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a true dof vector x to a true dof vector y.
MemoryType
Memory types supported by MFEM.