diff --git a/src/language/modules/SchemeModule.cpp b/src/language/modules/SchemeModule.cpp index 62744a6d96802f21fc25dd793a9f3f4efcc7ba87..fafb9a4b3c2314f26bd9bb1f3a423efd34cb31b4 100644 --- a/src/language/modules/SchemeModule.cpp +++ b/src/language/modules/SchemeModule.cpp @@ -1531,7 +1531,7 @@ SchemeModule::SchemeModule() )); - this->_addBuiltinFunction("hyperplastic_glace_solver_serraille", + this->_addBuiltinFunction("apply_hyperplastic_fluxes_cauchygreen_serraille", std::function( [](const std::shared_ptr<const DiscreteFunctionVariant>& rho, @@ -1540,9 +1540,9 @@ SchemeModule::SchemeModule() const std::shared_ptr<const DiscreteFunctionVariant>& Fe, const std::shared_ptr<const DiscreteFunctionVariant>& zeta, const std::shared_ptr<const DiscreteFunctionVariant>& yield, - const std::shared_ptr<const DiscreteFunctionVariant>& aL, - const std::shared_ptr<const DiscreteFunctionVariant>& aT, const std::shared_ptr<const DiscreteFunctionVariant>& sigma, + const std::shared_ptr<const ItemValueVariant>& ur, + const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr, const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>, @@ -1550,15 +1550,40 @@ SchemeModule::SchemeModule() std::shared_ptr<const DiscreteFunctionVariant>, std::shared_ptr<const DiscreteFunctionVariant>, std::shared_ptr<const DiscreteFunctionVariant>> { - return HyperplasticSolverSerrailleHandler{ - getCommonMesh({rho, u, E, Fe, zeta, yield, aL, aT, sigma})} - .solver() - .apply(HyperplasticSolverSerrailleHandler::SolverType::Glace, - HyperplasticSolverSerrailleHandler::RelaxationType::CauchyGreen, dt, rho, u, E, Fe, - zeta, yield, aL, aT, sigma, bc_descriptor_list); - } + return HyperplasticSolverSerrailleHandler{getCommonMesh({rho, u, E, Fe})} .solver() + .apply_fluxes(dt, rho, u, E, Fe, zeta, yield, sigma, ur, Fjr, bc_descriptor_list, + HyperplasticSolverSerrailleHandler::RelaxationType::CauchyGreen + ); + })); - )); + this->_addBuiltinFunction("hyperplastic_glace_solver_serraille", + std::function( + + [](const std::shared_ptr<const DiscreteFunctionVariant>& rho, + const std::shared_ptr<const DiscreteFunctionVariant>& u, + const std::shared_ptr<const DiscreteFunctionVariant>& E, + const std::shared_ptr<const DiscreteFunctionVariant>& Fe, + const std::shared_ptr<const DiscreteFunctionVariant>& zeta, + const std::shared_ptr<const DiscreteFunctionVariant>& yield, + const std::shared_ptr<const DiscreteFunctionVariant>& aL, + const std::shared_ptr<const DiscreteFunctionVariant>& aT, + const std::shared_ptr<const DiscreteFunctionVariant>& sigma, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& + bc_descriptor_list, + const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>, + std::shared_ptr<const DiscreteFunctionVariant>, + std::shared_ptr<const DiscreteFunctionVariant>, + std::shared_ptr<const DiscreteFunctionVariant>, + std::shared_ptr<const DiscreteFunctionVariant>> { + return HyperplasticSolverSerrailleHandler{ + getCommonMesh({rho, u, E, Fe, zeta, yield, aL, aT, sigma})} + .solver() + .apply(HyperplasticSolverSerrailleHandler::SolverType::Glace, + HyperplasticSolverSerrailleHandler::RelaxationType::CauchyGreen, dt, rho, u, E, Fe, + zeta, yield, aL, aT, sigma, bc_descriptor_list); + } + + )); diff --git a/src/scheme/SerrailleSolver.cpp b/src/scheme/SerrailleSolver.cpp index 717a4f093f58857beda9428e664011f0d0ac67b6..27912a1ea5cf35cc17ea69b6473e5bb81c0629ab 100644 --- a/src/scheme/SerrailleSolver.cpp +++ b/src/scheme/SerrailleSolver.cpp @@ -434,6 +434,29 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final this->_applySymmetryBC(bc_list, Ar, br); this->_applyVelocityBC(bc_list, Ar, br); } + void + _projectOnBoundary(const BoundaryConditionList& bc_list, NodeValue<Rd>& xr) const + { + for (const auto& boundary_condition : bc_list) { + std::visit( + [&](auto&& bc) { + using T = std::decay_t<decltype(bc)>; + if constexpr (std::is_same_v<SymmetryBoundaryCondition, T>) { + const Rd& n = bc.outgoingNormal(); + const Rd& o = bc.origin(); + + const Array<const NodeId>& node_list = bc.nodeList(); + parallel_for( + bc.numberOfNodes(), PUGS_LAMBDA(int r_number) { + const NodeId r = node_list[r_number]; + Rd& x = xr[node_list[r]]; + x -= dot(x - o, n) * n; + }); + } + }, + boundary_condition); + } + } NodeValue<const Rd> _computeUr(const MeshType& mesh, const NodeValue<Rdxd>& Ar, const NodeValue<Rd>& br) const @@ -536,7 +559,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final const DiscreteScalarFunction& yield, const DiscreteTensorFunction3d& sigma, const NodeValue<const Rd>& ur, - const NodeValuePerCell<const Rd>& Fjr) const + const NodeValuePerCell<const Rd>& Fjr, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list + ) const { const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix(); @@ -623,7 +648,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final const DiscreteScalarFunction& yield, const DiscreteTensorFunction3d& sigma, const NodeValue<const Rd>& ur, - const NodeValuePerCell<const Rd>& Fjr) const + const NodeValuePerCell<const Rd>& Fjr, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list + ) const { const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix(); @@ -711,7 +738,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final const DiscreteScalarFunction& yield, const DiscreteTensorFunction3d& sigma, const NodeValue<const Rd>& ur, - const NodeValuePerCell<const Rd>& Fjr) const + const NodeValuePerCell<const Rd>& Fjr, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list + ) const { const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix(); @@ -723,6 +752,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final NodeValue<Rd> new_xr = copy(mesh.xr()); parallel_for(mesh.numberOfNodes(), PUGS_LAMBDA(NodeId r) { new_xr[r] += dt * ur[r]; }); + const BoundaryConditionList bc_list = this->_getBCList(mesh, bc_descriptor_list); + this->_projectOnBoundary(bc_list, new_xr); + std::shared_ptr<const MeshType> new_mesh = std::make_shared<MeshType>(mesh.shared_connectivity(), new_xr); CellValue<const double> Vj = MeshDataManager::instance().getMeshData(mesh).Vj(); @@ -778,7 +810,6 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final // const R3x3 expMT = EigenvalueSolver{}.computeExpForTinyMatrix(transpose(M)); // calcul de S trial - const R3x3 Strial = (new_mu[j] / J) * (Be_barre - 1./3*trace(Be_barre)*I); double limit2 = 2. / 3 * yield[j] * yield[j]; double normS2 = frobeniusNorm(Strial); @@ -800,7 +831,6 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final // const double fenorm1 = frobeniusNorm(new_Be[j]); // rationorm[j] = fenorm1 / fenorm0; }); - std::cout << "sum_chi " << sum(chi) << "\n"; CellValue<const double> new_Vj = MeshDataManager::instance().getMeshData(*new_mesh).Vj(); @@ -828,7 +858,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final const std::shared_ptr<const DiscreteFunctionVariant>& sigma, const std::shared_ptr<const ItemValueVariant>& ur, const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr, - const RelaxationType& relaxation_type) const + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, + const RelaxationType& relaxation_type + ) const { std::shared_ptr mesh_v = getCommonMesh({rho, u, E}); if (not mesh_v) { @@ -850,7 +882,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final yield->get<DiscreteScalarFunction>(), // sigma->get<DiscreteTensorFunction3d>(), // ur->get<NodeValue<const Rd>>(), // - Fjr->get<NodeValuePerCell<const Rd>>()); + Fjr->get<NodeValuePerCell<const Rd>>(), + bc_descriptor_list + ); } case RelaxationType::Implicit: { return this->apply_fluxes2(dt, // @@ -863,7 +897,9 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final yield->get<DiscreteScalarFunction>(), // sigma->get<DiscreteTensorFunction3d>(), // ur->get<NodeValue<const Rd>>(), // - Fjr->get<NodeValuePerCell<const Rd>>()); + Fjr->get<NodeValuePerCell<const Rd>>(), + bc_descriptor_list + ); } case RelaxationType::CauchyGreen: { return this->apply_fluxes_B(dt, // @@ -876,7 +912,10 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final yield->get<DiscreteScalarFunction>(), // sigma->get<DiscreteTensorFunction3d>(), // ur->get<NodeValue<const Rd>>(), // - Fjr->get<NodeValuePerCell<const Rd>>()); + Fjr->get<NodeValuePerCell<const Rd>>(), + bc_descriptor_list + ); + } default: { throw UnexpectedError("invalid relaxation type"); @@ -1088,7 +1127,7 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver final const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) const { auto [ur, Fjr] = compute_fluxes(solver_type, rho, aL, aT, u, sigma, bc_descriptor_list); - return apply_fluxes(dt, rho, u, E, Fe, zeta, yield, sigma, ur, Fjr, relaxation_type); + return apply_fluxes(dt, rho, u, E, Fe, zeta, yield, sigma, ur, Fjr, bc_descriptor_list, relaxation_type); } std::tuple<std::shared_ptr<const MeshVariant>, @@ -1485,6 +1524,12 @@ class HyperplasticSolverSerrailleHandler::SerrailleSolver<MeshType>::SymmetryBou { return m_mesh_flat_node_boundary.outgoingNormal(); } + + const Rd& + origin() const + { + return m_mesh_flat_node_boundary.origin(); + } size_t numberOfNodes() const diff --git a/src/scheme/SerrailleSolver.hpp b/src/scheme/SerrailleSolver.hpp index 8b61a4c1910340e7ffaa155f81f51ed6cd86e960..a5fb58acd8ec6317a030bf7d2f0d4f5de84d7096 100644 --- a/src/scheme/SerrailleSolver.hpp +++ b/src/scheme/SerrailleSolver.hpp @@ -59,7 +59,9 @@ class HyperplasticSolverSerrailleHandler const std::shared_ptr<const DiscreteFunctionVariant>& sigma, const std::shared_ptr<const ItemValueVariant>& ur, const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr, - const RelaxationType& relaxation_type) const = 0; + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, + const RelaxationType& relaxation_type + ) const = 0; virtual std::tuple<std::shared_ptr<const MeshVariant>, std::shared_ptr<const DiscreteFunctionVariant>,