Skip to content
Snippets Groups Projects
Commit e5e0f5e7 authored by chantrait's avatar chantrait
Browse files

couplage hyperelasticity q=1 ok

parent 9911f348
No related branches found
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
src/scheme/HyperelasticSolver.cpp
+ 144
− 3
View file @ e5e0f5e7
......@@ -421,6 +421,10 @@ class HyperelasticSolverHandler::HyperelasticSolver final : public HyperelasticS
void _applySymmetryBC(const BoundaryConditionList& bc_list, NodeValue<Rdxd>& Ar, NodeValue<Rd>& br) const;
void _applyVelocityBC(const BoundaryConditionList& bc_list, NodeValue<Rdxd>& Ar, NodeValue<Rd>& br) const;
void _applyCouplingBC(const BoundaryConditionList& bc_list, NodeValue<Rdxd>& Ar, NodeValue<Rd>& br) const;
void _applyCouplingBC2(const MeshType& mesh,
const BoundaryConditionList& bc_list,
NodeValue<Rd>& ur,
NodeValuePerCell<Rd>& Fjr) const;
void
_applyBoundaryConditions(const BoundaryConditionList& bc_list,
const MeshType& mesh,
......@@ -436,7 +440,7 @@ class HyperelasticSolverHandler::HyperelasticSolver final : public HyperelasticS
this->_applyVelocityBC(bc_list, Ar, br);
}
NodeValue<const Rd>
NodeValue<Rd>
_computeUr(const MeshType& mesh, const NodeValue<Rdxd>& Ar, const NodeValue<Rd>& br) const
{
NodeValue<Rd> u{mesh.connectivity()};
......@@ -510,8 +514,10 @@ class HyperelasticSolverHandler::HyperelasticSolver final : public HyperelasticS
synchronize(Ar);
synchronize(br);
NodeValue<const Rd> ur = this->_computeUr(mesh, Ar, br);
NodeValuePerCell<const Rd> Fjr = this->_computeFjr(mesh, Ajr, ur, u, sigma);
NodeValue<Rd> ur = this->_computeUr(mesh, Ar, br);
NodeValuePerCell<Rd> Fjr = this->_computeFjr(mesh, Ajr, ur, u, sigma);
this->_applyCouplingBC2(mesh, bc_list, ur, Fjr);
return std::make_tuple(std::make_shared<const ItemValueVariant>(ur),
std::make_shared<const SubItemValuePerItemVariant>(Fjr));
......@@ -946,6 +952,141 @@ HyperelasticSolverHandler::HyperelasticSolver<Dimension>::_applyCouplingBC(const
}
}
template <size_t Dimension>
void
HyperelasticSolverHandler::HyperelasticSolver<Dimension>::_applyCouplingBC2(const MeshType& mesh,
const BoundaryConditionList& bc_list,
NodeValue<Rd>& ur,
NodeValuePerCell<Rd>& Fjr) 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<CouplingBoundaryCondition, T>) {
const auto& node_list = bc.nodeList();
/* std::cout << "\033[01;31m" */
/* << "un applyCoupling" */
/* << "Dimension: " << Dimension << "\033[00;00m" << std::endl; */
/* std::cout << "\033[01;31m" */
/* << "node_list.size()" << node_list.size() << "\033[00;00m" << std::endl; */
#ifdef PUGS_HAS_COSTO
auto costo = parallel::Messenger::getInstance().myCoupling;
auto rank = costo->myGlobalRank();
if (rank == 0) {
const size_t& n = node_list.size();
const auto& node_to_cell_matrix = mesh.connectivity().nodeToCellMatrix();
const auto& node_local_numbers_in_their_cells = mesh.connectivity().nodeLocalNumbersInTheirCells();
for (size_t i = 0; i < n; i++) {
const NodeId node_id = node_list[i];
const auto& node_to_cell = node_to_cell_matrix[node_id];
const auto& node_local_number_in_its_cells = node_local_numbers_in_their_cells.itemArray(node_id);
/* std::cout << "\033[01;33m" */
/* << "rank_0: " << node_id << " " << node_to_cell.size() << "\033[00;00m" << std::endl; */
for (size_t j = 0; j < node_to_cell.size(); j++) {
const CellId J = node_to_cell[j];
const unsigned int R = node_local_number_in_its_cells[j];
/* Fjr(J, R) = CR_Fjr(J, R); */
}
/* ur[node_id] = CR_ur[node_id]; */
}
}
if (rank == 1) {
/* std::vector<double> recv; */
/* tag = 310; */
/* costo->recvData(recv, dest, tag); */
/* CR_Fjr = ... */
/* CR_ur = */
const size_t& n = node_list.size();
const auto& node_to_cell_matrix = mesh.connectivity().nodeToCellMatrix();
const auto& node_local_numbers_in_their_cells = mesh.connectivity().nodeLocalNumbersInTheirCells();
for (size_t i = 0; i < n; i++) {
const NodeId node_id = node_list[i];
const auto& node_to_cell = node_to_cell_matrix[node_id];
const auto& node_local_number_in_its_cells = node_local_numbers_in_their_cells.itemArray(node_id);
/* std::cout << "\033[01;32m" */
/* << "rank_1: " << node_id << " " << node_to_cell.size() << "\033[00;00m" << std::endl; */
for (size_t j = 0; j < node_to_cell.size(); j++) {
const CellId J = node_to_cell[j];
const unsigned int R = node_local_number_in_its_cells[j];
/* Fjr(J, R) = CR_Fjr(J, R); */
}
/* ur[node_id] = CR_ur[node_id]; */
}
}
/* const int dest = costo->myGlobalSize() - 1; */
/* int tag = 300; */
/* std::vector<int> shape; */
/* shape.resize(3); */
/* shape[0] = node_list.size(); */
/* shape[1] = Dimension + Dimension * Dimension; */
/* /\* shape[1] = Dimension; *\/ */
/* shape[2] = 0; */
/* std::vector<double> data; */
/* data.resize(shape[0] * shape[1] + shape[2]); */
/* Array<TinyVector<Dimension>> br_extracted(node_list.size()); */
/* Array<TinyMatrix<Dimension>> Ar_extracted(node_list.size()); */
/* parallel_for( */
/* node_list.size(), PUGS_LAMBDA(size_t i_node) { */
/* NodeId node_id = node_list[i_node]; */
/* for (size_t i_dim = 0; i_dim < Dimension; i_dim++) { */
/* br_extracted[i_node] = br[node_id]; */
/* Ar_extracted[i_node] = Ar[node_id]; */
/* } */
/* }); */
/* /\*TODO: serializer directement Ar et br pour eviter une copie*\/ */
/* for (size_t i_node = 0; i_node < node_list.size(); ++i_node) { */
/* for (unsigned short i_dim = 0; i_dim < Dimension; ++i_dim) { */
/* data[i_node * shape[1] + i_dim] = br_extracted[i_node][i_dim]; */
/* for (size_t j_dim = 0; j_dim < Dimension; j_dim++) { */
/* data[i_node * shape[1] + Dimension + i_dim * Dimension + j_dim] = Ar_extracted[i_node](i_dim,
j_dim);
*/
/* } */
/* } */
/* } */
/* costo->sendData(shape, data, dest, tag); */
/* std::vector<double> recv; */
/* tag = 310; */
/* costo->recvData(recv, dest, tag); */
/* for (size_t i_node = 0; i_node < node_list.size(); ++i_node) { */
/* for (unsigned short i_dim = 0; i_dim < Dimension; ++i_dim) { */
/* br_extracted[i_node][i_dim] = recv[i_node * shape[1] + i_dim]; */
/* for (size_t j_dim = 0; j_dim < Dimension; j_dim++) { */
/* Ar_extracted[i_node](i_dim, j_dim) = recv[i_node * shape[1] + Dimension + i_dim * Dimension +
j_dim];
*/
/* } */
/* } */
/* } */
/* parallel_for( */
/* node_list.size(), PUGS_LAMBDA(size_t i_node) { */
/* NodeId node_id = node_list[i_node]; */
/* for (size_t i_dim = 0; i_dim < Dimension; i_dim++) { */
/* br[node_id] = br_extracted[i_node]; */
/* Ar[node_id] = Ar_extracted[i_node]; */
/* } */
/* }); */
#endif PUGS_HAS_COSTO
}
},
boundary_condition);
}
}
template <size_t Dimension>
class HyperelasticSolverHandler::HyperelasticSolver<Dimension>::FixedBoundaryCondition
{
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please to comment