From d935de5d8618569416095c7b94913486435b7144 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Mon, 19 Jul 2021 16:15:14 +0200
Subject: [PATCH] Reorganize MeshNodeBoundary: objects can no longer be built
manually
One now must use getMeshNodeBoundary which checks if the boundary
exists.
This also concerns
- MeshFlatNodeBoundary (getMeshFlatNodeBoundary) and
- MeshLineNodeBoundary (getMeshLineNodeBoundary)
---
src/language/modules/SchemeModule.cpp | 2 +-
src/{scheme => mesh}/IBoundaryDescriptor.hpp | 0
src/mesh/MeshNodeBoundary.hpp | 122 ++++++++++++++++-
src/mesh/MeshRandomizer.hpp | 68 ++--------
src/scheme/AcousticSolver.cpp | 123 +++++++++---------
.../AxisBoundaryConditionDescriptor.hpp | 4 +-
.../DirichletBoundaryConditionDescriptor.hpp | 4 +-
.../FixedBoundaryConditionDescriptor.hpp | 4 +-
.../FourierBoundaryConditionDescriptor.hpp | 4 +-
.../FreeBoundaryConditionDescriptor.hpp | 4 +-
src/scheme/IBoundaryConditionDescriptor.hpp | 4 +
src/scheme/NamedBoundaryDescriptor.hpp | 2 +-
.../NeumannBoundaryConditionDescriptor.hpp | 4 +-
src/scheme/NumberedBoundaryDescriptor.hpp | 2 +-
.../SymmetryBoundaryConditionDescriptor.hpp | 4 +-
15 files changed, 213 insertions(+), 138 deletions(-)
rename src/{scheme => mesh}/IBoundaryDescriptor.hpp (100%)
diff --git a/src/language/modules/SchemeModule.cpp b/src/language/modules/SchemeModule.cpp
index 634ffe04f..81812defd 100644
--- a/src/language/modules/SchemeModule.cpp
+++ b/src/language/modules/SchemeModule.cpp
@@ -7,6 +7,7 @@
#include <language/utils/BuiltinFunctionEmbedder.hpp>
#include <language/utils/TypeDescriptor.hpp>
#include <mesh/Connectivity.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <mesh/Mesh.hpp>
#include <mesh/MeshData.hpp>
#include <mesh/MeshDataManager.hpp>
@@ -23,7 +24,6 @@
#include <scheme/FourierBoundaryConditionDescriptor.hpp>
#include <scheme/FreeBoundaryConditionDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <scheme/IDiscreteFunction.hpp>
#include <scheme/IDiscreteFunctionDescriptor.hpp>
#include <scheme/NamedBoundaryDescriptor.hpp>
diff --git a/src/scheme/IBoundaryDescriptor.hpp b/src/mesh/IBoundaryDescriptor.hpp
similarity index 100%
rename from src/scheme/IBoundaryDescriptor.hpp
rename to src/mesh/IBoundaryDescriptor.hpp
diff --git a/src/mesh/MeshNodeBoundary.hpp b/src/mesh/MeshNodeBoundary.hpp
index 0743d2d95..ab98ded0e 100644
--- a/src/mesh/MeshNodeBoundary.hpp
+++ b/src/mesh/MeshNodeBoundary.hpp
@@ -4,9 +4,11 @@
#include <Kokkos_Vector.hpp>
#include <algebra/TinyMatrix.hpp>
#include <algebra/TinyVector.hpp>
+#include <mesh/Connectivity.hpp>
#include <mesh/ConnectivityMatrix.hpp>
-#include <mesh/IConnectivity.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <mesh/ItemValue.hpp>
+#include <mesh/Mesh.hpp>
#include <mesh/RefItemList.hpp>
#include <utils/Array.hpp>
#include <utils/Exceptions.hpp>
@@ -23,6 +25,10 @@ class MeshNodeBoundary // clazy:exclude=copyable-polymorphic
std::array<TinyVector<Dimension>, Dimension*(Dimension - 1)> _getBounds(const MeshType& mesh) const;
public:
+ template <typename MeshType>
+ friend MeshNodeBoundary<MeshType::Dimension> getMeshNodeBoundary(const MeshType& mesh,
+ const IBoundaryDescriptor& boundary_descriptor);
+
MeshNodeBoundary& operator=(const MeshNodeBoundary&) = default;
MeshNodeBoundary& operator=(MeshNodeBoundary&&) = default;
@@ -32,6 +38,7 @@ class MeshNodeBoundary // clazy:exclude=copyable-polymorphic
return m_node_list;
}
+ protected:
template <typename MeshType>
MeshNodeBoundary(const MeshType& mesh, const RefFaceList& ref_face_list)
: m_boundary_name(ref_face_list.refId().tagName())
@@ -146,6 +153,7 @@ class MeshNodeBoundary // clazy:exclude=copyable-polymorphic
static_assert(Dimension == MeshType::Dimension);
}
+ public:
MeshNodeBoundary(const MeshNodeBoundary&) = default;
MeshNodeBoundary(MeshNodeBoundary&&) = default;
@@ -153,6 +161,41 @@ class MeshNodeBoundary // clazy:exclude=copyable-polymorphic
virtual ~MeshNodeBoundary() = default;
};
+template <typename MeshType>
+MeshNodeBoundary<MeshType::Dimension>
+getMeshNodeBoundary(const MeshType& mesh, const IBoundaryDescriptor& boundary_descriptor)
+{
+ for (size_t i_ref_node_list = 0; i_ref_node_list < mesh.connectivity().template numberOfRefItemList<ItemType::node>();
+ ++i_ref_node_list) {
+ const auto& ref_node_list = mesh.connectivity().template refItemList<ItemType::node>(i_ref_node_list);
+ const RefId& ref = ref_node_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshNodeBoundary<MeshType::Dimension>{mesh, ref_node_list};
+ }
+ }
+ for (size_t i_ref_edge_list = 0; i_ref_edge_list < mesh.connectivity().template numberOfRefItemList<ItemType::edge>();
+ ++i_ref_edge_list) {
+ const auto& ref_edge_list = mesh.connectivity().template refItemList<ItemType::edge>(i_ref_edge_list);
+ const RefId& ref = ref_edge_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshNodeBoundary<MeshType::Dimension>{mesh, ref_edge_list};
+ }
+ }
+ for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<ItemType::face>();
+ ++i_ref_face_list) {
+ const auto& ref_face_list = mesh.connectivity().template refItemList<ItemType::face>(i_ref_face_list);
+ const RefId& ref = ref_face_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshNodeBoundary<MeshType::Dimension>{mesh, ref_face_list};
+ }
+ }
+
+ std::ostringstream ost;
+ ost << "cannot find surface with name " << rang::fgB::red << boundary_descriptor << rang::style::reset;
+
+ throw NormalError(ost.str());
+}
+
template <>
template <typename MeshType>
std::array<TinyVector<2>, 2>
@@ -372,6 +415,12 @@ class MeshFlatNodeBoundary : public MeshNodeBoundary<Dimension> // clazy:exclu
MeshFlatNodeBoundary& operator=(const MeshFlatNodeBoundary&) = default;
MeshFlatNodeBoundary& operator=(MeshFlatNodeBoundary&&) = default;
+ template <typename MeshType>
+ friend MeshFlatNodeBoundary<MeshType::Dimension> getMeshFlatNodeBoundary(
+ const MeshType& mesh,
+ const IBoundaryDescriptor& boundary_descriptor);
+
+ private:
template <typename MeshType>
MeshFlatNodeBoundary(const MeshType& mesh, const RefFaceList& ref_face_list)
: MeshNodeBoundary<Dimension>(mesh, ref_face_list), m_outgoing_normal(_getOutgoingNormal(mesh))
@@ -386,12 +435,40 @@ class MeshFlatNodeBoundary : public MeshNodeBoundary<Dimension> // clazy:exclu
;
}
+ public:
MeshFlatNodeBoundary() = default;
MeshFlatNodeBoundary(const MeshFlatNodeBoundary&) = default;
MeshFlatNodeBoundary(MeshFlatNodeBoundary&&) = default;
virtual ~MeshFlatNodeBoundary() = default;
};
+template <typename MeshType>
+MeshFlatNodeBoundary<MeshType::Dimension>
+getMeshFlatNodeBoundary(const MeshType& mesh, const IBoundaryDescriptor& boundary_descriptor)
+{
+ for (size_t i_ref_node_list = 0; i_ref_node_list < mesh.connectivity().template numberOfRefItemList<ItemType::node>();
+ ++i_ref_node_list) {
+ const auto& ref_node_list = mesh.connectivity().template refItemList<ItemType::node>(i_ref_node_list);
+ const RefId& ref = ref_node_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshFlatNodeBoundary<MeshType::Dimension>{mesh, ref_node_list};
+ }
+ }
+ for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<ItemType::face>();
+ ++i_ref_face_list) {
+ const auto& ref_face_list = mesh.connectivity().template refItemList<ItemType::face>(i_ref_face_list);
+ const RefId& ref = ref_face_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshFlatNodeBoundary<MeshType::Dimension>{mesh, ref_face_list};
+ }
+ }
+
+ std::ostringstream ost;
+ ost << "cannot find surface with name " << rang::fgB::red << boundary_descriptor << rang::style::reset;
+
+ throw NormalError(ost.str());
+}
+
template <>
template <typename MeshType>
PUGS_INLINE TinyVector<1, double>
@@ -687,6 +764,12 @@ class MeshLineNodeBoundary : public MeshNodeBoundary<Dimension> // clazy:exclu
}
public:
+ template <typename MeshType>
+ friend MeshLineNodeBoundary<MeshType::Dimension> getMeshLineNodeBoundary(
+ const MeshType& mesh,
+ const IBoundaryDescriptor& boundary_descriptor);
+
+ PUGS_INLINE
const Rd&
direction() const
{
@@ -696,6 +779,7 @@ class MeshLineNodeBoundary : public MeshNodeBoundary<Dimension> // clazy:exclu
MeshLineNodeBoundary& operator=(const MeshLineNodeBoundary&) = default;
MeshLineNodeBoundary& operator=(MeshLineNodeBoundary&&) = default;
+ private:
template <typename MeshType>
MeshLineNodeBoundary(const MeshType& mesh, const RefEdgeList& ref_edge_list)
: MeshNodeBoundary<Dimension>(mesh, ref_edge_list), m_direction(_getDirection(mesh))
@@ -717,12 +801,48 @@ class MeshLineNodeBoundary : public MeshNodeBoundary<Dimension> // clazy:exclu
;
}
+ public:
MeshLineNodeBoundary() = default;
MeshLineNodeBoundary(const MeshLineNodeBoundary&) = default;
MeshLineNodeBoundary(MeshLineNodeBoundary&&) = default;
virtual ~MeshLineNodeBoundary() = default;
};
+template <typename MeshType>
+MeshLineNodeBoundary<MeshType::Dimension>
+getMeshLineNodeBoundary(const MeshType& mesh, const IBoundaryDescriptor& boundary_descriptor)
+{
+ for (size_t i_ref_node_list = 0; i_ref_node_list < mesh.connectivity().template numberOfRefItemList<ItemType::node>();
+ ++i_ref_node_list) {
+ const auto& ref_node_list = mesh.connectivity().template refItemList<ItemType::node>(i_ref_node_list);
+ const RefId& ref = ref_node_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshLineNodeBoundary<MeshType::Dimension>{mesh, ref_node_list};
+ }
+ }
+ for (size_t i_ref_edge_list = 0; i_ref_edge_list < mesh.connectivity().template numberOfRefItemList<ItemType::edge>();
+ ++i_ref_edge_list) {
+ const auto& ref_edge_list = mesh.connectivity().template refItemList<ItemType::edge>(i_ref_edge_list);
+ const RefId& ref = ref_edge_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshLineNodeBoundary<MeshType::Dimension>{mesh, ref_edge_list};
+ }
+ }
+ for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<ItemType::face>();
+ ++i_ref_face_list) {
+ const auto& ref_face_list = mesh.connectivity().template refItemList<ItemType::face>(i_ref_face_list);
+ const RefId& ref = ref_face_list.refId();
+ if (ref == boundary_descriptor) {
+ return MeshLineNodeBoundary<MeshType::Dimension>{mesh, ref_face_list};
+ }
+ }
+
+ std::ostringstream ost;
+ ost << "cannot find surface with name " << rang::fgB::red << boundary_descriptor << rang::style::reset;
+
+ throw NormalError(ost.str());
+}
+
template <>
template <typename MeshType>
PUGS_INLINE TinyVector<1>
diff --git a/src/mesh/MeshRandomizer.hpp b/src/mesh/MeshRandomizer.hpp
index b6bfe625b..aa2258c59 100644
--- a/src/mesh/MeshRandomizer.hpp
+++ b/src/mesh/MeshRandomizer.hpp
@@ -43,82 +43,32 @@ class MeshRandomizer
{
BoundaryConditionList bc_list;
- constexpr ItemType EdgeType = [] {
- if constexpr (Dimension == 3) {
- return ItemType::edge;
- } else if constexpr (Dimension == 2) {
- return ItemType::face;
- } else {
- return ItemType::node;
- }
- }();
-
- constexpr ItemType FaceType = [] {
- if constexpr (Dimension > 1) {
- return ItemType::face;
- } else {
- return ItemType::node;
- }
- }();
-
for (const auto& bc_descriptor : bc_descriptor_list) {
- bool is_valid_boundary_condition = true;
-
switch (bc_descriptor->type()) {
case IBoundaryConditionDescriptor::Type::axis: {
- const AxisBoundaryConditionDescriptor& axis_bc_descriptor =
- dynamic_cast<const AxisBoundaryConditionDescriptor&>(*bc_descriptor);
- for (size_t i_ref_edge_list = 0; i_ref_edge_list < mesh.connectivity().template numberOfRefItemList<EdgeType>();
- ++i_ref_edge_list) {
- const auto& ref_edge_list = mesh.connectivity().template refItemList<EdgeType>(i_ref_edge_list);
- const RefId& ref = ref_edge_list.refId();
- if (ref == axis_bc_descriptor.boundaryDescriptor()) {
- if constexpr (Dimension == 1) {
- bc_list.emplace_back(FixedBoundaryCondition{MeshNodeBoundary<Dimension>{mesh, ref_edge_list}});
- } else {
- bc_list.emplace_back(AxisBoundaryCondition{MeshLineNodeBoundary<Dimension>(mesh, ref_edge_list)});
- }
- }
+ if constexpr (Dimension == 1) {
+ bc_list.emplace_back(FixedBoundaryCondition{getMeshNodeBoundary(mesh, bc_descriptor->boundaryDescriptor())});
+ } else {
+ bc_list.emplace_back(
+ AxisBoundaryCondition{getMeshLineNodeBoundary(mesh, bc_descriptor->boundaryDescriptor())});
}
- is_valid_boundary_condition = true;
break;
}
case IBoundaryConditionDescriptor::Type::symmetry: {
- const SymmetryBoundaryConditionDescriptor& sym_bc_descriptor =
- dynamic_cast<const SymmetryBoundaryConditionDescriptor&>(*bc_descriptor);
- for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<FaceType>();
- ++i_ref_face_list) {
- const auto& ref_face_list = mesh.connectivity().template refItemList<FaceType>(i_ref_face_list);
- const RefId& ref = ref_face_list.refId();
- if (ref == sym_bc_descriptor.boundaryDescriptor()) {
- bc_list.emplace_back(SymmetryBoundaryCondition{MeshFlatNodeBoundary<Dimension>(mesh, ref_face_list)});
- }
- }
- is_valid_boundary_condition = true;
+ bc_list.emplace_back(
+ SymmetryBoundaryCondition{getMeshFlatNodeBoundary(mesh, bc_descriptor->boundaryDescriptor())});
break;
}
case IBoundaryConditionDescriptor::Type::fixed: {
- const FixedBoundaryConditionDescriptor& fixed_bc_descriptor =
- dynamic_cast<const FixedBoundaryConditionDescriptor&>(*bc_descriptor);
- for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<FaceType>();
- ++i_ref_face_list) {
- const auto& ref_face_list = mesh.connectivity().template refItemList<FaceType>(i_ref_face_list);
- const RefId& ref = ref_face_list.refId();
- if (ref == fixed_bc_descriptor.boundaryDescriptor()) {
- bc_list.emplace_back(FixedBoundaryCondition{MeshNodeBoundary<Dimension>{mesh, ref_face_list}});
- }
- }
+ bc_list.emplace_back(FixedBoundaryCondition{getMeshNodeBoundary(mesh, bc_descriptor->boundaryDescriptor())});
break;
}
default: {
- is_valid_boundary_condition = false;
- }
- }
- if (not is_valid_boundary_condition) {
std::ostringstream error_msg;
error_msg << *bc_descriptor << " is an invalid boundary condition for mesh randomizer";
throw NormalError(error_msg.str());
}
+ }
}
return bc_list;
diff --git a/src/scheme/AcousticSolver.cpp b/src/scheme/AcousticSolver.cpp
index 0ad9b7901..7bd97aef4 100644
--- a/src/scheme/AcousticSolver.cpp
+++ b/src/scheme/AcousticSolver.cpp
@@ -6,6 +6,7 @@
#include <scheme/DirichletBoundaryConditionDescriptor.hpp>
#include <scheme/DiscreteFunctionP0.hpp>
#include <scheme/DiscreteFunctionUtils.hpp>
+#include <scheme/FixedBoundaryConditionDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
#include <scheme/IDiscreteFunction.hpp>
#include <scheme/IDiscreteFunctionDescriptor.hpp>
@@ -68,12 +69,13 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
using DiscreteScalarFunction = DiscreteFunctionP0<Dimension, double>;
using DiscreteVectorFunction = DiscreteFunctionP0<Dimension, Rd>;
+ class FixedBoundaryCondition;
class PressureBoundaryCondition;
class SymmetryBoundaryCondition;
class VelocityBoundaryCondition;
- using BoundaryCondition =
- std::variant<PressureBoundaryCondition, SymmetryBoundaryCondition, VelocityBoundaryCondition>;
+ using BoundaryCondition = std::
+ variant<FixedBoundaryCondition, PressureBoundaryCondition, SymmetryBoundaryCondition, VelocityBoundaryCondition>;
using BoundaryConditionList = std::vector<BoundaryCondition>;
@@ -83,22 +85,8 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
NodeValue<const Rd> m_ur;
NodeValuePerCell<const Rd> m_Fjr;
- CellValue<const double>
- _getRhoC(const DiscreteScalarFunction& rho, const DiscreteScalarFunction& c)
- {
- Assert(rho.mesh() == c.mesh());
-
- const MeshType& mesh = dynamic_cast<const MeshType&>(*rho.mesh());
-
- CellValue<double> rhoc{mesh.connectivity()};
- parallel_for(
- mesh.numberOfCells(), PUGS_LAMBDA(CellId cell_id) { rhoc[cell_id] = rho[cell_id] * c[cell_id]; });
-
- return rhoc;
- }
-
NodeValuePerCell<const Rdxd>
- _computeGlaceAjr(const MeshType& mesh, const CellValue<const double>& rhoc)
+ _computeGlaceAjr(const MeshType& mesh, const DiscreteScalarFunction& rhoc)
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(mesh);
@@ -120,7 +108,7 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
}
NodeValuePerCell<const Rdxd>
- _computeEucclhydAjr(const MeshType& mesh, const CellValue<const double>& rhoc)
+ _computeEucclhydAjr(const MeshType& mesh, const DiscreteScalarFunction& rhoc)
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(mesh);
@@ -168,7 +156,7 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
}
NodeValuePerCell<const Rdxd>
- _computeAjr(const MeshType& mesh, const CellValue<const double>& rhoc)
+ _computeAjr(const MeshType& mesh, const DiscreteScalarFunction& rhoc)
{
if constexpr (Dimension == 1) {
return _computeGlaceAjr(mesh, rhoc);
@@ -264,42 +252,26 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
switch (bc_descriptor->type()) {
case IBoundaryConditionDescriptor::Type::symmetry: {
- const SymmetryBoundaryConditionDescriptor& sym_bc_descriptor =
- dynamic_cast<const SymmetryBoundaryConditionDescriptor&>(*bc_descriptor);
- for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<FaceType>();
- ++i_ref_face_list) {
- const auto& ref_face_list = mesh.connectivity().template refItemList<FaceType>(i_ref_face_list);
- const RefId& ref = ref_face_list.refId();
- if (ref == sym_bc_descriptor.boundaryDescriptor()) {
- SymmetryBoundaryCondition{MeshFlatNodeBoundary<Dimension>(mesh, ref_face_list)};
- bc_list.push_back(SymmetryBoundaryCondition{MeshFlatNodeBoundary<Dimension>(mesh, ref_face_list)});
- }
- }
- is_valid_boundary_condition = true;
+ bc_list.push_back(
+ SymmetryBoundaryCondition(getMeshFlatNodeBoundary(mesh, bc_descriptor->boundaryDescriptor())));
+ break;
+ }
+ case IBoundaryConditionDescriptor::Type::fixed: {
+ bc_list.push_back(FixedBoundaryCondition(getMeshNodeBoundary(mesh, bc_descriptor->boundaryDescriptor())));
break;
}
-
case IBoundaryConditionDescriptor::Type::dirichlet: {
const DirichletBoundaryConditionDescriptor& dirichlet_bc_descriptor =
dynamic_cast<const DirichletBoundaryConditionDescriptor&>(*bc_descriptor);
if (dirichlet_bc_descriptor.name() == "velocity") {
- for (size_t i_ref_face_list = 0;
- i_ref_face_list < mesh.connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
- const auto& ref_face_list = mesh.connectivity().template refItemList<FaceType>(i_ref_face_list);
- const RefId& ref = ref_face_list.refId();
- if (ref == dirichlet_bc_descriptor.boundaryDescriptor()) {
- MeshNodeBoundary<Dimension> mesh_node_boundary{mesh, ref_face_list};
-
- const FunctionSymbolId velocity_id = dirichlet_bc_descriptor.rhsSymbolId();
-
- const auto& node_list = mesh_node_boundary.nodeList();
+ MeshNodeBoundary<Dimension> mesh_node_boundary =
+ getMeshNodeBoundary(mesh, dirichlet_bc_descriptor.boundaryDescriptor());
- Array<const TinyVector<Dimension>> value_list = InterpolateItemValue<TinyVector<Dimension>(
- TinyVector<Dimension>)>::template interpolate<ItemType::node>(velocity_id, mesh.xr(), node_list);
+ Array<const TinyVector<Dimension>> value_list = InterpolateItemValue<TinyVector<Dimension>(
+ TinyVector<Dimension>)>::template interpolate<ItemType::node>(dirichlet_bc_descriptor.rhsSymbolId(),
+ mesh.xr(), mesh_node_boundary.nodeList());
- bc_list.push_back(VelocityBoundaryCondition{node_list, value_list});
- }
- }
+ bc_list.push_back(VelocityBoundaryCondition{mesh_node_boundary, value_list});
} else if (dirichlet_bc_descriptor.name() == "pressure") {
for (size_t i_ref_face_list = 0;
i_ref_face_list < mesh.connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
@@ -401,8 +373,7 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list)
: m_solver_type{solver_type}, m_boundary_condition_list{this->_getBCList(mesh, bc_descriptor_list)}
{
- CellValue<const double> rhoc = this->_getRhoC(rho, c);
- NodeValuePerCell<const Rdxd> Ajr = this->_computeAjr(mesh, rhoc);
+ NodeValuePerCell<const Rdxd> Ajr = this->_computeAjr(mesh, rho * c);
NodeValue<Rdxd> Ar = this->_computeAr(mesh, Ajr);
NodeValue<Rd> br = this->_computeBr(mesh, Ajr, u, p);
@@ -485,10 +456,9 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
throw NormalError("acoustic solver expects P0 functions");
}
- return this->apply(dt, dynamic_cast<const MeshType&>(*i_mesh),
- *std::dynamic_pointer_cast<const DiscreteScalarFunction>(rho),
- *std::dynamic_pointer_cast<const DiscreteVectorFunction>(u),
- *std::dynamic_pointer_cast<const DiscreteScalarFunction>(E));
+ return this->apply(dt, dynamic_cast<const MeshType&>(*i_mesh), dynamic_cast<const DiscreteScalarFunction&>(*rho),
+ dynamic_cast<const DiscreteVectorFunction&>(*u),
+ dynamic_cast<const DiscreteScalarFunction&>(*E));
}
AcousticSolver(SolverType solver_type,
@@ -500,10 +470,10 @@ class AcousticSolverHandler::AcousticSolver final : public AcousticSolverHandler
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list)
: AcousticSolver{solver_type,
dynamic_cast<const MeshType&>(*mesh),
- *std::dynamic_pointer_cast<const DiscreteScalarFunction>(rho),
- *std::dynamic_pointer_cast<const DiscreteScalarFunction>(c),
- *std::dynamic_pointer_cast<const DiscreteVectorFunction>(u),
- *std::dynamic_pointer_cast<const DiscreteScalarFunction>(p),
+ dynamic_cast<const DiscreteScalarFunction&>(*rho),
+ dynamic_cast<const DiscreteScalarFunction&>(*c),
+ dynamic_cast<const DiscreteVectorFunction&>(*u),
+ dynamic_cast<const DiscreteScalarFunction&>(*p),
bc_descriptor_list}
{}
@@ -624,12 +594,41 @@ AcousticSolverHandler::AcousticSolver<Dimension>::_applyVelocityBC(NodeValue<Rdx
Ar[node_id] = identity;
br[node_id] = value;
});
+ } else if constexpr (std::is_same_v<FixedBoundaryCondition, T>) {
+ const auto& node_list = bc.nodeList();
+ parallel_for(
+ node_list.size(), PUGS_LAMBDA(size_t i_node) {
+ NodeId node_id = node_list[i_node];
+
+ Ar[node_id] = identity;
+ br[node_id] = zero;
+ });
}
},
boundary_condition);
}
}
+template <size_t Dimension>
+class AcousticSolverHandler::AcousticSolver<Dimension>::FixedBoundaryCondition
+{
+ private:
+ const MeshNodeBoundary<Dimension> m_mesh_node_boundary;
+
+ public:
+ const Array<const NodeId>&
+ nodeList() const
+ {
+ return m_mesh_node_boundary.nodeList();
+ }
+
+ FixedBoundaryCondition(const MeshNodeBoundary<Dimension> mesh_node_boundary)
+ : m_mesh_node_boundary{mesh_node_boundary}
+ {}
+
+ ~FixedBoundaryCondition() = default;
+};
+
template <size_t Dimension>
class AcousticSolverHandler::AcousticSolver<Dimension>::PressureBoundaryCondition
{
@@ -688,14 +687,15 @@ template <size_t Dimension>
class AcousticSolverHandler::AcousticSolver<Dimension>::VelocityBoundaryCondition
{
private:
+ const MeshNodeBoundary<Dimension> m_mesh_node_boundary;
+
const Array<const TinyVector<Dimension>> m_value_list;
- const Array<const NodeId> m_node_list;
public:
const Array<const NodeId>&
nodeList() const
{
- return m_node_list;
+ return m_mesh_node_boundary.nodeList();
}
const Array<const TinyVector<Dimension>>&
@@ -704,8 +704,9 @@ class AcousticSolverHandler::AcousticSolver<Dimension>::VelocityBoundaryConditio
return m_value_list;
}
- VelocityBoundaryCondition(const Array<const NodeId>& node_list, const Array<const TinyVector<Dimension>>& value_list)
- : m_value_list{value_list}, m_node_list{node_list}
+ VelocityBoundaryCondition(const MeshNodeBoundary<Dimension>& mesh_node_boundary,
+ const Array<const TinyVector<Dimension>>& value_list)
+ : m_mesh_node_boundary{mesh_node_boundary}, m_value_list{value_list}
{}
~VelocityBoundaryCondition() = default;
diff --git a/src/scheme/AxisBoundaryConditionDescriptor.hpp b/src/scheme/AxisBoundaryConditionDescriptor.hpp
index a03a4bcdd..07ab3be48 100644
--- a/src/scheme/AxisBoundaryConditionDescriptor.hpp
+++ b/src/scheme/AxisBoundaryConditionDescriptor.hpp
@@ -1,8 +1,8 @@
#ifndef AXIS_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#define AXIS_BOUNDARY_CONDITION_DESCRIPTOR_HPP
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -20,7 +20,7 @@ class AxisBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
public:
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
diff --git a/src/scheme/DirichletBoundaryConditionDescriptor.hpp b/src/scheme/DirichletBoundaryConditionDescriptor.hpp
index 5fa6c42f9..36f95df0d 100644
--- a/src/scheme/DirichletBoundaryConditionDescriptor.hpp
+++ b/src/scheme/DirichletBoundaryConditionDescriptor.hpp
@@ -2,8 +2,8 @@
#define DIRICHLET_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#include <language/utils/FunctionSymbolId.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -36,7 +36,7 @@ class DirichletBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
}
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
diff --git a/src/scheme/FixedBoundaryConditionDescriptor.hpp b/src/scheme/FixedBoundaryConditionDescriptor.hpp
index 9167c24be..b69648e38 100644
--- a/src/scheme/FixedBoundaryConditionDescriptor.hpp
+++ b/src/scheme/FixedBoundaryConditionDescriptor.hpp
@@ -2,8 +2,8 @@
#define FIXED_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#include <language/utils/FunctionSymbolId.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -23,7 +23,7 @@ class FixedBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
public:
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
diff --git a/src/scheme/FourierBoundaryConditionDescriptor.hpp b/src/scheme/FourierBoundaryConditionDescriptor.hpp
index e7aeb830c..c9d75883d 100644
--- a/src/scheme/FourierBoundaryConditionDescriptor.hpp
+++ b/src/scheme/FourierBoundaryConditionDescriptor.hpp
@@ -2,8 +2,8 @@
#define FOURIER_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#include <language/utils/FunctionSymbolId.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -43,7 +43,7 @@ class FourierBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
}
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
diff --git a/src/scheme/FreeBoundaryConditionDescriptor.hpp b/src/scheme/FreeBoundaryConditionDescriptor.hpp
index b0eb2ceb2..2d46f2a36 100644
--- a/src/scheme/FreeBoundaryConditionDescriptor.hpp
+++ b/src/scheme/FreeBoundaryConditionDescriptor.hpp
@@ -1,8 +1,8 @@
#ifndef FREE_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#define FREE_BOUNDARY_CONDITION_DESCRIPTOR_HPP
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -20,7 +20,7 @@ class FreeBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
public:
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
diff --git a/src/scheme/IBoundaryConditionDescriptor.hpp b/src/scheme/IBoundaryConditionDescriptor.hpp
index bf92cd450..c0ea161db 100644
--- a/src/scheme/IBoundaryConditionDescriptor.hpp
+++ b/src/scheme/IBoundaryConditionDescriptor.hpp
@@ -3,6 +3,8 @@
#include <iostream>
+class IBoundaryDescriptor;
+
class IBoundaryConditionDescriptor
{
public:
@@ -27,6 +29,8 @@ class IBoundaryConditionDescriptor
return bcd._write(os);
}
+ virtual const IBoundaryDescriptor& boundaryDescriptor() const = 0;
+
virtual Type type() const = 0;
IBoundaryConditionDescriptor() = default;
diff --git a/src/scheme/NamedBoundaryDescriptor.hpp b/src/scheme/NamedBoundaryDescriptor.hpp
index db88ee4a0..5d8959261 100644
--- a/src/scheme/NamedBoundaryDescriptor.hpp
+++ b/src/scheme/NamedBoundaryDescriptor.hpp
@@ -1,7 +1,7 @@
#ifndef NAMED_BOUNDARY_DESCRIPTOR_HPP
#define NAMED_BOUNDARY_DESCRIPTOR_HPP
-#include <scheme/IBoundaryDescriptor.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <iostream>
#include <string>
diff --git a/src/scheme/NeumannBoundaryConditionDescriptor.hpp b/src/scheme/NeumannBoundaryConditionDescriptor.hpp
index edc1b1fb3..597ab95e7 100644
--- a/src/scheme/NeumannBoundaryConditionDescriptor.hpp
+++ b/src/scheme/NeumannBoundaryConditionDescriptor.hpp
@@ -2,8 +2,8 @@
#define NEUMANN_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#include <language/utils/FunctionSymbolId.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -36,7 +36,7 @@ class NeumannBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
}
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
diff --git a/src/scheme/NumberedBoundaryDescriptor.hpp b/src/scheme/NumberedBoundaryDescriptor.hpp
index 4ed9a2628..47aa7fa53 100644
--- a/src/scheme/NumberedBoundaryDescriptor.hpp
+++ b/src/scheme/NumberedBoundaryDescriptor.hpp
@@ -1,7 +1,7 @@
#ifndef NUMBERED_BOUNDARY_DESCRIPTOR_HPP
#define NUMBERED_BOUNDARY_DESCRIPTOR_HPP
-#include <scheme/IBoundaryDescriptor.hpp>
+#include <mesh/IBoundaryDescriptor.hpp>
#include <iostream>
diff --git a/src/scheme/SymmetryBoundaryConditionDescriptor.hpp b/src/scheme/SymmetryBoundaryConditionDescriptor.hpp
index f5657ad58..9364090dd 100644
--- a/src/scheme/SymmetryBoundaryConditionDescriptor.hpp
+++ b/src/scheme/SymmetryBoundaryConditionDescriptor.hpp
@@ -1,8 +1,8 @@
#ifndef SYMMETRY_BOUNDARY_CONDITION_DESCRIPTOR_HPP
#define SYMMETRY_BOUNDARY_CONDITION_DESCRIPTOR_HPP
+#include <mesh/IBoundaryDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
-#include <scheme/IBoundaryDescriptor.hpp>
#include <memory>
@@ -20,7 +20,7 @@ class SymmetryBoundaryConditionDescriptor : public IBoundaryConditionDescriptor
public:
const IBoundaryDescriptor&
- boundaryDescriptor() const
+ boundaryDescriptor() const final
{
return *m_boundary_descriptor;
}
--
GitLab