Skip to content
Snippets Groups Projects
Commit 033af9cb authored by Stéphane Del Pino's avatar Stéphane Del Pino
Browse files

Continue clean-up

Face dispatching is now performed in a separate function
parent 603ade90
No related branches found
No related tags found
1 merge request!11Feature/mpi
...@@ -118,224 +118,15 @@ MeshDispatcher<Dimension>::_buildNbCellToSend() ...@@ -118,224 +118,15 @@ MeshDispatcher<Dimension>::_buildNbCellToSend()
} }
template <int Dimension> template <int Dimension>
MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) void
: m_mesh(mesh), MeshDispatcher<Dimension>::_dispatchFaces()
m_cell_new_owner(_getCellNewOwner()),
m_face_new_owner(_getFaceNewOwner()),
m_node_new_owner(_getNodeNewOwner()),
m_cell_list_to_send_by_proc(_buildCellListToSend()),
m_nb_cell_to_send_by_proc(_buildNbCellToSend()),
m_nb_cell_to_recv_by_proc(parallel::allToAll(m_nb_cell_to_send_by_proc))
{
std::vector<Array<int>> recv_cell_number_by_proc
= this->exchange(mesh.connectivity().cellNumber());
const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix();
std::vector<Array<int>> cell_node_number_to_send_by_proc =
[&] () {
const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber();
std::vector<Array<int>> cell_node_number_to_send_by_proc(parallel::size());
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
std::vector<int> node_number_by_cell_vector;
for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) {
const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j];
const auto& cell_node_list = cell_to_node_matrix[cell_id];
for (size_t r=0; r<cell_node_list.size(); ++r) {
const NodeId& node_id = cell_node_list[r];
node_number_by_cell_vector.push_back(node_number[node_id]);
}
}
cell_node_number_to_send_by_proc[i_rank] = convert_to_array(node_number_by_cell_vector);
}
return cell_node_number_to_send_by_proc;
} ();
std::vector<Array<int>> recv_number_of_node_per_cell_by_proc
= [&] () {
CellValue<int> number_of_node_per_cell(mesh.connectivity());
parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){
number_of_node_per_cell[j] = cell_to_node_matrix[j].size();
});
return this->exchange(number_of_node_per_cell);
} ();
std::vector<Array<int>> recv_cell_node_number_by_proc(parallel::size());
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
recv_cell_node_number_by_proc[i_rank]
= Array<int>(sum(recv_number_of_node_per_cell_by_proc[i_rank]));
}
parallel::exchange(cell_node_number_to_send_by_proc, recv_cell_node_number_by_proc);
const std::unordered_map<int, int> node_number_id_map
= [&] () {
std::unordered_map<int, int> node_number_id_map;
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
int cpt=0;
for (size_t i=0; i<recv_cell_node_number_by_proc[i_rank].size(); ++i) {
int node_number = recv_cell_node_number_by_proc[i_rank][i];
auto [iterator, inserted] = node_number_id_map.insert(std::make_pair(node_number, cpt));
if (inserted) cpt++;
}
}
return node_number_id_map;
} ();
m_new_descriptor.node_number_vector.resize(node_number_id_map.size());
for (const auto& [number, id] : node_number_id_map) {
m_new_descriptor.node_number_vector[id] = number;
}
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
int l=0;
for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) {
std::vector<unsigned int> node_vector;
for (int k=0; k<recv_number_of_node_per_cell_by_proc[i_rank][i]; ++k) {
const auto& searched_node_id = node_number_id_map.find(recv_cell_node_number_by_proc[i_rank][l++]);
Assert(searched_node_id != node_number_id_map.end());
node_vector.push_back(searched_node_id->second);
}
m_new_descriptor.cell_by_node_vector.emplace_back(node_vector);
}
}
const std::unordered_map<int, int> cell_number_id_map
= [&] () {
std::unordered_map<int, int> cell_number_id_map;
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
int cpt=0;
for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) {
int cell_number = recv_cell_number_by_proc[i_rank][i];
auto [iterator, inserted] = cell_number_id_map.insert(std::make_pair(cell_number, cpt));
if (inserted) cpt++;
}
}
return cell_number_id_map;
} ();
m_new_descriptor.cell_number_vector.resize(cell_number_id_map.size());
for (const auto& [number, id] : cell_number_id_map) {
m_new_descriptor.cell_number_vector[id] = number;
}
{
std::vector<Array<CellType>> recv_cell_type_by_proc
= this->exchange(mesh.connectivity().cellType());
m_new_descriptor.cell_type_vector.resize(cell_number_id_map.size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) {
int cell_number = recv_cell_number_by_proc[i_rank][i];
const auto& searched_cell_id = cell_number_id_map.find(cell_number);
Assert(searched_cell_id != cell_number_id_map.end());
m_new_descriptor.cell_type_vector[searched_cell_id->second] = recv_cell_type_by_proc[i_rank][i];
}
}
}
std::vector<Array<const NodeId>> send_node_id_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<bool> tag(mesh.numberOfNodes());
tag.fill(false);
std::vector<NodeId> node_id_vector;
for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) {
const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j];
const auto& cell_node_list = cell_to_node_matrix[cell_id];
for (size_t r=0; r<cell_node_list.size(); ++r) {
const NodeId& node_id = cell_node_list[r];
if (not tag[node_id]) {
node_id_vector.push_back(node_id);
tag[node_id] = true;
}
}
}
send_node_id_by_proc[i_rank] = convert_to_array(node_id_vector);
}
Array<unsigned int> nb_node_to_send_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
nb_node_to_send_by_proc[i_rank] = send_node_id_by_proc[i_rank].size();
}
Array<const unsigned int> nb_node_to_recv_by_proc
= parallel::allToAll(nb_node_to_send_by_proc);
std::vector<Array<const NodeId>> recv_node_id_correspondance_by_proc(parallel::size());
{ {
const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber(); if constexpr (Dimension>1) {
std::vector<Array<const int>> send_node_number_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<int> send_node_number(send_node_id_by_proc[i_rank].size());
const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank];
parallel_for(send_node_number.size(), PASTIS_LAMBDA(const size_t& j){
send_node_number[j] = node_number[send_node_id[j]];
});
send_node_number_by_proc[i_rank] = send_node_number;
}
std::vector<Array<int>> recv_node_number_by_proc(parallel::size());
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
recv_node_number_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]);
}
parallel::exchange(send_node_number_by_proc, recv_node_number_by_proc);
for (size_t i_rank=0; i_rank<nb_node_to_recv_by_proc.size(); ++i_rank) {
Array<NodeId> node_id_correspondace(nb_node_to_recv_by_proc[i_rank]);
for (size_t l=0; l<nb_node_to_recv_by_proc[i_rank]; ++l) {
const int& node_number = recv_node_number_by_proc[i_rank][l];
const auto& searched_node_id = node_number_id_map.find(node_number);
Assert(searched_node_id != node_number_id_map.end());
node_id_correspondace[l] = searched_node_id->second;
}
recv_node_id_correspondance_by_proc[i_rank] = node_id_correspondace;
}
}
{
std::vector<Array<int>> recv_cell_new_owner_by_proc
= this->exchange(m_cell_new_owner);
m_new_descriptor.cell_owner_vector.resize(cell_number_id_map.size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) {
int cell_number = recv_cell_number_by_proc[i_rank][i];
const auto& searched_cell_id = cell_number_id_map.find(cell_number);
Assert(searched_cell_id != cell_number_id_map.end());
m_new_descriptor.cell_owner_vector[searched_cell_id->second] = recv_cell_new_owner_by_proc[i_rank][i];
}
}
}
{
std::vector<Array<const int>> send_node_owner_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<int> send_node_owner(nb_node_to_send_by_proc[i_rank]);
const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank];
parallel_for(send_node_id.size(), PASTIS_LAMBDA(const size_t& r) {
const NodeId& node_id = send_node_id[r];
send_node_owner[r] = m_node_new_owner[node_id];
});
send_node_owner_by_proc[i_rank] = send_node_owner;
}
std::vector<Array<int>> recv_node_owner_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
recv_node_owner_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]);
}
parallel::exchange(send_node_owner_by_proc, recv_node_owner_by_proc);
m_new_descriptor.node_owner_vector.resize(m_new_descriptor.node_number_vector.size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
for (size_t r=0; r<recv_node_owner_by_proc[i_rank].size(); ++r) {
const NodeId& node_id = recv_node_id_correspondance_by_proc[i_rank][r];
m_new_descriptor.node_owner_vector[node_id] = recv_node_owner_by_proc[i_rank][r];
}
}
}
if constexpr(Dimension>1) { // Faces
std::vector<Array<int>> recv_number_of_face_per_cell_by_proc = std::vector<Array<int>> recv_number_of_face_per_cell_by_proc =
[&] () { [&] () {
CellValue<int> number_of_face_per_cell(mesh.connectivity()); CellValue<int> number_of_face_per_cell(m_mesh.connectivity());
const auto& cell_to_face_matrix = mesh.connectivity().cellToFaceMatrix(); const auto& cell_to_face_matrix = m_mesh.connectivity().cellToFaceMatrix();
parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){ parallel_for(m_mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){
number_of_face_per_cell[j] = cell_to_face_matrix[j].size(); number_of_face_per_cell[j] = cell_to_face_matrix[j].size();
}); });
...@@ -344,8 +135,8 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -344,8 +135,8 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
std::vector<Array<int>> recv_cell_face_number_by_proc(parallel::size()); std::vector<Array<int>> recv_cell_face_number_by_proc(parallel::size());
{ {
const auto& cell_to_face_matrix = mesh.connectivity().cellToFaceMatrix(); const auto& cell_to_face_matrix = m_mesh.connectivity().cellToFaceMatrix();
const FaceValue<const int>& face_number = mesh.connectivity().faceNumber(); const FaceValue<const int>& face_number = m_mesh.connectivity().faceNumber();
std::vector<Array<int>> cell_face_number_to_send_by_proc(parallel::size()); std::vector<Array<int>> cell_face_number_to_send_by_proc(parallel::size());
{ {
...@@ -393,7 +184,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -393,7 +184,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
int l=0; int l=0;
for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { for (size_t i=0; i<m_recv_cell_number_by_proc[i_rank].size(); ++i) {
std::vector<unsigned int> face_vector; std::vector<unsigned int> face_vector;
for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) { for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) {
const auto& searched_face_id = face_number_id_map.find(recv_cell_face_number_by_proc[i_rank][l++]); const auto& searched_face_id = face_number_id_map.find(recv_cell_face_number_by_proc[i_rank][l++]);
...@@ -407,7 +198,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -407,7 +198,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
{ {
std::vector<Array<bool>> cell_face_is_reversed_to_send_by_proc(parallel::size()); std::vector<Array<bool>> cell_face_is_reversed_to_send_by_proc(parallel::size());
{ {
const auto& cell_face_is_reversed = mesh.connectivity().cellFaceIsReversed(); const auto& cell_face_is_reversed = m_mesh.connectivity().cellFaceIsReversed();
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
std::vector<bool> face_is_reversed_by_cell_vector; std::vector<bool> face_is_reversed_by_cell_vector;
for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) {
...@@ -431,7 +222,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -431,7 +222,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
int l=0; int l=0;
for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { for (size_t i=0; i<m_recv_cell_number_by_proc[i_rank].size(); ++i) {
std::vector<bool> face_is_reversed_vector; std::vector<bool> face_is_reversed_vector;
for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) { for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) {
face_is_reversed_vector.push_back(recv_cell_face_is_reversed_by_proc[i_rank][l++]); face_is_reversed_vector.push_back(recv_cell_face_is_reversed_by_proc[i_rank][l++]);
...@@ -443,10 +234,10 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -443,10 +234,10 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
std::vector<Array<const FaceId>> send_face_id_by_proc(parallel::size()); std::vector<Array<const FaceId>> send_face_id_by_proc(parallel::size());
{ {
const auto& cell_to_face_matrix = mesh.connectivity().cellToFaceMatrix(); const auto& cell_to_face_matrix = m_mesh.connectivity().cellToFaceMatrix();
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<bool> tag(mesh.numberOfFaces()); Array<bool> tag(m_mesh.numberOfFaces());
tag.fill(false); tag.fill(false);
std::vector<FaceId> face_id_vector; std::vector<FaceId> face_id_vector;
for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) {
...@@ -466,7 +257,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -466,7 +257,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
std::vector<Array<const int>> send_face_number_by_proc(parallel::size()); std::vector<Array<const int>> send_face_number_by_proc(parallel::size());
{ {
const FaceValue<const int>& face_number = mesh.connectivity().faceNumber(); const FaceValue<const int>& face_number = m_mesh.connectivity().faceNumber();
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<int> send_face_number(send_face_id_by_proc[i_rank].size()); Array<int> send_face_number(send_face_id_by_proc[i_rank].size());
...@@ -537,9 +328,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -537,9 +328,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
std::vector<Array<int>> recv_face_number_of_node_by_proc(parallel::size()); std::vector<Array<int>> recv_face_number_of_node_by_proc(parallel::size());
{ {
const auto& face_to_node_matrix = mesh.connectivity().faceToNodeMatrix(); const auto& face_to_node_matrix = m_mesh.connectivity().faceToNodeMatrix();
FaceValue<int> number_of_node_per_face(mesh.connectivity()); FaceValue<int> number_of_node_per_face(m_mesh.connectivity());
parallel_for(mesh.numberOfFaces(), PASTIS_LAMBDA(const FaceId& j){ parallel_for(m_mesh.numberOfFaces(), PASTIS_LAMBDA(const FaceId& j){
number_of_node_per_face[j] = face_to_node_matrix[j].size(); number_of_node_per_face[j] = face_to_node_matrix[j].size();
}); });
...@@ -562,8 +353,8 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -562,8 +353,8 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
std::vector<Array<int>> face_node_number_to_send_by_proc(parallel::size()); std::vector<Array<int>> face_node_number_to_send_by_proc(parallel::size());
{ {
const auto& face_to_node_matrix = mesh.connectivity().faceToNodeMatrix(); const auto& face_to_node_matrix = m_mesh.connectivity().faceToNodeMatrix();
const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber(); const NodeValue<const int>& node_number = m_mesh.connectivity().nodeNumber();
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
std::vector<int> node_number_by_face_vector; std::vector<int> node_number_by_face_vector;
for (size_t l=0; l<send_face_id_by_proc[i_rank].size(); ++l) { for (size_t l=0; l<send_face_id_by_proc[i_rank].size(); ++l) {
...@@ -592,8 +383,8 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -592,8 +383,8 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
for (size_t i=0; i<recv_face_number_by_proc[i_rank].size(); ++i) { for (size_t i=0; i<recv_face_number_by_proc[i_rank].size(); ++i) {
std::vector<unsigned int> node_vector; std::vector<unsigned int> node_vector;
for (int k=0; k<recv_face_number_of_node_by_proc[i_rank][i]; ++k) { for (int k=0; k<recv_face_number_of_node_by_proc[i_rank][i]; ++k) {
const auto& searched_node_id = node_number_id_map.find(recv_face_node_number_by_proc[i_rank][l++]); const auto& searched_node_id = m_node_number_id_map.find(recv_face_node_number_by_proc[i_rank][l++]);
Assert(searched_node_id != node_number_id_map.end()); Assert(searched_node_id != m_node_number_id_map.end());
node_vector.push_back(searched_node_id->second); node_vector.push_back(searched_node_id->second);
} }
m_new_descriptor.face_to_node_vector.emplace_back(node_vector); m_new_descriptor.face_to_node_vector.emplace_back(node_vector);
...@@ -603,7 +394,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -603,7 +394,7 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
// Getting references // Getting references
Array<const size_t> number_of_ref_face_list_per_proc Array<const size_t> number_of_ref_face_list_per_proc
= parallel::allGather(mesh.connectivity().numberOfRefFaceList()); = parallel::allGather(m_mesh.connectivity().numberOfRefFaceList());
const size_t number_of_face_list_sender const size_t number_of_face_list_sender
= [&] () { = [&] () {
...@@ -640,9 +431,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -640,9 +431,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
// sending references tags // sending references tags
Array<RefId::TagNumberType> ref_tag_list{number_of_ref_face_list_per_proc[sender_rank]}; Array<RefId::TagNumberType> ref_tag_list{number_of_ref_face_list_per_proc[sender_rank]};
if (parallel::rank() == sender_rank){ if (parallel::rank() == sender_rank){
for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); for (size_t i_ref_face_list=0; i_ref_face_list<m_mesh.connectivity().numberOfRefFaceList();
++i_ref_face_list) { ++i_ref_face_list) {
auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); auto ref_face_list = m_mesh.connectivity().refFaceList(i_ref_face_list);
ref_tag_list[i_ref_face_list] = ref_face_list.refId().tagNumber(); ref_tag_list[i_ref_face_list] = ref_face_list.refId().tagNumber();
} }
} }
...@@ -651,9 +442,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -651,9 +442,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
// sending references name size // sending references name size
Array<size_t> ref_name_size_list{number_of_ref_face_list_per_proc[sender_rank]}; Array<size_t> ref_name_size_list{number_of_ref_face_list_per_proc[sender_rank]};
if (parallel::rank() == sender_rank){ if (parallel::rank() == sender_rank){
for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); for (size_t i_ref_face_list=0; i_ref_face_list<m_mesh.connectivity().numberOfRefFaceList();
++i_ref_face_list) { ++i_ref_face_list) {
auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); auto ref_face_list = m_mesh.connectivity().refFaceList(i_ref_face_list);
ref_name_size_list[i_ref_face_list] = ref_face_list.refId().tagName().size(); ref_name_size_list[i_ref_face_list] = ref_face_list.refId().tagName().size();
} }
} }
...@@ -663,9 +454,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -663,9 +454,9 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
Array<RefId::TagNameType::value_type> ref_name_cat{sum(ref_name_size_list)}; Array<RefId::TagNameType::value_type> ref_name_cat{sum(ref_name_size_list)};
if (parallel::rank() == sender_rank){ if (parallel::rank() == sender_rank){
size_t i_char=0; size_t i_char=0;
for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); for (size_t i_ref_face_list=0; i_ref_face_list<m_mesh.connectivity().numberOfRefFaceList();
++i_ref_face_list) { ++i_ref_face_list) {
auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); auto ref_face_list = m_mesh.connectivity().refFaceList(i_ref_face_list);
for (auto c : ref_face_list.refId().tagName()) { for (auto c : ref_face_list.refId().tagName()) {
ref_name_cat[i_char++] = c; ref_name_cat[i_char++] = c;
} }
...@@ -691,14 +482,14 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -691,14 +482,14 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
constexpr size_t block_size = sizeof(block_type); constexpr size_t block_size = sizeof(block_type);
const size_t nb_block = ref_id_list.size()/block_size + (ref_id_list.size()%block_size != 0); const size_t nb_block = ref_id_list.size()/block_size + (ref_id_list.size()%block_size != 0);
for (size_t i_block=0; i_block<nb_block; ++i_block) { for (size_t i_block=0; i_block<nb_block; ++i_block) {
FaceValue<block_type> face_references(mesh.connectivity()); FaceValue<block_type> face_references(m_mesh.connectivity());
face_references.fill(0); face_references.fill(0);
if (mesh.connectivity().numberOfRefFaceList() > 0) { if (m_mesh.connectivity().numberOfRefFaceList() > 0) {
const size_t max_i_ref = std::min(ref_id_list.size(), block_size*(i_block+1)); const size_t max_i_ref = std::min(ref_id_list.size(), block_size*(i_block+1));
for (size_t i_ref=block_size*i_block, i=0; i_ref<max_i_ref; ++i_ref, ++i) { for (size_t i_ref=block_size*i_block, i=0; i_ref<max_i_ref; ++i_ref, ++i) {
block_type ref_bit{1<<i}; block_type ref_bit{1<<i};
auto ref_face_list = mesh.connectivity().refFaceList(i_ref); auto ref_face_list = m_mesh.connectivity().refFaceList(i_ref);
const auto& face_list = ref_face_list.faceList(); const auto& face_list = ref_face_list.faceList();
for (size_t i_face=0; i_face<face_list.size(); ++i_face) { for (size_t i_face=0; i_face<face_list.size(); ++i_face) {
...@@ -755,8 +546,223 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) ...@@ -755,8 +546,223 @@ MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
} }
} }
} }
}
}
template <int Dimension>
MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh)
: m_mesh(mesh),
m_cell_new_owner(_getCellNewOwner()),
m_face_new_owner(_getFaceNewOwner()),
m_node_new_owner(_getNodeNewOwner()),
m_cell_list_to_send_by_proc(_buildCellListToSend()),
m_nb_cell_to_send_by_proc(_buildNbCellToSend()),
m_nb_cell_to_recv_by_proc(parallel::allToAll(m_nb_cell_to_send_by_proc))
{
m_recv_cell_number_by_proc = this->exchange(mesh.connectivity().cellNumber());
const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix();
std::vector<Array<int>> cell_node_number_to_send_by_proc =
[&] () {
const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber();
std::vector<Array<int>> cell_node_number_to_send_by_proc(parallel::size());
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
std::vector<int> node_number_by_cell_vector;
for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) {
const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j];
const auto& cell_node_list = cell_to_node_matrix[cell_id];
for (size_t r=0; r<cell_node_list.size(); ++r) {
const NodeId& node_id = cell_node_list[r];
node_number_by_cell_vector.push_back(node_number[node_id]);
}
}
cell_node_number_to_send_by_proc[i_rank] = convert_to_array(node_number_by_cell_vector);
}
return cell_node_number_to_send_by_proc;
} ();
std::vector<Array<int>> recv_number_of_node_per_cell_by_proc
= [&] () {
CellValue<int> number_of_node_per_cell(mesh.connectivity());
parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){
number_of_node_per_cell[j] = cell_to_node_matrix[j].size();
});
return this->exchange(number_of_node_per_cell);
} ();
std::vector<Array<int>> recv_cell_node_number_by_proc(parallel::size());
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
recv_cell_node_number_by_proc[i_rank]
= Array<int>(sum(recv_number_of_node_per_cell_by_proc[i_rank]));
}
parallel::exchange(cell_node_number_to_send_by_proc, recv_cell_node_number_by_proc);
m_node_number_id_map =
[&] () {
std::unordered_map<int, int> node_number_id_map;
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
int cpt=0;
for (size_t i=0; i<recv_cell_node_number_by_proc[i_rank].size(); ++i) {
int node_number = recv_cell_node_number_by_proc[i_rank][i];
auto [iterator, inserted] = node_number_id_map.insert(std::make_pair(node_number, cpt));
if (inserted) cpt++;
}
}
return node_number_id_map;
} ();
m_new_descriptor.node_number_vector.resize(m_node_number_id_map.size());
for (const auto& [number, id] : m_node_number_id_map) {
m_new_descriptor.node_number_vector[id] = number;
}
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
int l=0;
for (size_t i=0; i<m_recv_cell_number_by_proc[i_rank].size(); ++i) {
std::vector<unsigned int> node_vector;
for (int k=0; k<recv_number_of_node_per_cell_by_proc[i_rank][i]; ++k) {
const auto& searched_node_id = m_node_number_id_map.find(recv_cell_node_number_by_proc[i_rank][l++]);
Assert(searched_node_id != m_node_number_id_map.end());
node_vector.push_back(searched_node_id->second);
}
m_new_descriptor.cell_by_node_vector.emplace_back(node_vector);
}
}
const std::unordered_map<int, int> cell_number_id_map
= [&] () {
std::unordered_map<int, int> cell_number_id_map;
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
int cpt=0;
for (size_t i=0; i<m_recv_cell_number_by_proc[i_rank].size(); ++i) {
int cell_number = m_recv_cell_number_by_proc[i_rank][i];
auto [iterator, inserted] = cell_number_id_map.insert(std::make_pair(cell_number, cpt));
if (inserted) cpt++;
}
}
return cell_number_id_map;
} ();
m_new_descriptor.cell_number_vector.resize(cell_number_id_map.size());
for (const auto& [number, id] : cell_number_id_map) {
m_new_descriptor.cell_number_vector[id] = number;
}
{
std::vector<Array<CellType>> recv_cell_type_by_proc
= this->exchange(mesh.connectivity().cellType());
m_new_descriptor.cell_type_vector.resize(cell_number_id_map.size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
for (size_t i=0; i<m_recv_cell_number_by_proc[i_rank].size(); ++i) {
int cell_number = m_recv_cell_number_by_proc[i_rank][i];
const auto& searched_cell_id = cell_number_id_map.find(cell_number);
Assert(searched_cell_id != cell_number_id_map.end());
m_new_descriptor.cell_type_vector[searched_cell_id->second] = recv_cell_type_by_proc[i_rank][i];
}
}
}
std::vector<Array<const NodeId>> send_node_id_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<bool> tag(mesh.numberOfNodes());
tag.fill(false);
std::vector<NodeId> node_id_vector;
for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) {
const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j];
const auto& cell_node_list = cell_to_node_matrix[cell_id];
for (size_t r=0; r<cell_node_list.size(); ++r) {
const NodeId& node_id = cell_node_list[r];
if (not tag[node_id]) {
node_id_vector.push_back(node_id);
tag[node_id] = true;
}
}
}
send_node_id_by_proc[i_rank] = convert_to_array(node_id_vector);
}
Array<unsigned int> nb_node_to_send_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
nb_node_to_send_by_proc[i_rank] = send_node_id_by_proc[i_rank].size();
} }
Array<const unsigned int> nb_node_to_recv_by_proc
= parallel::allToAll(nb_node_to_send_by_proc);
std::vector<Array<const NodeId>> recv_node_id_correspondance_by_proc(parallel::size());
{
const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber();
std::vector<Array<const int>> send_node_number_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<int> send_node_number(send_node_id_by_proc[i_rank].size());
const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank];
parallel_for(send_node_number.size(), PASTIS_LAMBDA(const size_t& j){
send_node_number[j] = node_number[send_node_id[j]];
});
send_node_number_by_proc[i_rank] = send_node_number;
}
std::vector<Array<int>> recv_node_number_by_proc(parallel::size());
for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) {
recv_node_number_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]);
}
parallel::exchange(send_node_number_by_proc, recv_node_number_by_proc);
for (size_t i_rank=0; i_rank<nb_node_to_recv_by_proc.size(); ++i_rank) {
Array<NodeId> node_id_correspondace(nb_node_to_recv_by_proc[i_rank]);
for (size_t l=0; l<nb_node_to_recv_by_proc[i_rank]; ++l) {
const int& node_number = recv_node_number_by_proc[i_rank][l];
const auto& searched_node_id = m_node_number_id_map.find(node_number);
Assert(searched_node_id != m_node_number_id_map.end());
node_id_correspondace[l] = searched_node_id->second;
}
recv_node_id_correspondance_by_proc[i_rank] = node_id_correspondace;
}
}
{
std::vector<Array<int>> recv_cell_new_owner_by_proc
= this->exchange(m_cell_new_owner);
m_new_descriptor.cell_owner_vector.resize(cell_number_id_map.size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
for (size_t i=0; i<m_recv_cell_number_by_proc[i_rank].size(); ++i) {
int cell_number = m_recv_cell_number_by_proc[i_rank][i];
const auto& searched_cell_id = cell_number_id_map.find(cell_number);
Assert(searched_cell_id != cell_number_id_map.end());
m_new_descriptor.cell_owner_vector[searched_cell_id->second] = recv_cell_new_owner_by_proc[i_rank][i];
}
}
}
{
std::vector<Array<const int>> send_node_owner_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
Array<int> send_node_owner(nb_node_to_send_by_proc[i_rank]);
const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank];
parallel_for(send_node_id.size(), PASTIS_LAMBDA(const size_t& r) {
const NodeId& node_id = send_node_id[r];
send_node_owner[r] = m_node_new_owner[node_id];
});
send_node_owner_by_proc[i_rank] = send_node_owner;
}
std::vector<Array<int>> recv_node_owner_by_proc(parallel::size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
recv_node_owner_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]);
}
parallel::exchange(send_node_owner_by_proc, recv_node_owner_by_proc);
m_new_descriptor.node_owner_vector.resize(m_new_descriptor.node_number_vector.size());
for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) {
for (size_t r=0; r<recv_node_owner_by_proc[i_rank].size(); ++r) {
const NodeId& node_id = recv_node_id_correspondance_by_proc[i_rank][r];
m_new_descriptor.node_owner_vector[node_id] = recv_node_owner_by_proc[i_rank][r];
}
}
}
this->_dispatchFaces();
using ConnectivityType = Connectivity<Dimension>; using ConnectivityType = Connectivity<Dimension>;
std::shared_ptr p_connectivity std::shared_ptr p_connectivity
......
...@@ -5,6 +5,8 @@ ...@@ -5,6 +5,8 @@
#include <ItemValue.hpp> #include <ItemValue.hpp>
#include <ItemValueUtils.hpp> #include <ItemValueUtils.hpp>
#include <unordered_map>
template <int Dimension> template <int Dimension>
class MeshDispatcher class MeshDispatcher
{ {
...@@ -28,6 +30,9 @@ class MeshDispatcher ...@@ -28,6 +30,9 @@ class MeshDispatcher
Array<int> m_nb_cell_to_send_by_proc; Array<int> m_nb_cell_to_send_by_proc;
Array<int> m_nb_cell_to_recv_by_proc; Array<int> m_nb_cell_to_recv_by_proc;
std::vector<Array<int>> m_recv_cell_number_by_proc;
std::unordered_map<int, int> m_node_number_id_map;
CellValue<int> _getCellNewOwner(); CellValue<int> _getCellNewOwner();
FaceValue<int> _getFaceNewOwner(); FaceValue<int> _getFaceNewOwner();
NodeValue<int> _getNodeNewOwner(); NodeValue<int> _getNodeNewOwner();
...@@ -35,6 +40,8 @@ class MeshDispatcher ...@@ -35,6 +40,8 @@ class MeshDispatcher
const CellListToSendByProc _buildCellListToSend() const; const CellListToSendByProc _buildCellListToSend() const;
Array<int> _buildNbCellToSend(); Array<int> _buildNbCellToSend();
void _dispatchFaces();
public: public:
const std::shared_ptr<MeshType> dispatchedMesh() const const std::shared_ptr<MeshType> dispatchedMesh() const
...@@ -84,12 +91,6 @@ class MeshDispatcher ...@@ -84,12 +91,6 @@ class MeshDispatcher
return recv_cell_value_by_proc; return recv_cell_value_by_proc;
} }
[[deprecated]]
const CellListToSendByProc& cell_list_to_send_by_proc() const
{
return m_cell_list_to_send_by_proc;
}
MeshDispatcher(const MeshType& mesh); MeshDispatcher(const MeshType& mesh);
MeshDispatcher(const MeshDispatcher&) = delete; MeshDispatcher(const MeshDispatcher&) = delete;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment