phynexis-deps/mfem
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge pull request #4851 from mfem/najlkin/warn-gridfunc

Browse Source 
Vector dimension for bdr/face elements
This commit is contained in:
Andrew Ho
2025-06-02 11:45:47 -07:00
committed by GitHub
parent 1730b05078 8a522f5e7d
commit a443f003bb
7 changed files with 119 additions and 31 deletions
Download Patch File Download Diff File Expand all files Collapse all files
fem/fespace.cpp
+42 -4
View File
@@ -1490,10 +1490,8 @@ int FiniteElementSpace::GetNConformingDofs() const
int FiniteElementSpace::GetVectorDim() const
{
const FiniteElement *fe = GetTypicalFE();
if (fe == nullptr)
{
fe = GetTypicalTraceElement();
}
MFEM_VERIFY(fe, "A typical finite element does not exist!");
if (fe->GetRangeType() == FiniteElement::SCALAR)
{
return GetVDim();
@@ -1501,6 +1499,31 @@ int FiniteElementSpace::GetVectorDim() const
return GetVDim()*std::max(GetMesh()->SpaceDimension(), fe->GetRangeDim());
}
int FiniteElementSpace::GetBdrVectorDim() const
{
const FiniteElement *be = GetTypicalBE();
MFEM_VERIFY(be, "A typical boundary finite element does not exist!");
if (be->GetRangeType() == FiniteElement::SCALAR)
{
return GetVDim();
}
return GetVDim()*std::max(GetMesh()->SpaceDimension()-1, be->GetRangeDim());
}
int FiniteElementSpace::GetFaceVectorDim() const
{
const FiniteElement *face_el = GetTypicalFaceElement();
MFEM_VERIFY(face_el, "A typical face finite element does not exist!");
if (face_el->GetRangeType() == FiniteElement::SCALAR)
{
return GetVDim();
}
return GetVDim()*std::max(GetMesh()->SpaceDimension()-1,
face_el->GetRangeDim());
}
int FiniteElementSpace::GetCurlDim() const
{
const FiniteElement *fe = GetTypicalFE();
@@ -3915,6 +3938,16 @@ const FiniteElement *FiniteElementSpace::GetBE(int i) const
return BE;
}
const FiniteElement *FiniteElementSpace::GetTypicalBE() const
{
if (mesh->GetNBE() > 0) { return GetBE(0); }
Geometry::Type geom = mesh->GetTypicalFaceGeometry();
const FiniteElement *be = fec->FiniteElementForGeometry(geom);
MFEM_VERIFY(be != nullptr, "Could not determine a typical BE!");
return be;
}
const FiniteElement *FiniteElementSpace::GetFaceElement(int i) const
{
MFEM_VERIFY(!IsVariableOrder(), "not implemented");
@@ -3945,6 +3978,11 @@ const FiniteElement *FiniteElementSpace::GetFaceElement(int i) const
return fe;
}
const FiniteElement *FiniteElementSpace::GetTypicalFaceElement() const
{
return fec->FiniteElementForGeometry(mesh->GetTypicalFaceGeometry());
}
const FiniteElement *FiniteElementSpace::GetEdgeElement(int i,
int variant) const
{
fem/fespace.hpp
+20
View File
@@ -868,6 +868,14 @@ public:
*/
int GetVectorDim() const;
/// Return the total dimension of a vector on boundary
/** @see GetVectorDim() */
int GetBdrVectorDim() const;
/// Return the total dimension of a vector on a face
/** @see GetVectorDim() */
int GetFaceVectorDim() const;
/// Return the dimension of the curl of a GridFunction defined on this space.
/** Note: This assumes a space dimension of 2 or 3 only. */
int GetCurlDim() const;
@@ -1333,12 +1341,24 @@ public:
associated with i'th boundary face in the mesh object. */
const FiniteElement *GetBE(int i) const;
/// @brief Return a "typical" boundary element.
///
/// This can be used in situations where the local mesh partition may be
/// empty.
const FiniteElement *GetTypicalBE() const;
/** @brief Returns pointer to the FiniteElement in the FiniteElementCollection
associated with i'th face in the mesh object. Faces in this case refer
to the MESHDIM-1 primitive so in 2D they are segments and in 1D they are
points.*/
const FiniteElement *GetFaceElement(int i) const;
/// @brief Return a "typical" face element.
///
/// This can be used in situations where the local mesh partition may be
/// empty.
const FiniteElement *GetTypicalFaceElement() const;
/** @brief Returns pointer to the FiniteElement in the FiniteElementCollection
associated with i'th edge in the mesh object. */
const FiniteElement *GetEdgeElement(int i, int variant = 0) const;
fem/gridfunc.cpp
+7 -15
View File
@@ -351,16 +351,6 @@ void GridFunction::ComputeFlux(BilinearFormIntegrator &blfi,
}
}
int GridFunction::VectorDim() const
{
return fes->GetVectorDim();
}
int GridFunction::CurlDim() const
{
return fes->GetCurlDim();
}
void GridFunction::GetTrueDofs(Vector &tv) const
{
const SparseMatrix *R = fes->GetRestrictionMatrix();
@@ -2660,7 +2650,8 @@ void GridFunction::ProjectDiscCoefficient(VectorCoefficient &coeff,
void GridFunction::ProjectBdrCoefficient(VectorCoefficient &vcoeff,
const Array<int> &attr)
{
MFEM_VERIFY(VectorDim() == vcoeff.GetVDim(), "vcoeff vdim != VectorDim()");
MFEM_VERIFY(BdrVectorDim() == vcoeff.GetVDim(),
"vcoeff vdim != BdrVectorDim()");
Array<int> values_counter;
AccumulateAndCountBdrValues(NULL, &vcoeff, attr, values_counter);
ComputeMeans(ARITHMETIC, values_counter);
@@ -2709,9 +2700,9 @@ void GridFunction::ProjectBdrCoefficient(Coefficient *coeff[],
void GridFunction::ProjectBdrCoefficientNormal(
VectorCoefficient &vcoeff, const Array<int> &bdr_attr)
{
MFEM_VERIFY(fes->GetVDim() * fes->GetMesh()->SpaceDimension() ==
vcoeff.GetVDim(),
"vcoeff vdim mismatch");
MFEM_VERIFY(BdrVectorDim() == 1, "BdrVectorDim() != 1");
MFEM_VERIFY(vcoeff.GetVDim() == fes->GetMesh()->SpaceDimension(),
"vcoeff vdim != space dim");
#if 0
// implementation for the case when the face dofs are integrals of the
// normal component.
@@ -2787,7 +2778,8 @@ void GridFunction::ProjectBdrCoefficientNormal(
void GridFunction::ProjectBdrCoefficientTangent(
VectorCoefficient &vcoeff, const Array<int> &bdr_attr)
{
MFEM_VERIFY(VectorDim() == vcoeff.GetVDim(), "vcoeff vdim != VectorDim()");
MFEM_VERIFY(BdrVectorDim()+1 == vcoeff.GetVDim(),
"vcoeff vdim != BdrVectorDim()+1");
Array<int> values_counter;
AccumulateAndCountBdrTangentValues(vcoeff, bdr_attr, values_counter);
ComputeMeans(ARITHMETIC, values_counter);
fem/gridfunc.hpp
+14 -4
View File
@@ -123,11 +123,21 @@ public:
FiniteElementCollection *OwnFEC() { return fec_owned; }
/// Shortcut for calling FiniteElementSpace::GetVectorDim() on the underlying #fes
int VectorDim() const;
/** @brief Shortcut for calling FiniteElementSpace::GetVectorDim() on the
underlying #fes */
int VectorDim() const { return fes->GetVectorDim(); }
/// Shortcut for calling FiniteElementSpace::GetCurlDim() on the underlying #fes
int CurlDim() const;
/** @brief Shortcut for calling FiniteElementSpace::GetBdrVectorDim() on the
underlying #fes */
int BdrVectorDim() const { return fes->GetBdrVectorDim(); }
/** @brief Shortcut for calling FiniteElementSpace::GetFaceVectorDim() on
the underlying #fes */
int FaceVectorDim() const { return fes->GetFaceVectorDim(); }
/** @brief Shortcut for calling FiniteElementSpace::GetCurlDim() on the
underlying #fes */
int CurlDim() const { return fes->GetCurlDim(); }
/// Read only access to the (optional) internal true-dof Vector.
const Vector &GetTrueVector() const
fem/pgridfunc.cpp
+4 -2
View File
@@ -740,14 +740,16 @@ void ParGridFunction::ProjectBdrCoefficient(
void ParGridFunction::ProjectBdrCoefficient(VectorCoefficient &vcoeff,
const Array<int> &attr)
{
MFEM_VERIFY(VectorDim() == vcoeff.GetVDim(), "vcoeff vdim != VectorDim()");
MFEM_VERIFY(BdrVectorDim() == vcoeff.GetVDim(),
"vcoeff vdim != BdrVectorDim()");
ProjectBdrCoefficient(NULL, &vcoeff, attr);
}
void ParGridFunction::ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff,
const Array<int> &bdr_attr)
{
MFEM_VERIFY(VectorDim() == vcoeff.GetVDim(), "vcoeff vdim != VectorDim()");
MFEM_VERIFY(BdrVectorDim()+1 == vcoeff.GetVDim(),
"vcoeff vdim != BdrVectorDim()+1");
Array<int> values_counter;
AccumulateAndCountBdrTangentValues(vcoeff, bdr_attr, values_counter);
tests/unit/fem/test_var_order.cpp
+16 -2
View File
@@ -696,7 +696,14 @@ static void TestSolveVec(FiniteElementSpace &fespace)
GridFunction x(&fespace);
x = 0.0;
x.ProjectBdrCoefficient(exsol, ess_attr);
if (x.BdrVectorDim() == 1)
{
x.ProjectBdrCoefficientNormal(exsol, ess_attr);
}
else
{
x.ProjectBdrCoefficientTangent(exsol, ess_attr);
}
// Assemble the linear form
LinearForm lf(&fespace);
@@ -1081,7 +1088,14 @@ static void TestSolveParVec(ParFiniteElementSpace &fespace)
ParGridFunction x(&fespace);
x = 0.0;
x.ProjectBdrCoefficient(exsol, ess_attr);
if (x.BdrVectorDim() == 1)
{
x.ProjectBdrCoefficientNormal(exsol, ess_attr);
}
else
{
x.ProjectBdrCoefficientTangent(exsol, ess_attr);
}
// Assemble the linear form
ParLinearForm lf(&fespace);
tests/unit/mesh/test_psubmesh.cpp
+16 -4
View File
@@ -486,8 +486,14 @@ void multidomain_test_3d(FECType fec_type)
{
cylinder_gf.ProjectCoefficient(vcoeff);
outer_gf.ProjectCoefficient(vcoeff);
outer_gf.ProjectBdrCoefficient(vzerocoeff,
outer_cyl_surf_marker);
if (fec_type == FECType::RT)
{
outer_gf.ProjectBdrCoefficientNormal(vzerocoeff, outer_cyl_surf_marker);
}
else
{
outer_gf.ProjectBdrCoefficientTangent(vzerocoeff, outer_cyl_surf_marker);
}
outer_gf_ex.ProjectCoefficient(vcoeff);
}
ParSubMesh::Transfer(cylinder_gf, outer_gf);
@@ -507,8 +513,14 @@ void multidomain_test_3d(FECType fec_type)
{
outer_gf.ProjectCoefficient(vcoeff);
cylinder_gf.ProjectCoefficient(vcoeff);
cylinder_gf.ProjectBdrCoefficient(vzerocoeff,
cylinder_cyl_surf_marker);
if (fec_type == FECType::RT)
{
cylinder_gf.ProjectBdrCoefficientNormal(vzerocoeff, cylinder_cyl_surf_marker);
}
else
{
cylinder_gf.ProjectBdrCoefficientTangent(vzerocoeff, cylinder_cyl_surf_marker);
}
cylinder_gf_ex.ProjectCoefficient(vcoeff);
}
ParSubMesh::Transfer(outer_gf, cylinder_gf);