diff --git a/src/language/algorithms/ElasticityDiamondAlgorithm.cpp b/src/language/algorithms/ElasticityDiamondAlgorithm.cpp index 0dd1739e0cbcdcd6c85d3e1817fa5ed58de833c7..81fbcea6a803d13bcb85e9b016732a94de7b6895 100644 --- a/src/language/algorithms/ElasticityDiamondAlgorithm.cpp +++ b/src/language/algorithms/ElasticityDiamondAlgorithm.cpp @@ -68,10 +68,24 @@ ElasticityDiamondScheme<Dimension>::ElasticityDiamondScheme( switch (bc_descriptor->type()) { case IBoundaryConditionDescriptor::Type::symmetry: { - // const SymmetryBoundaryConditionDescriptor& sym_bc_descriptor = - // dynamic_cast<const SymmetryBoundaryConditionDescriptor&>(*bc_descriptor); - throw NotImplementedError("NIY"); - break; + const SymmetryBoundaryConditionDescriptor& sym_bc_descriptor = + dynamic_cast<const SymmetryBoundaryConditionDescriptor&>(*bc_descriptor); + if (sym_bc_descriptor.name() == "symmetry") { + 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()) { + if constexpr (Dimension > 1) { + Array<const FaceId> face_list = ref_face_list.list(); + boundary_condition_list.push_back(SymmetryBoundaryCondition{face_list}); + } else { + throw NotImplementedError("Symmetry conditions are not supported in 1d"); + break; + } + } + } + } } case IBoundaryConditionDescriptor::Type::dirichlet: { const DirichletBoundaryConditionDescriptor& dirichlet_bc_descriptor = @@ -195,6 +209,28 @@ ElasticityDiamondScheme<Dimension>::ElasticityDiamondScheme( return face_is_neumann; }(); + const FaceValue<const bool> primal_face_is_symmetry = [&] { + FaceValue<bool> face_is_symmetry{mesh->connectivity()}; + face_is_symmetry.fill(false); + for (const auto& boundary_condition : boundary_condition_list) { + std::visit( + [&](auto&& bc) { + using T = std::decay_t<decltype(bc)>; + if constexpr ((std::is_same_v<T, SymmetryBoundaryCondition>)) { + const auto& face_list = bc.faceList(); + + for (size_t i_face = 0; i_face < face_list.size(); ++i_face) { + const FaceId face_id = face_list[i_face]; + face_is_symmetry[face_id] = true; + } + } + }, + boundary_condition); + } + + return face_is_symmetry; + }(); + NodeValue<bool> primal_node_is_on_boundary(mesh->connectivity()); if (parallel::size() > 1) { throw NotImplementedError("Calculation of node_is_on_boundary is incorrect"); @@ -754,6 +790,7 @@ template <size_t Dimension> class ElasticityDiamondScheme<Dimension>::SymmetryBoundaryCondition { private: + const Array<const TinyVector<Dimension>> m_value_list; const Array<const FaceId> m_face_list; public: diff --git a/src/scheme/SymmetryBoundaryConditionDescriptor.hpp b/src/scheme/SymmetryBoundaryConditionDescriptor.hpp index f5657ad58856f38c5d645501a14b62e77039b4c1..f1b9fd81ea81b8e39efb2be3eaa33b258aacda2b 100644 --- a/src/scheme/SymmetryBoundaryConditionDescriptor.hpp +++ b/src/scheme/SymmetryBoundaryConditionDescriptor.hpp @@ -16,9 +16,16 @@ class SymmetryBoundaryConditionDescriptor : public IBoundaryConditionDescriptor return os; } + const std::string_view m_name; + std::shared_ptr<const IBoundaryDescriptor> m_boundary_descriptor; public: + std::string_view + name() const + { + return m_name; + } const IBoundaryDescriptor& boundaryDescriptor() const {