10MoFEMErrorCode MatPartitioningApply_Parmetis_MoFEM(MatPartitioning part,
24 : cOre(const_cast<
MoFEM::
Core &>(core)), dEbug(false) {}
28 ParallelComm *pcomm = ParallelComm::get_pcomm(
32 auto get_pstatus = [&](
const auto ent) {
33 unsigned char pstatus;
36 "can not get pstatus");
40 auto get_sharing_procs = [&](
const auto ent,
const auto pstatus) {
41 std::vector<int> sharing_procs(MAX_SHARING_PROCS, -1);
42 if (pstatus & PSTATUS_MULTISHARED) {
45 pcomm->sharedps_tag(), &ent, 1, &sharing_procs[0]),
46 "can not ger sharing_procs_ptr");
47 }
else if (pstatus & PSTATUS_SHARED) {
50 1, &sharing_procs[0]),
51 "can not get sharing proc");
56 auto get_sharing_handles = [&](
const auto ent,
const auto pstatus) {
57 std::vector<EntityHandle> sharing_handles(MAX_SHARING_PROCS, 0);
58 if (pstatus & PSTATUS_MULTISHARED) {
61 pcomm->sharedhs_tag(), &ent, 1, &sharing_handles[0]),
62 "get shared handles");
63 }
else if (pstatus & PSTATUS_SHARED) {
66 1, &sharing_handles[0]),
67 "get sharing handle");
69 return sharing_handles;
73 std::vector<std::vector<EntityHandle>> sbuffer(m_field.
get_comm_size());
75 for (
auto ent : ents) {
77 auto pstatus = get_pstatus(ent);
79 auto sharing_procs = get_sharing_procs(ent, pstatus);
80 auto sharing_handles = get_sharing_handles(ent, pstatus);
83 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"pstatus " << std::bitset<8>(pstatus);
86 for (
int proc = 0; proc < MAX_SHARING_PROCS && -1 != sharing_procs[proc];
88 if (sharing_procs[proc] == -1)
90 "sharing processor not set");
95 const auto handle_on_sharing_proc = sharing_handles[proc];
96 sbuffer[sharing_procs[proc]].push_back(handle_on_sharing_proc);
98 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"send %lu (%lu) to %d at %d\n", ent,
99 handle_on_sharing_proc, sharing_procs[proc],
102 if (!(pstatus & PSTATUS_MULTISHARED))
109 std::vector<int> sbuffer_lengths(
111 const size_t block_size =
sizeof(
EntityHandle) /
sizeof(
int);
114 if (!sbuffer[proc].empty()) {
116 sbuffer_lengths[proc] = sbuffer[proc].size() * block_size;
121 sbuffer_lengths[proc] = 0;
133 CHKERR PetscGatherNumberOfMessages(comm, NULL, &sbuffer_lengths[0], &nrecvs);
139 CHKERR PetscGatherMessageLengths(comm, nsends, nrecvs, &sbuffer_lengths[0],
147 CHKERR PetscCommGetNewTag(comm, &tag);
152 MPI_Request *r_waits;
156 CHKERR PetscPostIrecvInt(comm, tag, nrecvs, onodes, olengths, &rbuf,
159 MPI_Request *s_waits;
160 CHKERR PetscMalloc1(nsends, &s_waits);
163 for (
int proc = 0, kk = 0; proc < m_field.
get_comm_size(); proc++) {
164 if (!sbuffer_lengths[proc])
166 CHKERR MPI_Isend(&(sbuffer[proc])[0],
167 sbuffer_lengths[proc],
169 tag, comm, s_waits + kk);
175 CHKERR MPI_Waitall(nrecvs, r_waits, &status[0]);
179 CHKERR MPI_Waitall(nsends, s_waits, &status[0]);
182 MOFEM_LOG_C(
"SYNC", Sev::verbose,
"Rank %d nb. before ents %u\n",
187 for (
int kk = 0; kk < nrecvs; kk++) {
189 int len = olengths[kk];
190 int *data_from_proc = rbuf[kk];
192 for (
int ee = 0; ee < len;) {
199 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"received %lu from %d at %d\n", ent,
205 MOFEM_LOG_C(
"SYNC", Sev::verbose,
"Rank %d nb. after ents %u",
209 CHKERR PetscFree(s_waits);
210 CHKERR PetscFree(rbuf[0]);
212 CHKERR PetscFree(r_waits);
214 CHKERR PetscFree(olengths);
215 CHKERR PetscCommDestroy(&comm);
229 CHKERR m_field.
get_moab().get_entities_by_handle(idm, ents,
false);
238 ParallelComm *pcomm = ParallelComm::get_pcomm(
244 CHKERR pcomm->filter_pstatus(shared, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1,
246 CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
247 PSTATUS_OR, -1,
nullptr);
252 auto get_pstatus = [&](
const auto ent) {
253 unsigned char pstatus;
256 "can not get pstatus");
260 auto get_sharing_procs = [&](
const auto ent,
const auto pstatus) {
261 std::vector<int> sharing_procs(MAX_SHARING_PROCS, -1);
262 if (pstatus & PSTATUS_MULTISHARED) {
265 pcomm->sharedps_tag(), &ent, 1, &sharing_procs[0]),
266 "can not ger sharing_procs_ptr");
267 }
else if (pstatus & PSTATUS_SHARED) {
270 1, &sharing_procs[0]),
271 "can not get sharing proc");
273 return sharing_procs;
276 auto get_sharing_handles = [&](
const auto ent,
const auto pstatus) {
277 std::vector<EntityHandle> sharing_handles(MAX_SHARING_PROCS, 0);
278 if (pstatus & PSTATUS_MULTISHARED) {
281 pcomm->sharedhs_tag(), &ent, 1, &sharing_handles[0]),
282 "get shared handles");
283 }
else if (pstatus & PSTATUS_SHARED) {
286 1, &sharing_handles[0]),
287 "get sharing handle");
289 return sharing_handles;
292 auto get_parent_and_bit = [&](
const auto ent) {
295 m_field.
get_moab().tag_get_data(th_RefParentHandle, &ent, 1, &parent),
299 m_field.
get_moab().tag_get_data(th_RefBitLevel, &ent, 1, &
bit),
301 return std::make_pair(parent,
bit);
304 auto set_parent = [&](
const auto ent,
const auto parent) {
305 return m_field.
get_moab().tag_set_data(th_RefParentHandle, &ent, 1,
309 auto set_bit = [&](
const auto ent,
const auto bit) {
310 return m_field.
get_moab().tag_set_data(th_RefBitLevel, &ent, 1, &
bit);
314 std::vector<std::vector<unsigned long long>> sbuffer(m_field.
get_comm_size());
316 for (
auto ent : shared) {
318 auto pstatus = get_pstatus(ent);
319 auto sharing_procs = get_sharing_procs(ent, pstatus);
320 auto sharing_handles = get_sharing_handles(ent, pstatus);
321 auto [parent,
bit] = get_parent_and_bit(ent);
324 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"pstatus " << std::bitset<8>(pstatus);
327 auto pstatus_parent = get_pstatus(parent);
328 auto sharing_procs_parent = get_sharing_procs(parent, pstatus_parent);
329 auto sharing_handles_parent = get_sharing_handles(parent, pstatus_parent);
331 for (
int proc = 0; proc < MAX_SHARING_PROCS && -1 != sharing_procs[proc];
333 if (sharing_procs[proc] == -1)
335 "sharing processor not set");
339 auto it = std::find(sharing_procs_parent.begin(),
340 sharing_procs_parent.end(), sharing_procs[proc]);
341 if (it == sharing_procs_parent.end()) {
344 "Sharing proc for parent entity can not be found proc = %u",
345 sharing_procs[proc]);
348 auto handle_on_sharing_proc = sharing_handles[proc];
349 auto parent_handle_on_sharing_proc =
350 sharing_handles_parent[std::distance(sharing_procs_parent.begin(),
352 sbuffer[sharing_procs[proc]].push_back(handle_on_sharing_proc);
353 sbuffer[sharing_procs[proc]].push_back(parent_handle_on_sharing_proc);
355 sbuffer[sharing_procs[proc]].push_back(
bit.to_ullong());
356 }
catch (std::exception &ex) {
357 MOFEM_LOG(
"SELF", Sev::warning) << ex.what();
359 <<
"On " << ent <<
" "
361 MOFEM_LOG(
"SELF", Sev::warning) <<
"For bit ref " <<
bit;
364 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"send %lu (%lu) to %d at %d\n", ent,
365 handle_on_sharing_proc, sharing_procs[proc],
368 if (!(pstatus & PSTATUS_MULTISHARED))
373 for (
int proc = 0; proc < MAX_SHARING_PROCS && -1 != sharing_procs[proc];
375 if (sharing_procs[proc] == -1)
377 "sharing processor not set");
380 auto handle_on_sharing_proc = sharing_handles[proc];
381 sbuffer[sharing_procs[proc]].push_back(handle_on_sharing_proc);
382 sbuffer[sharing_procs[proc]].push_back(parent);
385 sbuffer[sharing_procs[proc]].push_back(
bit.to_ullong());
386 }
catch (std::exception &ex) {
387 MOFEM_LOG(
"SELF", Sev::warning) << ex.what();
389 <<
"On " << ent <<
" "
391 MOFEM_LOG(
"SELF", Sev::warning) <<
"For bit ref " <<
bit;
395 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"send %lu (%lu) to %d at %d\n", ent,
396 handle_on_sharing_proc, sharing_procs[proc],
399 if (!(pstatus & PSTATUS_MULTISHARED))
407 std::vector<int> sbuffer_lengths(
410 const size_t block_size =
sizeof(
unsigned long long) /
sizeof(
int);
413 if (!sbuffer[proc].empty()) {
415 sbuffer_lengths[proc] = sbuffer[proc].size() * block_size;
420 sbuffer_lengths[proc] = 0;
432 CHKERR PetscGatherNumberOfMessages(comm, NULL, &sbuffer_lengths[0], &nrecvs);
438 CHKERR PetscGatherMessageLengths(comm, nsends, nrecvs, &sbuffer_lengths[0],
446 CHKERR PetscCommGetNewTag(comm, &tag);
451 MPI_Request *r_waits;
455 CHKERR PetscPostIrecvInt(comm, tag, nrecvs, onodes, olengths, &rbuf,
458 MPI_Request *s_waits;
459 CHKERR PetscMalloc1(nsends, &s_waits);
462 for (
int proc = 0, kk = 0; proc < m_field.
get_comm_size(); proc++) {
463 if (!sbuffer_lengths[proc])
465 CHKERR MPI_Isend(&(sbuffer[proc])[0],
466 sbuffer_lengths[proc],
468 tag, comm, s_waits + kk);
474 CHKERR MPI_Waitall(nrecvs, r_waits, &status[0]);
478 CHKERR MPI_Waitall(nsends, s_waits, &status[0]);
482 "Rank %d nb. shared to synchronise parents ents %u\n",
487 for (
int kk = 0; kk < nrecvs; kk++) {
489 int len = olengths[kk];
490 int *data_from_proc = rbuf[kk];
492 for (
int ee = 0; ee < len;) {
495 unsigned long long uulong_bit;
498 bcopy(&data_from_proc[ee], &parent,
sizeof(
EntityHandle));
500 bcopy(&data_from_proc[ee], &uulong_bit,
sizeof(
unsigned long long));
503 CHKERR set_parent(ent, parent);
507 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"received %lu (%lu) from %d at %d\n",
515 CHKERR PetscFree(s_waits);
516 CHKERR PetscFree(rbuf[0]);
518 CHKERR PetscFree(r_waits);
520 CHKERR PetscFree(olengths);
521 CHKERR PetscCommDestroy(&comm);
530 const Problem *problem_ptr,
const std::string &fe_name,
int verb) {
533 ParallelComm *pcomm = ParallelComm::get_pcomm(
535 std::vector<int> shprocs(MAX_SHARING_PROCS, 0);
536 std::vector<EntityHandle> shhandles(MAX_SHARING_PROCS, 0);
538 Tag th_gid = m_field.
get_moab().globalId_tag();
543 CHKERR PetscLayoutGetRange(layout, &gid, &last_gid);
544 CHKERR PetscLayoutDestroy(&layout);
554 0, (*fe_miit)->getFEUId()));
558 get_id_for_max_type<MBENTITYSET>(), (*fe_miit)->getFEUId()));
559 for (; fit != hi_fe_it; ++fit) {
561 auto ent = (*fit)->getEnt();
562 auto part = (*fit)->getPart();
565 CHKERR m_field.
get_moab().tag_set_data(pcomm->part_tag(), &ent, 1, &part);
566 if (part == pcomm->rank()) {
573 if (pcomm->size() > 1) {
575 unsigned char pstatus = 0;
576 if (pcomm->rank() != part) {
577 pstatus = PSTATUS_NOT_OWNED;
578 pstatus |= PSTATUS_GHOST;
581 if (pcomm->size() > 2) {
582 pstatus |= PSTATUS_SHARED;
583 pstatus |= PSTATUS_MULTISHARED;
585 pstatus |= PSTATUS_SHARED;
589 for (
size_t rr = 0; rr < pcomm->size(); ++rr) {
590 if (rr != pcomm->rank()) {
591 shhandles[rrr] = ent;
596 for (; rrr != pcomm->size(); ++rrr)
599 if (pstatus & PSTATUS_SHARED) {
600 CHKERR m_field.
get_moab().tag_set_data(pcomm->sharedp_tag(), &ent, 1,
602 CHKERR m_field.
get_moab().tag_set_data(pcomm->sharedh_tag(), &ent, 1,
606 if (pstatus & PSTATUS_MULTISHARED) {
607 CHKERR m_field.
get_moab().tag_set_data(pcomm->sharedps_tag(), &ent, 1,
609 CHKERR m_field.
get_moab().tag_set_data(pcomm->sharedhs_tag(), &ent, 1,
612 CHKERR m_field.
get_moab().tag_set_data(pcomm->pstatus_tag(), &ent, 1,
618 CHKERR pcomm->exchange_tags(th_gid, ents);
623 const std::string name,
const std::string &fe_name,
int verb) {
634 const int num_entities,
635 const int owner_proc,
int verb) {
641 ParallelComm *pcomm = ParallelComm::get_pcomm(
644 Range all_ents_range;
645 all_ents_range.insert_list(entities, entities + num_entities);
647 auto get_tag = [&]() {
return m_field.
get_moab().globalId_tag(); };
649 auto delete_tag = [&](
auto &&th_gid) {
655 auto resolve_shared_ents = [&](
auto &&th_gid,
auto &all_ents_range) {
656 auto set_gid = [&](
auto &th_gid) {
657 std::vector<int> gids(num_entities);
658 for (
size_t g = 0;
g != all_ents_range.size(); ++
g)
666 auto get_skin_ents = [&](
auto &all_ents_range) {
667 std::array<Range, 4> proc_ents_skin;
668 proc_ents_skin[3] = all_ents_range.subset_by_dimension(3);
669 proc_ents_skin[2] = all_ents_range.subset_by_dimension(2);
670 proc_ents_skin[1] = all_ents_range.subset_by_dimension(1);
671 proc_ents_skin[0] = all_ents_range.subset_by_dimension(0);
672 return proc_ents_skin;
675 auto resolve_dim = [&](
auto &all_ents_range) {
676 for (
int resolve_dim = 3; resolve_dim >= 0; --resolve_dim) {
677 if (all_ents_range.num_of_dimension(resolve_dim))
683 auto get_proc_ent = [&](
auto &all_ents_range) {
686 proc_ent = all_ents_range;
690 auto resolve_shared_ents = [&](
auto &&proc_ents,
auto &&skin_ents) {
691 return pcomm->resolve_shared_ents(
692 0, proc_ents, resolve_dim(all_ents_range),
693 resolve_dim(all_ents_range), skin_ents.data(), set_gid(th_gid));
696 CHKERR resolve_shared_ents(get_proc_ent(all_ents_range),
697 get_skin_ents(all_ents_range));
702 CHKERR delete_tag(resolve_shared_ents(get_tag(), all_ents_range));
710 CHKERR pcomm->get_owner_handle(e, moab_owner_proc, moab_owner_handle);
712 unsigned char pstatus = 0;
714 CHKERR m_field.
get_moab().tag_get_data(pcomm->pstatus_tag(), &e, 1,
717 std::vector<int> shprocs(MAX_SHARING_PROCS, 0);
718 std::vector<EntityHandle> shhandles(MAX_SHARING_PROCS, 0);
720 CHKERR m_field.
get_moab().tag_get_data(pcomm->sharedp_tag(), &e, 1,
722 CHKERR m_field.
get_moab().tag_get_data(pcomm->sharedh_tag(), &e, 1,
724 if (pstatus & PSTATUS_MULTISHARED) {
725 CHKERR m_field.
get_moab().tag_get_data(pcomm->sharedps_tag(), &e, 1,
727 CHKERR m_field.
get_moab().tag_get_data(pcomm->sharedhs_tag(), &e, 1,
731 std::ostringstream ss;
734 if (!(pstatus & PSTATUS_NOT_OWNED))
736 if (pstatus & PSTATUS_SHARED)
737 ss <<
"PSTATUS_SHARED ";
738 if (pstatus & PSTATUS_MULTISHARED)
739 ss <<
"PSTATUS_MULTISHARED ";
741 ss <<
"owner " << moab_owner_proc <<
" (" << owner_proc <<
") ";
745 ss << shprocs[r] <<
" ";
749 ss << shhandles[r] <<
" ";
752 MOFEM_LOG(
"SYNC", Sev::noisy) << ss.str();
758 for (
auto e : all_ents_range)
767 const int owner_proc,
772 const int num_ents = entities.size();
773 std::vector<EntityHandle> vec_ents(num_ents);
774 std::copy(entities.begin(), entities.end(), vec_ents.begin());
783 const int owner_proc,
int verb) {
788 std::vector<EntityHandle> field_ents;
789 CHKERR m_field.
get_moab().get_entities_by_handle(field_meshset, field_ents,
803 Range exchange_ents_data_verts, exchange_ents_data;
812 for (
auto it = lo; it != hi; ++it)
815 ((*it)->getPStatus()) &&
817 (*it)->getNbDofsOnEnt()
820 if ((*it)->getEntType() == MBVERTEX)
821 exchange_ents_data_verts.insert((*it)->getEnt());
823 exchange_ents_data.insert((*it)->getEnt());
827 ParallelComm *pcomm = ParallelComm::get_pcomm(
830 auto exchange = [&](
const Range &ents, Tag
th) {
833 std::vector<Tag> tags;
835 CHKERR pcomm->exchange_tags(tags, tags, ents);
840 CHKERR exchange(exchange_ents_data_verts, field_ptr->th_FieldDataVerts);
841 CHKERR exchange(exchange_ents_data, field_ptr->th_FieldData);
848 const int adj_dim,
const int n_parts,
849 Tag *th_vertex_weights, Tag *th_edge_weights,
850 Tag *th_part_weights,
int verb,
const bool debug) {
855 int rstart, rend, nb_elems;
859 CHKERR PetscLayoutSetBlockSize(layout, 1);
860 CHKERR PetscLayoutSetSize(layout, ents.size());
861 CHKERR PetscLayoutSetUp(layout);
862 CHKERR PetscLayoutGetSize(layout, &nb_elems);
863 CHKERR PetscLayoutGetRange(layout, &rstart, &rend);
864 CHKERR PetscLayoutDestroy(&layout);
867 <<
"Finite elements in problem: row lower " << rstart <<
" row upper "
868 << rend <<
" nb. elems " << nb_elems <<
" ( " << ents.size() <<
" )";
873 std::vector<EntityHandle> weight_ents;
874 weight_ents.reserve(rend - rstart + 1);
878 std::vector<int> adj;
879 AdjBridge(
const EntityHandle ent, std::vector<int> &adj)
880 : ent(ent), adj(adj) {}
883 typedef multi_index_container<
887 hashed_unique<member<AdjBridge, EntityHandle, &AdjBridge::ent>>
892 auto get_it = [&](
auto i) {
893 auto it = ents.begin();
895 if (it == ents.end())
903 proc_ents.insert(get_it(rstart), get_it(rend));
904 if (proc_ents.size() != rend - rstart)
906 "Wrong number of elements in range %d != %d", proc_ents.size(),
911 proc_ents, adj_dim,
true, all_dim_ents, moab::Interface::UNION);
913 AdjBridgeMap adj_bridge_map;
914 auto hint = adj_bridge_map.begin();
915 std::vector<int> adj;
916 for (
auto ent : all_dim_ents) {
919 adj_ents = intersect(adj_ents, ents);
921 adj.reserve(adj_ents.size());
922 for (
auto a : adj_ents)
923 adj.emplace_back(ents.index(
a));
924 hint = adj_bridge_map.emplace_hint(hint, ent, adj);
930 const int nb_loc_elements = rend - rstart;
931 std::vector<int>
i(nb_loc_elements + 1, 0),
j;
933 std::vector<int> row_adj;
934 Range::iterator fe_it;
937 for (fe_it = proc_ents.begin(), ii = rstart, jj = 0, max_row_size = 0;
938 fe_it != proc_ents.end(); ++fe_it, ++ii) {
942 "not yet implemented, don't know what to do for meshset "
947 CHKERR m_field.
get_moab().get_adjacencies(&*fe_it, 1, adj_dim,
false,
949 dim_ents = intersect(dim_ents, all_dim_ents);
952 for (
auto e : dim_ents) {
953 auto adj_it = adj_bridge_map.find(e);
954 if (adj_it != adj_bridge_map.end()) {
956 for (
const auto idx : adj_it->adj)
957 row_adj.push_back(idx);
964 std::sort(row_adj.begin(), row_adj.end());
965 auto end = std::unique(row_adj.begin(), row_adj.end());
967 size_t row_size = std::distance(row_adj.begin(), end);
968 max_row_size = std::max(max_row_size, row_size);
969 if (
j.capacity() < (nb_loc_elements - jj) * max_row_size)
970 j.reserve(nb_loc_elements * max_row_size);
973 auto diag = ents.index(*fe_it);
974 for (
auto it = row_adj.begin(); it != end; ++it)
981 if (th_vertex_weights != NULL)
982 weight_ents.push_back(*fe_it);
988 CHKERR PetscMalloc(
i.size() *
sizeof(
int), &_i);
989 CHKERR PetscMalloc(
j.size() *
sizeof(
int), &_j);
990 copy(
i.begin(),
i.end(), _i);
991 copy(
j.begin(),
j.end(), _j);
995 int *vertex_weights = NULL;
996 if (th_vertex_weights != NULL) {
997 CHKERR PetscMalloc(weight_ents.size() *
sizeof(
int), &vertex_weights);
999 &*weight_ents.begin(),
1000 weight_ents.size(), vertex_weights);
1006 CHKERR MatCreateMPIAdj(m_field.
get_comm(), rend - rstart, nb_elems, _i, _j,
1008 CHKERR MatSetOption(Adj, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
1012 CHKERR MatConvert(Adj, MATMPIAIJ, MAT_INITIAL_MATRIX, &
A);
1013 CHKERR MatView(
A, PETSC_VIEWER_DRAW_WORLD);
1020 MOFEM_LOG(
"WORLD", Sev::verbose) <<
"Start";
1022 MatPartitioning part;
1026 CHKERR MatPartitioningSetAdjacency(part, Adj);
1027 CHKERR MatPartitioningSetFromOptions(part);
1028 CHKERR MatPartitioningSetNParts(part, n_parts);
1029 if (th_vertex_weights != NULL) {
1030 CHKERR MatPartitioningSetVertexWeights(part, vertex_weights);
1033 PetscObjectTypeCompare((PetscObject)part, MATPARTITIONINGPARMETIS, &same);
1036 CHKERR MatPartitioningApply_Parmetis_MoFEM(part, &is);
1039 CHKERR MatPartitioningApply(part, &is);
1042 MOFEM_LOG(
"WORLD", Sev::verbose) <<
"End";
1045 IS is_gather, is_num, is_gather_num;
1046 CHKERR ISAllGather(is, &is_gather);
1047 CHKERR ISPartitioningToNumbering(is, &is_num);
1048 CHKERR ISAllGather(is_num, &is_gather_num);
1050 const int *part_number, *gids;
1051 CHKERR ISGetIndices(is_gather, &part_number);
1052 CHKERR ISGetIndices(is_gather_num, &gids);
1055 ParallelComm *pcomm = ParallelComm::get_pcomm(
1057 Tag part_tag = pcomm->part_tag();
1058 CHKERR m_field.
get_moab().tag_set_data(part_tag, ents, part_number);
1059 Tag gid_tag = m_field.
get_moab().globalId_tag();
1061 std::map<int, Range> parts_ents;
1064 Range::iterator eit = ents.begin();
1065 for (
int ii = 0; eit != ents.end(); eit++, ii++) {
1066 parts_ents[part_number[ii]].insert(*eit);
1070 0, MBENTITYSET, &part_tag, NULL, 1, tagged_sets,
1071 moab::Interface::UNION);
1072 if (!tagged_sets.empty())
1073 CHKERR m_field.
get_moab().tag_delete_data(part_tag, tagged_sets);
1075 if (n_parts > (
int)tagged_sets.size()) {
1077 int num_new = n_parts - tagged_sets.size();
1078 for (
int i = 0;
i < num_new;
i++) {
1081 tagged_sets.insert(new_set);
1083 }
else if (n_parts < (
int)tagged_sets.size()) {
1085 int num_del = tagged_sets.size() - n_parts;
1086 for (
int i = 0;
i < num_del;
i++) {
1093 std::vector<int> dum_ids(n_parts);
1094 for (
int i = 0;
i < n_parts;
i++)
1101 for (
int pp = 0; pp != n_parts; pp++) {
1102 Range dim_ents = parts_ents[pp].subset_by_dimension(
dim);
1103 for (
int dd =
dim - 1; dd >= 0; dd--) {
1107 dim_ents, dd,
false, adj_ents, moab::Interface::UNION);
1112 parts_ents[pp].merge(adj_ents);
1115 for (
int pp = 1; pp != n_parts; pp++) {
1116 for (
int ppp = 0; ppp != pp; ppp++) {
1117 parts_ents[pp] = subtract(parts_ents[pp], parts_ents[ppp]);
1122 for (
int rr = 0; rr != n_parts; rr++) {
1124 ss <<
"out_part_" << rr <<
".vtk";
1126 <<
"Save debug part mesh " << ss.str();
1130 CHKERR m_field.
get_moab().write_file(ss.str().c_str(),
"VTK",
"",
1136 for (
int pp = 0; pp != n_parts; pp++) {
1137 CHKERR m_field.
get_moab().add_entities(tagged_sets[pp], parts_ents[pp]);
1147 for (
int pp = 0; pp != n_parts; pp++) {
1148 Range type_ents = parts_ents[pp].subset_by_type(
t);
1149 if (
t != MBVERTEX) {
1150 CHKERR m_field.
get_moab().tag_clear_data(part_tag, type_ents, &pp);
1153 auto eit = type_ents.begin();
1154 for (; eit != type_ents.end();) {
1155 CHKERR m_field.
get_moab().tag_iterate(gid_tag, eit, type_ents.end(),
1157 auto gid_tag_ptr =
static_cast<int *
>(ptr);
1158 for (; count > 0; --count) {
1169 CHKERR ISRestoreIndices(is_gather, &part_number);
1170 CHKERR ISRestoreIndices(is_gather_num, &gids);
1171 CHKERR ISDestroy(&is_num);
1172 CHKERR ISDestroy(&is_gather_num);
1173 CHKERR ISDestroy(&is_gather);
1175 CHKERR MatPartitioningDestroy(&part);
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
#define MOFEM_LOG_C(channel, severity, format,...)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
virtual const FieldEntity_multiIndex * get_field_ents() const =0
Get the field ents object.
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
virtual const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const =0
get field structure
MoFEMErrorCode synchroniseFieldEntities(const std::string name, int verb=DEFAULT_VERBOSITY)
#define MOFEM_LOG(channel, severity)
Log.
MoFEMErrorCode partitionMesh(const Range &ents, const int dim, const int adj_dim, const int n_parts, Tag *th_vertex_weights=nullptr, Tag *th_edge_weights=nullptr, Tag *th_part_weights=nullptr, int verb=VERBOSE, const bool debug=false)
Set partition tag to each finite element in the problem.
MoFEMErrorCode resolveSharedFiniteElements(const Problem *problem_ptr, const std::string &fe_name, int verb=DEFAULT_VERBOSITY)
resolve shared entities for finite elements in the problem
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
implementation of Data Operators for Forces and Sources
auto type_from_handle(const EntityHandle h)
get type from entity handle
constexpr double t
plate stiffness
constexpr auto field_name
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
CommInterface(const MoFEM::Core &core)
MoFEMErrorCode makeEntitiesMultishared(const EntityHandle *entities, const int num_entities, const int owner_proc=0, int verb=DEFAULT_VERBOSITY)
make entities from proc 0 shared on all proc
MoFEMErrorCode synchroniseEntities(Range &ent, int verb=DEFAULT_VERBOSITY)
MoFEMErrorCode makeFieldEntitiesMultishared(const std::string field_name, const int owner_proc=0, int verb=DEFAULT_VERBOSITY)
make field entities multi shared
MoFEMErrorCode resolveParentEntities(const Range &ent, int verb=DEFAULT_VERBOSITY)
Synchronise parent entities.
MoFEMErrorCode exchangeFieldData(const std::string field_name, int verb=DEFAULT_VERBOSITY)
Exchange field data.
virtual int get_comm_size() const =0
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
virtual moab::Interface & get_moab()=0
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
virtual boost::shared_ptr< BasicEntityData > & get_basic_entity_data_ptr()=0
Get pointer to basic entity data.
Tag get_th_RefBitLevel() const
Tag get_th_RefParentHandle() const
Deprecated interface functions.
UId getLocalUniqueIdCalculate() const
Generate UId for finite element entity.
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
keeps basic data about problem
MoFEMErrorCode getNumberOfElementsByNameAndPart(MPI_Comm comm, const std::string name, PetscLayout *layout) const
Get number of finite elements by name on processors.
boost::shared_ptr< NumeredEntFiniteElement_multiIndex > numeredFiniteElementsPtr
store finite elements
base class for all interface classes