22#define DOCTEST_CONFIG_IMPLEMENT
56 {
"UniformSquare_10",
"UniformSquare_10.cgns", 2,
false, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, 100, 40},
57 {
"IV10_10",
"IV10_10.cgns", 2,
true, {10, 0, 0}, {0, 10, 0}, {0, 0, 10}, 100, -1},
58 {
"NACA0012_H2",
"NACA0012_H2.cgns", 2,
false, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, 20816, 484},
59 {
"IV10U_10",
"IV10U_10.cgns", 2,
true, {10, 0, 0}, {0, 10, 0}, {0, 0, 10}, 322, -1},
60 {
"Ball2",
"Ball2.cgns", 3,
false, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, 958994, 16402},
62static constexpr int N_CONFIGS =
sizeof(g_configs) /
sizeof(g_configs[0]);
76static std::string meshPath(
const std::string &name)
78 std::string f(__FILE__);
79 for (
int i = 0;
i < 4;
i++)
81 auto pos = f.rfind(
'/');
82 if (pos == std::string::npos)
84 if (pos != std::string::npos)
87 return f +
"/data/mesh/" + name;
91 int meshId,
const MeshConfig &cfg,
int nReorderParts = 1)
94 std::cout <<
"[buildMesh " << meshId <<
" " << cfg.
name <<
"] START" << std::endl;
96 auto mesh = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
102 mesh->SetPeriodicGeometry(
108 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
109 reader.Deduplicate1to1Periodic(1e-8);
110 reader.BuildCell2Cell();
117 reader.MeshPartitionCell2Cell(pOpt);
118 reader.PartitionReorderToMeshCell2Cell();
120 mesh->RecoverNode2CellAndNode2Bnd();
121 mesh->RecoverCell2CellAndBnd2Cell();
122 mesh->BuildGhostPrimary();
124 mesh->AdjGlobal2LocalPrimary();
125 mesh->AdjGlobal2LocalN2CB();
127 if (nReorderParts > 1)
128 mesh->ReorderLocalCells(nReorderParts);
130 mesh->InterpolateFace();
131 mesh->AssertOnFaces();
133 mesh->AdjLocal2GlobalN2CB();
134 mesh->BuildGhostN2CB();
135 mesh->AdjGlobal2LocalN2CB();
136 mesh->AssertOnN2CB();
138 mesh->BuildCell2CellFace();
139 mesh->AdjGlobal2LocalC2CFace();
142 mesh->RecreatePeriodicNodes();
143 mesh->BuildVTKConnectivity();
146 std::cout <<
"[buildMesh " << meshId <<
" " << cfg.
name <<
"] DONE" << std::endl;
151static std::vector<std::vector<DNDS::index>> snapshotAdj(
154 std::vector<std::vector<DNDS::index>> out(nRows);
157 out[
i].resize(
adj.RowSize(
i));
164static std::vector<std::pair<DNDS::index, DNDS::index>> snapshotAdj2(
167 std::vector<std::pair<DNDS::index, DNDS::index>> out(nRows);
173static std::vector<DNDS::index> snapshotAdj1(
176 std::vector<DNDS::index> out(nRows);
182static void checkAdjMatchesSnapshot(
184 const std::vector<std::vector<DNDS::index>> &snap)
197 MPI_Init(&argc, &argv);
201 for (
int ic = 0; ic < N_CONFIGS; ic++)
206 if (ic == 4 && g_mpi.
size != 8)
208 g_full[ic] = buildFullMesh(ic, g_configs[ic]);
211 doctest::Context ctx;
212 ctx.applyCommandLine(argc, argv);
215 for (
auto &
m : g_full)
225#define FOR_EACH_MESH_CONFIG(body) \
226 for (int _ci = 0; _ci < N_CONFIGS; _ci++) \
229 if (_ci == 4 && g_mpi.size != 8) \
232 const auto &cfg = g_configs[_ci]; \
233 auto m = g_full[_ci]; \
234 SUBCASE(cfg.name) { body } \
250TEST_CASE(
"Pipeline: global cell count matches expected"){
255 CHECK(globalCells > 0);
260TEST_CASE(
"Pipeline: global bnd count matches expected"){
268 CHECK(globalBnds > 0);
278TEST_CASE(
"InterpolateFace: face2cell has exactly 2 entries per face"){
281 CHECK(
m->face2cell.RowSize(iF) == 2);
284TEST_CASE(
"InterpolateFace: face2cell owner is a valid local cell"){
291 CHECK(owner < totalCells);
295TEST_CASE(
"InterpolateFace: face2node indices in valid range"){
303 CHECK(iN < totalNodes);
307TEST_CASE(
"InterpolateFace: cell2face row sizes match expected face count"){
311 auto elem =
m->GetCellElement(iC);
312 int nFaceExpected = elem.GetNumFaces();
313 CHECK(
m->cell2face.RowSize(iC) == nFaceExpected);
317TEST_CASE(
"InterpolateFace: cell2face entries are valid face indices"){
325 CHECK(iF < totalFaces);
329TEST_CASE(
"InterpolateFace: bnd2face maps to valid faces"){
336 CHECK(iFace < m->NumFace() +
m->NumFaceGhost());
340TEST_CASE(
"InterpolateFace: face element types are valid")
344 CHECK(
m->GetFaceElement(iF).type != Elem::ElemType::UnknownElem);
355 auto mesh = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
361 mesh->SetPeriodicGeometry(
367 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
368 reader.Deduplicate1to1Periodic(1e-8);
369 reader.BuildCell2Cell();
376 reader.MeshPartitionCell2Cell(pOpt);
377 reader.PartitionReorderToMeshCell2Cell();
379 mesh->RecoverNode2CellAndNode2Bnd();
380 mesh->RecoverCell2CellAndBnd2Cell();
381 mesh->BuildGhostPrimary();
382 mesh->AdjGlobal2LocalPrimary();
383 mesh->AdjGlobal2LocalN2CB();
386 mesh->InterpolateFaceLegacy();
388 mesh->InterpolateFace();
394static std::vector<std::set<std::set<DNDS::index>>>
398 std::vector<std::set<std::set<DNDS::index>>>
result(nCells);
399 for (
DNDS::index iCell = 0; iCell < nCells; iCell++)
401 auto eCell =
Elem::Element{
m->cellElemInfo[iCell]->getElemType()};
407 std::set<DNDS::index> vs;
408 for (
int k = 0; k < nVerts; k++)
409 vs.insert(
m->face2node[iFace][k]);
416TEST_CASE(
"InterpolateFace: DSL matches Legacy on all mesh configs")
418 for (
int ci = 0; ci < N_CONFIGS; ci++)
420 if (ci == 4 && g_mpi.
size != 8)
423 const auto &cfg = g_configs[ci];
426 auto mDSL = buildMeshUpToFace(cfg,
false);
427 auto mLeg = buildMeshUpToFace(cfg,
true);
430 CHECK(mDSL->NumFace() == mLeg->NumFace());
431 CHECK(mDSL->NumFaceGhost() == mLeg->NumFaceGhost());
434 auto dslSets = collectCellFaceVertexSets(mDSL);
435 auto legSets = collectCellFaceVertexSets(mLeg);
436 DNDS::index nCells = mDSL->cell2node.father->Size();
439 for (
DNDS::index iCell = 0; iCell < nCells; iCell++)
442 CHECK(dslSets[iCell] == legSets[iCell]);
448 std::map<std::set<DNDS::index>, std::pair<DNDS::index, DNDS::index>> dslF2C, legF2C;
449 for (
DNDS::index iF = 0; iF < mDSL->NumFace(); iF++)
451 auto eFace =
Elem::Element{mDSL->faceElemInfo(iF, 0).getElemType()};
452 std::set<DNDS::index> vs;
453 for (
int k = 0; k < eFace.GetNumVertices(); k++)
454 vs.insert(mDSL->face2node[iF][k]);
455 dslF2C[vs] = {mDSL->face2cell(iF, 0), mDSL->face2cell(iF, 1)};
457 for (
DNDS::index iF = 0; iF < mLeg->NumFace(); iF++)
459 auto eFace =
Elem::Element{mLeg->faceElemInfo(iF, 0).getElemType()};
460 std::set<DNDS::index> vs;
461 for (
int k = 0; k < eFace.GetNumVertices(); k++)
462 vs.insert(mLeg->face2node[iF][k]);
463 legF2C[vs] = {mLeg->face2cell(iF, 0), mLeg->face2cell(iF, 1)};
465 CHECK(dslF2C.size() == legF2C.size());
467 for (
auto &[vs, dslPair] : dslF2C)
469 auto it = legF2C.find(vs);
472 std::set<DNDS::index> dslCells{dslPair.first, dslPair.second};
473 std::set<DNDS::index> legCells{it->second.first, it->second.second};
474 CHECK(dslCells == legCells);
478 std::map<std::set<DNDS::index>,
t_index> dslZones, legZones;
479 for (
DNDS::index iF = 0; iF < mDSL->NumFace(); iF++)
481 auto eFace =
Elem::Element{mDSL->faceElemInfo(iF, 0).getElemType()};
482 std::set<DNDS::index> vs;
483 for (
int k = 0; k < eFace.GetNumVertices(); k++)
484 vs.insert(mDSL->face2node[iF][k]);
485 dslZones[vs] = mDSL->faceElemInfo(iF, 0).zone;
487 for (
DNDS::index iF = 0; iF < mLeg->NumFace(); iF++)
489 auto eFace =
Elem::Element{mLeg->faceElemInfo(iF, 0).getElemType()};
490 std::set<DNDS::index> vs;
491 for (
int k = 0; k < eFace.GetNumVertices(); k++)
492 vs.insert(mLeg->face2node[iF][k]);
493 legZones[vs] = mLeg->faceElemInfo(iF, 0).zone;
495 for (
auto &[vs, dslZone] : dslZones)
497 auto it = legZones.find(vs);
499 CHECK(dslZone == it->second);
509TEST_CASE(
"Periodic: face pbi has uniform XOR between owner and neighbor (collaborating check)")
511 for (
int _ci = 0; _ci < N_CONFIGS; _ci++)
513 if (_ci == 4 && g_mpi.
size != 8)
515 const auto &cfg = g_configs[_ci];
518 auto m = g_full[_ci];
522 for (
DNDS::index iFace = 0; iFace <
m->NumFace(); iFace++)
529 auto eCellL =
Elem::Element{
m->cellElemInfo[iCellL]->getElemType()};
530 auto eCellR =
Elem::Element{
m->cellElemInfo[iCellR]->getElemType()};
532 int slotL = -1, slotR = -1;
534 if (
m->cell2face(iCellL,
j) == iFace)
540 if (
m->cell2face(iCellR,
j) == iFace)
548 auto eFaceL = eCellL.ObtainFace(slotL);
549 int nFN = eFaceL.GetNumNodes();
551 std::vector<DNDS::index> nodesL(nFN), nodesR(nFN);
552 std::vector<NodePeriodicBits> pbiL(nFN), pbiR(nFN);
553 eCellL.ExtractFaceNodes(slotL,
m->cell2node[iCellL], nodesL);
554 eCellL.ExtractFaceNodes(slotL,
m->cell2nodePbi[iCellL], pbiL);
555 eCellR.ExtractFaceNodes(slotR,
m->cell2node[iCellR], nodesR);
556 eCellR.ExtractFaceNodes(slotR,
m->cell2nodePbi[iCellR], pbiR);
558 std::vector<idx_pbi_t> pairsL(nFN), pairsR(nFN);
559 for (
int k = 0; k < nFN; k++)
561 pairsL[k] = {nodesL[k], uint8_t(pbiL[k])};
562 pairsR[k] = {nodesR[k], uint8_t(pbiR[k])};
564 std::sort(pairsL.begin(), pairsL.end());
565 std::sort(pairsR.begin(), pairsR.end());
567 for (
int k = 0; k < nFN; k++)
571 CHECK(pairsL[k].first == pairsR[k].first);
574 uint8_t v0 = pairsL[0].second ^ pairsR[0].second;
575 for (
int k = 1; k < nFN; k++)
577 uint8_t vk = pairsL[k].second ^ pairsR[k].second;
589TEST_CASE(
"Periodic: face2nodePbi is consistent with owner cell's cell2nodePbi")
591 for (
int _ci = 0; _ci < N_CONFIGS; _ci++)
593 if (_ci == 4 && g_mpi.
size != 8)
595 const auto &cfg = g_configs[_ci];
598 auto m = g_full[_ci];
602 for (
DNDS::index iFace = 0; iFace <
m->NumFace(); iFace++)
605 auto eCell =
Elem::Element{
m->cellElemInfo[iCellOwner]->getElemType()};
609 if (
m->cell2face(iCellOwner,
j) == iFace)
616 auto eFace = eCell.ObtainFace(iSlot);
617 int nFN = eFace.GetNumNodes();
619 std::vector<DNDS::index> cellFaceNodes(nFN);
620 std::vector<NodePeriodicBits> cellFacePbi(nFN);
621 eCell.ExtractFaceNodes(iSlot,
m->cell2node[iCellOwner], cellFaceNodes);
622 eCell.ExtractFaceNodes(iSlot,
m->cell2nodePbi[iCellOwner], cellFacePbi);
624 std::set<idx_pbi_t> fromCell, fromFace;
625 for (
int k = 0; k < nFN; k++)
627 fromCell.insert({cellFaceNodes[k], uint8_t(cellFacePbi[k])});
628 fromFace.insert({
m->face2node(iFace, k), uint8_t(
m->face2nodePbi(iFace, k))});
631 CHECK(fromCell == fromFace);
637TEST_CASE(
"Periodic: RecreatePeriodicNodes produces consistent counts DSL vs Legacy")
639 for (
int ci = 0; ci < N_CONFIGS; ci++)
641 if (ci == 4 && g_mpi.
size != 8)
643 const auto &cfg = g_configs[ci];
649 auto mDSL = buildMeshUpToFace(cfg,
false);
650 auto mLeg = buildMeshUpToFace(cfg,
true);
652 mDSL->AdjLocal2GlobalN2CB();
653 mDSL->BuildGhostN2CB();
654 mDSL->AdjGlobal2LocalN2CB();
655 mDSL->RecreatePeriodicNodes();
657 mLeg->AdjLocal2GlobalN2CB();
658 mLeg->BuildGhostN2CB();
659 mLeg->AdjGlobal2LocalN2CB();
660 mLeg->RecreatePeriodicNodes();
662 CHECK(mDSL->coordsPeriodicRecreated.father->Size() ==
663 mLeg->coordsPeriodicRecreated.father->Size());
665 auto dslGlobal = mDSL->coordsPeriodicRecreated.father->globalSize();
666 auto legGlobal = mLeg->coordsPeriodicRecreated.father->globalSize();
667 CHECK(dslGlobal == legGlobal);
672TEST_CASE(
"Periodic: face2nodePbi matches DSL vs Legacy")
674 for (
int ci = 0; ci < N_CONFIGS; ci++)
676 if (ci == 4 && g_mpi.
size != 8)
678 const auto &cfg = g_configs[ci];
684 auto mDSL = buildMeshUpToFace(cfg,
false);
685 auto mLeg = buildMeshUpToFace(cfg,
true);
687 REQUIRE(mDSL->NumFace() == mLeg->NumFace());
690 using FaceKey = std::set<DNDS::index>;
691 using NodePbiSet = std::set<std::pair<DNDS::index, uint8_t>>;
692 std::map<FaceKey, NodePbiSet> dslMap, legMap;
694 for (
DNDS::index iF = 0; iF < mDSL->NumFace(); iF++)
696 auto eF =
Elem::Element{mDSL->faceElemInfo(iF, 0).getElemType()};
698 for (
int k = 0; k < eF.GetNumVertices(); k++)
699 vs.insert(mDSL->face2node(iF, k));
701 for (
int k = 0; k < eF.GetNumNodes(); k++)
702 nps.insert({mDSL->face2node(iF, k), uint8_t(mDSL->face2nodePbi(iF, k))});
705 for (
DNDS::index iF = 0; iF < mLeg->NumFace(); iF++)
707 auto eF =
Elem::Element{mLeg->faceElemInfo(iF, 0).getElemType()};
709 for (
int k = 0; k < eF.GetNumVertices(); k++)
710 vs.insert(mLeg->face2node(iF, k));
712 for (
int k = 0; k < eF.GetNumNodes(); k++)
713 nps.insert({mLeg->face2node(iF, k), uint8_t(mLeg->face2nodePbi(iF, k))});
717 CHECK(dslMap.size() == legMap.size());
718 for (
auto &[vs, dslNps] : dslMap)
720 auto it = legMap.find(vs);
722 CHECK(dslNps == it->second);
734 auto mesh = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
740 mesh->SetPeriodicGeometry(
746 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
747 reader.Deduplicate1to1Periodic(1e-8);
748 reader.BuildCell2Cell();
755 reader.MeshPartitionCell2Cell(pOpt);
756 reader.PartitionReorderToMeshCell2Cell();
760 mesh->RecoverNode2CellAndNode2BndLegacy();
761 mesh->RecoverCell2CellAndBnd2CellLegacy();
765 mesh->RecoverNode2CellAndNode2Bnd();
766 mesh->RecoverCell2CellAndBnd2Cell();
772TEST_CASE(
"bnd2cell: DSL matches Legacy on all mesh configs")
774 for (
int ci = 0; ci < N_CONFIGS; ci++)
776 if (ci == 4 && g_mpi.
size != 8)
778 const auto &cfg = g_configs[ci];
782 auto mDSL = buildMeshUpToBnd2Cell(cfg,
false);
783 auto mLeg = buildMeshUpToBnd2Cell(cfg,
true);
785 REQUIRE(mDSL->NumBnd() == mLeg->NumBnd());
788 for (
DNDS::index iBnd = 0; iBnd < mDSL->NumBnd(); iBnd++)
791 std::set<DNDS::index> dslCells{mDSL->bnd2cell(iBnd, 0), mDSL->bnd2cell(iBnd, 1)};
792 std::set<DNDS::index> legCells{mLeg->bnd2cell(iBnd, 0), mLeg->bnd2cell(iBnd, 1)};
793 CHECK(dslCells == legCells);
799 CHECK(mDSL->bnd2cell(iBnd, 0) == mLeg->bnd2cell(iBnd, 0));
800 CHECK(mDSL->bnd2cell(iBnd, 1) == mLeg->bnd2cell(iBnd, 1));
805 REQUIRE(mDSL->NumCell() == mLeg->NumCell());
806 for (
DNDS::index iCell = 0; iCell < mDSL->NumCell(); iCell++)
809 std::set<DNDS::index> dslNeighbors, legNeighbors;
811 dslNeighbors.insert(mDSL->cell2cell.father->operator()(iCell,
j));
813 legNeighbors.insert(mLeg->cell2cell.father->operator()(iCell,
j));
814 CHECK(dslNeighbors == legNeighbors);
826 auto mesh = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
832 mesh->SetPeriodicGeometry(
838 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
839 reader.Deduplicate1to1Periodic(1e-8);
840 reader.BuildCell2Cell();
847 reader.MeshPartitionCell2Cell(pOpt);
848 reader.PartitionReorderToMeshCell2Cell();
850 mesh->RecoverNode2CellAndNode2Bnd();
851 mesh->RecoverCell2CellAndBnd2Cell();
854 mesh->BuildGhostPrimaryLegacy();
856 mesh->BuildGhostPrimary();
861TEST_CASE(
"BuildGhostPrimary: DSL matches Legacy ghost sets")
863 for (
int ci = 0; ci < N_CONFIGS; ci++)
865 if (ci == 4 && g_mpi.
size != 8)
867 const auto &cfg = g_configs[ci];
871 auto mDSL = buildMeshUpToGhost(cfg,
false);
872 auto mLeg = buildMeshUpToGhost(cfg,
true);
877 CHECK(nCellDSL == nCellLeg);
881 auto &dslGhost = mDSL->cell2cell.trans.pLGhostMapping->ghostIndex;
882 auto &legGhost = mLeg->cell2cell.trans.pLGhostMapping->ghostIndex;
883 std::set<DNDS::index> dslSet(dslGhost.begin(), dslGhost.end());
884 std::set<DNDS::index> legSet(legGhost.begin(), legGhost.end());
885 CHECK(dslSet == legSet);
893 CHECK(nNodeDSL >= nNodeLeg);
897 auto &dslGhost = mDSL->coords.trans.pLGhostMapping->ghostIndex;
898 auto &legGhost = mLeg->coords.trans.pLGhostMapping->ghostIndex;
899 std::set<DNDS::index> dslSet(dslGhost.begin(), dslGhost.end());
900 std::set<DNDS::index> legSet(legGhost.begin(), legGhost.end());
901 CHECK(std::includes(dslSet.begin(), dslSet.end(),
902 legSet.begin(), legSet.end()));
911 CHECK(nBndDSL >= nBndLeg);
915 auto &dslGhost = mDSL->bnd2cell.trans.pLGhostMapping->ghostIndex;
916 auto &legGhost = mLeg->bnd2cell.trans.pLGhostMapping->ghostIndex;
917 std::set<DNDS::index> dslSet(dslGhost.begin(), dslGhost.end());
918 std::set<DNDS::index> legSet(legGhost.begin(), legGhost.end());
919 CHECK(std::includes(dslSet.begin(), dslSet.end(),
920 legSet.begin(), legSet.end()));
925 fmt::print(
" [{}] ghost: cells={}, nodes={}, bnds={}\n",
926 cfg.
name, nCellDSL, nNodeDSL, nBndDSL);
928 if (ci == 0 || ci == 2)
930 MPI_Barrier(g_mpi.
comm);
931 fmt::print(
" [{}] r{} DSL nBndGhost={} LEG nBndGhost={}\n",
932 cfg.
name, g_mpi.
rank, nBndDSL, nBndLeg);
934 for (
DNDS::index iBnd = mDSL->bnd2node.father->Size(); iBnd < mDSL->bnd2node.Size(); iBnd++)
936 DNDS::index iSon = iBnd - mDSL->bnd2node.father->Size();
937 DNDS::index bndGlobal = mDSL->bnd2node.trans.pLGhostMapping->ghostIndex[iSon];
938 fmt::print(
" r{} gBnd[{}]:", g_mpi.
rank, bndGlobal);
940 fmt::print(
" {}", mDSL->bnd2node(iBnd,
j));
954 auto mesh = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
960 mesh->SetPeriodicGeometry(
966 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
967 reader.Deduplicate1to1Periodic(1e-8);
968 reader.BuildCell2Cell();
975 reader.MeshPartitionCell2Cell(pOpt);
976 reader.PartitionReorderToMeshCell2Cell();
980TEST_CASE(
"Benchmark: DSL vs Legacy pipeline steps")
985 const auto &cfg = g_configs[2];
994 double tDSL = 0, tLeg = 0;
995 for (
int r = 0;
r < nRuns;
r++)
997 auto m1 = buildMeshForBenchmark(cfg);
998 MPI_Barrier(g_mpi.
comm);
999 double t0 = MPI_Wtime();
1000 m1->RecoverNode2CellAndNode2Bnd();
1001 MPI_Barrier(g_mpi.
comm);
1002 double t1 = MPI_Wtime();
1005 auto m2 = buildMeshForBenchmark(cfg);
1006 MPI_Barrier(g_mpi.
comm);
1007 double t2 = MPI_Wtime();
1008 m2->RecoverNode2CellAndNode2BndLegacy();
1009 MPI_Barrier(g_mpi.
comm);
1010 double t3 = MPI_Wtime();
1013 if (g_mpi.
rank == 0)
1014 fmt::print(
" Inverse (RecoverN2CB): DSL {:.4f}s Legacy {:.4f}s ({:.2f}x, {} runs)\n",
1015 tDSL / nRuns, tLeg / nRuns, tLeg / tDSL, nRuns);
1020 double tDSL = 0, tLeg = 0;
1021 for (
int r = 0;
r < nRuns;
r++)
1023 auto m1 = buildMeshForBenchmark(cfg);
1024 m1->RecoverNode2CellAndNode2Bnd();
1025 MPI_Barrier(g_mpi.
comm);
1026 double t0 = MPI_Wtime();
1027 m1->RecoverCell2CellAndBnd2Cell();
1028 MPI_Barrier(g_mpi.
comm);
1029 double t1 = MPI_Wtime();
1032 auto m2 = buildMeshForBenchmark(cfg);
1033 m2->RecoverNode2CellAndNode2BndLegacy();
1034 MPI_Barrier(g_mpi.
comm);
1035 double t2 = MPI_Wtime();
1036 m2->RecoverCell2CellAndBnd2CellLegacy();
1037 MPI_Barrier(g_mpi.
comm);
1038 double t3 = MPI_Wtime();
1041 if (g_mpi.
rank == 0)
1042 fmt::print(
" Compose (RecoverC2C): DSL {:.4f}s Legacy {:.4f}s ({:.2f}x, {} runs)\n",
1043 tDSL / nRuns, tLeg / nRuns, tLeg / tDSL, nRuns);
1048 double tDSL = 0, tLeg = 0;
1049 for (
int r = 0;
r < nRuns;
r++)
1051 auto m1 = buildMeshForBenchmark(cfg);
1052 m1->RecoverNode2CellAndNode2Bnd();
1053 m1->RecoverCell2CellAndBnd2Cell();
1054 MPI_Barrier(g_mpi.
comm);
1055 double t0 = MPI_Wtime();
1056 m1->BuildGhostPrimary();
1057 MPI_Barrier(g_mpi.
comm);
1058 double t1 = MPI_Wtime();
1061 auto m2 = buildMeshForBenchmark(cfg);
1062 m2->RecoverNode2CellAndNode2Bnd();
1063 m2->RecoverCell2CellAndBnd2Cell();
1064 MPI_Barrier(g_mpi.
comm);
1065 double t2 = MPI_Wtime();
1066 m2->BuildGhostPrimaryLegacy();
1067 MPI_Barrier(g_mpi.
comm);
1068 double t3 = MPI_Wtime();
1071 if (g_mpi.
rank == 0)
1072 fmt::print(
" Ghost (BuildGhostPri): DSL {:.4f}s Legacy {:.4f}s ({:.2f}x, {} runs)\n",
1073 tDSL / nRuns, tLeg / nRuns, tLeg / tDSL, nRuns);
1078 double tDSL = 0, tLeg = 0;
1079 for (
int r = 0;
r < nRuns;
r++)
1081 auto m1 = buildMeshForBenchmark(cfg);
1082 m1->RecoverNode2CellAndNode2Bnd();
1083 m1->RecoverCell2CellAndBnd2Cell();
1084 m1->BuildGhostPrimary();
1085 m1->AdjGlobal2LocalPrimary();
1086 m1->AdjGlobal2LocalN2CB();
1087 MPI_Barrier(g_mpi.
comm);
1088 double t0 = MPI_Wtime();
1089 m1->InterpolateFace();
1090 MPI_Barrier(g_mpi.
comm);
1091 double t1 = MPI_Wtime();
1094 auto m2 = buildMeshForBenchmark(cfg);
1095 m2->RecoverNode2CellAndNode2Bnd();
1096 m2->RecoverCell2CellAndBnd2Cell();
1097 m2->BuildGhostPrimary();
1098 m2->AdjGlobal2LocalPrimary();
1099 m2->AdjGlobal2LocalN2CB();
1100 MPI_Barrier(g_mpi.
comm);
1101 double t2 = MPI_Wtime();
1102 m2->InterpolateFaceLegacy();
1103 MPI_Barrier(g_mpi.
comm);
1104 double t3 = MPI_Wtime();
1107 if (g_mpi.
rank == 0)
1108 fmt::print(
" Interpolate (Face): DSL {:.4f}s Legacy {:.4f}s ({:.2f}x, {} runs)\n",
1109 tDSL / nRuns, tLeg / nRuns, tLeg / tDSL, nRuns);
1115TEST_CASE(
"N2CB: every local node has at least one adjacent cell"){
1117 for (
DNDS::index iNode = 0; iNode <
m->NumNode(); iNode++)
1119 bool hasCell =
false;
1121 if (
m->node2cell(iNode,
j) >= 0)
1127TEST_CASE(
"N2CB: node2cell entries are valid local indices"){
1130 for (
DNDS::index iNode = 0; iNode <
m->NumNodeProc(); iNode++)
1135 CHECK(iC < totalCells);
1141TEST_CASE(
"BuildCell2CellFace: row sizes match cell2face"){
1144 CHECK(
m->cell2cellFace.RowSize(iC) ==
m->cell2face.RowSize(iC));
1147TEST_CASE(
"BuildCell2CellFace: entries are valid cells or negative"){
1155 CHECK(nb < totalCells);
1161TEST_CASE(
"ConstructBndMesh: correct dimensions and state"){
1166 m->ConstructBndMesh(bMesh);
1167 CHECK(bMesh.dim ==
m->dim - 1);
1168 CHECK(bMesh.NumCell() ==
m->NumBnd());
1172TEST_CASE(
"ConstructBndMesh: node2parentNode valid"){
1177 m->ConstructBndMesh(bMesh);
1178 CHECK(bMesh.node2parentNode.size() ==
1179 static_cast<size_t>(bMesh.NumNode() + bMesh.NumNodeGhost()));
1184 CHECK(pn < m->NumNode() +
m->NumNodeGhost());
1190TEST_CASE(
"BuildVTKConnectivity: arrays populated correctly"){
1192 CHECK(
m->vtkCell2node.size() > 0);
1193 CHECK(
m->vtkCell2nodeOffsets.size() ==
static_cast<size_t>(
m->NumCell() + 1));
1194 CHECK(
m->vtkCellType.size() ==
static_cast<size_t>(
m->NumCell()));
1195 for (
size_t i = 1;
i <
m->vtkCell2nodeOffsets.size();
i++)
1196 CHECK(
m->vtkCell2nodeOffsets[
i] >=
m->vtkCell2nodeOffsets[
i - 1]);
1197 for (
auto ct :
m->vtkCellType)
1205TEST_CASE(
"ReorderLocalCells: all states preserved"){
1207 m->ReorderLocalCells(2);
1215TEST_CASE(
"ReorderLocalCells: cell count preserved")
1222 m->ReorderLocalCells(2);
1228TEST_CASE(
"ReorderLocalCells: partition starts valid"){
1230 m->ReorderLocalCells(2);
1231 int nParts =
m->NLocalParts();
1233 CHECK(
m->LocalPartStart(0) == 0);
1234 CHECK(
m->LocalPartEnd(nParts - 1) ==
m->NumCell());
1235 for (
int ip = 0; ip < nParts; ip++)
1236 CHECK(
m->LocalPartStart(ip) <=
m->LocalPartEnd(ip));
1239TEST_CASE(
"ReorderLocalCells: cell2node still valid"){
1241 m->ReorderLocalCells(2);
1248 CHECK(iN < totalNodes);
1252TEST_CASE(
"ReorderLocalCells: face count preserved")
1259 m->ReorderLocalCells(2);
1269TEST_CASE(
"AdjFacial round-trip: face2node and face2cell identity")
1274 auto f2nSnap = snapshotAdj(
m->face2node,
m->NumFace());
1275 auto f2cSnap = snapshotAdj2(
m->face2cell,
m->NumFace());
1277 m->AdjLocal2GlobalFacial();
1281 CHECK(
m->face2node(iF,
j) >= 0);
1283 m->AdjGlobal2LocalFacial();
1286 checkAdjMatchesSnapshot(
m->face2node,
m->NumFace(), f2nSnap);
1288 CHECK(
m->face2cell(iF, 0) == f2cSnap[iF].first);
1289 CHECK(
m->face2cell(iF, 1) == f2cSnap[iF].second);
1295TEST_CASE(
"AdjC2F round-trip: cell2face and bnd2face identity")
1300 auto c2fSnap = snapshotAdj(
m->cell2face,
m->NumCell());
1301 auto b2fSnap = snapshotAdj1(
m->bnd2face,
m->NumBnd());
1303 m->AdjLocal2GlobalC2F();
1306 m->AdjGlobal2LocalC2F();
1309 checkAdjMatchesSnapshot(
m->cell2face,
m->NumCell(), c2fSnap);
1311 CHECK(
m->bnd2face(ib, 0) == b2fSnap[ib]);)
1314TEST_CASE(
"AdjN2CB round-trip: node2cell and node2bnd identity")
1319 auto n2cSnap = snapshotAdj(
m->node2cell,
m->NumNodeProc());
1320 auto n2bSnap = snapshotAdj(
m->node2bnd,
m->NumNodeProc());
1322 m->AdjLocal2GlobalN2CB();
1325 m->AdjGlobal2LocalN2CB();
1328 checkAdjMatchesSnapshot(
m->node2cell,
m->NumNodeProc(), n2cSnap);
1329 checkAdjMatchesSnapshot(
m->node2bnd,
m->NumNodeProc(), n2bSnap);)
1332TEST_CASE(
"AdjC2CFace round-trip: cell2cellFace identity")
1337 auto c2cfSnap = snapshotAdj(
m->cell2cellFace,
m->NumCell());
1339 m->AdjLocal2GlobalC2CFace();
1342 m->AdjGlobal2LocalC2CFace();
1345 checkAdjMatchesSnapshot(
m->cell2cellFace,
m->NumCell(), c2cfSnap);)
1348TEST_CASE(
"AdjPrimary round-trip: serial-out pattern")
1353 auto c2nSnap = snapshotAdj(
m->cell2node,
m->NumCell());
1354 auto c2cSnap = snapshotAdj(
m->cell2cell,
m->NumCell());
1355 auto b2nSnap = snapshotAdj(
m->bnd2node,
m->NumBnd());
1356 auto b2cSnap = snapshotAdj2(
m->bnd2cell,
m->NumBnd());
1358 m->AdjLocal2GlobalPrimary();
1361 DNDS::index globalNodeSize =
m->coords.father->globalSize();
1365 CHECK(gN < globalNodeSize);
1368 m->AdjGlobal2LocalPrimary();
1371 checkAdjMatchesSnapshot(
m->cell2node,
m->NumCell(), c2nSnap);
1372 checkAdjMatchesSnapshot(
m->cell2cell,
m->NumCell(), c2cSnap);
1373 checkAdjMatchesSnapshot(
m->bnd2node,
m->NumBnd(), b2nSnap);
1375 CHECK(
m->bnd2cell(ib, 0) == b2cSnap[ib].first);
1376 CHECK(
m->bnd2cell(ib, 1) == b2cSnap[ib].second);
1380TEST_CASE(
"AdjForBnd round-trip: ConstructBndMesh + ForBnd")
1384 m->ConstructBndMesh(bMesh);
1387 auto c2nSnap = snapshotAdj(bMesh.cell2node, bMesh.NumCell());
1389 bMesh.AdjLocal2GlobalPrimaryForBnd();
1392 DNDS::index globalNodeSize = bMesh.coords.father->globalSize();
1396 CHECK(gN < globalNodeSize);
1399 bMesh.AdjGlobal2LocalPrimaryForBnd();
1402 checkAdjMatchesSnapshot(bMesh.cell2node, bMesh.NumCell(), c2nSnap);)
1409TEST_CASE(
"Ball2: mixed 3D element types")
1412 if (g_mpi.
size != 8)
1419 std::map<Elem::ElemType, int> typeCounts;
1422 auto elem =
m->GetCellElement(iC);
1423 typeCounts[elem.type]++;
1427 std::vector<int> localTypes;
1428 for (
const auto &[type, count] : typeCounts)
1429 localTypes.
push_back(static_cast<int>(type));
1431 int nLocalTypes = localTypes.size();
1432 std::vector<int> allRankSizes(g_mpi.
size);
1433 MPI_Gather(&nLocalTypes, 1, MPI_INT, allRankSizes.data(), 1, MPI_INT, 0, g_mpi.
comm);
1435 std::vector<int> allTypes;
1436 std::vector<int> displs;
1439 if (g_mpi.
rank == 0)
1441 displs.resize(g_mpi.
size);
1442 for (
int i = 0;
i < g_mpi.
size;
i++)
1444 displs[
i] = totalTypes;
1445 totalTypes += allRankSizes[
i];
1447 allTypes.resize(totalTypes);
1450 MPI_Gatherv(localTypes.data(), nLocalTypes, MPI_INT,
1451 allTypes.data(), allRankSizes.data(), displs.data(), MPI_INT, 0, g_mpi.
comm);
1454 std::set<Elem::ElemType> uniqueTypes;
1455 if (g_mpi.
rank == 0)
1457 for (
int typeInt : allTypes)
1458 uniqueTypes.insert(static_cast<Elem::
ElemType>(typeInt));
1462 int nUnique = uniqueTypes.size();
1463 MPI_Bcast(&nUnique, 1, MPI_INT, 0, g_mpi.
comm);
1464 std::vector<int> uniqueTypesVec;
1465 if (g_mpi.
rank == 0)
1467 for (
auto type : uniqueTypes)
1468 uniqueTypesVec.
push_back(static_cast<int>(type));
1470 uniqueTypesVec.resize(nUnique);
1471 MPI_Bcast(uniqueTypesVec.data(), nUnique, MPI_INT, 0, g_mpi.
comm);
1474 std::map<Elem::ElemType, int> globalCounts;
1475 for (
int typeInt : uniqueTypesVec)
1478 int localCount = typeCounts.count(type) ? typeCounts[type] : 0;
1480 MPI_Reduce(&localCount, &globalCount, 1, MPI_INT, MPI_SUM, 0, g_mpi.
comm);
1481 if (g_mpi.
rank == 0)
1482 globalCounts[type] = globalCount;
1485 if (g_mpi.
rank == 0)
1487 std::cout <<
"\nBall2 element type counts:" << std::endl;
1488 for (
const auto &[type, count] : globalCounts)
1490 std::cout <<
" Type " <<
static_cast<int>(type) <<
": " << count << std::endl;
1495 for (
const auto &[type, count] : globalCounts)
1497 CHECK(total == 958994);
1505TEST_CASE(
"Elevation: O2 mesh cell count equals O1 cell count")
1508 const auto &cfg = g_configs[0];
1511 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1520 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1521 reader.Deduplicate1to1Periodic(1e-8);
1522 reader.BuildCell2Cell();
1529 reader.MeshPartitionCell2Cell(pOpt);
1530 reader.PartitionReorderToMeshCell2Cell();
1532 mO1->RecoverNode2CellAndNode2Bnd();
1533 mO1->RecoverCell2CellAndBnd2Cell();
1534 mO1->BuildGhostPrimary();
1535 mO1->AdjGlobal2LocalPrimary();
1538 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1539 mO2->BuildO2FromO1Elevation(*mO1);
1542 mO2->RecoverNode2CellAndNode2Bnd();
1543 mO2->RecoverCell2CellAndBnd2Cell();
1544 mO2->BuildGhostPrimary();
1545 mO2->AdjGlobal2LocalPrimary();
1551 CHECK(nCellO2 == nCellO1);
1554 for (
DNDS::index iC = 0; iC < mO2->NumCell(); iC++)
1556 auto elem = mO2->GetCellElement(iC);
1557 CHECK(elem.GetOrder() == 2);
1561TEST_CASE(
"Elevation: O2 mesh has more nodes than O1")
1563 const auto &cfg = g_configs[0];
1566 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1575 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1576 reader.Deduplicate1to1Periodic(1e-8);
1577 reader.BuildCell2Cell();
1584 reader.MeshPartitionCell2Cell(pOpt);
1585 reader.PartitionReorderToMeshCell2Cell();
1587 mO1->RecoverNode2CellAndNode2Bnd();
1588 mO1->RecoverCell2CellAndBnd2Cell();
1589 mO1->BuildGhostPrimary();
1590 mO1->AdjGlobal2LocalPrimary();
1593 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1594 mO2->BuildO2FromO1Elevation(*mO1);
1597 mO2->RecoverNode2CellAndNode2Bnd();
1598 mO2->RecoverCell2CellAndBnd2Cell();
1599 mO2->BuildGhostPrimary();
1600 mO2->AdjGlobal2LocalPrimary();
1602 DNDS::index nNodeO1 = mO1->coords.father->globalSize();
1603 DNDS::index nNodeO2 = mO2->coords.father->globalSize();
1606 CHECK(nNodeO2 > nNodeO1);
1609TEST_CASE(
"Bisection: O1 mesh has more cells than O2")
1611 const auto &cfg = g_configs[0];
1614 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1623 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1624 reader.Deduplicate1to1Periodic(1e-8);
1625 reader.BuildCell2Cell();
1632 reader.MeshPartitionCell2Cell(pOpt);
1633 reader.PartitionReorderToMeshCell2Cell();
1635 mO1->RecoverNode2CellAndNode2Bnd();
1636 mO1->RecoverCell2CellAndBnd2Cell();
1637 mO1->BuildGhostPrimary();
1638 mO1->AdjGlobal2LocalPrimary();
1641 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1642 mO2->BuildO2FromO1Elevation(*mO1);
1643 mO2->RecoverNode2CellAndNode2Bnd();
1644 mO2->RecoverCell2CellAndBnd2Cell();
1645 mO2->BuildGhostPrimary();
1646 mO2->AdjGlobal2LocalPrimary();
1649 mO2->AdjLocal2GlobalPrimary();
1650 auto mO1B = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1651 mO1B->BuildBisectO1FormO2(*mO2);
1657 CHECK(nCellO1B > nCellO2);
1660 for (
DNDS::index iC = 0; iC < mO1B->NumCell(); iC++)
1662 auto elem = mO1B->GetCellElement(iC);
1663 CHECK(elem.GetOrder() == 1);
1667TEST_CASE(
"Bisection: global node count is preserved from O2 mesh")
1669 const auto &cfg = g_configs[0];
1672 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1681 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1682 reader.Deduplicate1to1Periodic(1e-8);
1683 reader.BuildCell2Cell();
1690 reader.MeshPartitionCell2Cell(pOpt);
1691 reader.PartitionReorderToMeshCell2Cell();
1693 mO1->RecoverNode2CellAndNode2Bnd();
1694 mO1->RecoverCell2CellAndBnd2Cell();
1695 mO1->BuildGhostPrimary();
1696 mO1->AdjGlobal2LocalPrimary();
1699 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1700 mO2->BuildO2FromO1Elevation(*mO1);
1701 mO2->RecoverNode2CellAndNode2Bnd();
1702 mO2->RecoverCell2CellAndBnd2Cell();
1703 mO2->BuildGhostPrimary();
1704 mO2->AdjGlobal2LocalPrimary();
1707 mO2->AdjLocal2GlobalPrimary();
1708 auto mO1B = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1709 mO1B->BuildBisectO1FormO2(*mO2);
1712 DNDS::index nNodeO2 = mO2->coords.father->globalSize();
1713 DNDS::index nNodeO1B = mO1B->coords.father->globalSize();
1714 CHECK(nNodeO1B == nNodeO2);
1717TEST_CASE(
"Elevation+Bisection: exact cell count progression")
1722 const auto &cfg = g_configs[0];
1725 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1734 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1735 reader.Deduplicate1to1Periodic(1e-8);
1736 reader.BuildCell2Cell();
1743 reader.MeshPartitionCell2Cell(pOpt);
1744 reader.PartitionReorderToMeshCell2Cell();
1746 mO1->RecoverNode2CellAndNode2Bnd();
1747 mO1->RecoverCell2CellAndBnd2Cell();
1748 mO1->BuildGhostPrimary();
1749 mO1->AdjGlobal2LocalPrimary();
1752 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1753 mO2->BuildO2FromO1Elevation(*mO1);
1754 mO2->RecoverNode2CellAndNode2Bnd();
1755 mO2->RecoverCell2CellAndBnd2Cell();
1756 mO2->BuildGhostPrimary();
1757 mO2->AdjGlobal2LocalPrimary();
1760 mO2->AdjLocal2GlobalPrimary();
1761 auto mO1B = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1762 mO1B->BuildBisectO1FormO2(*mO2);
1766 DNDS::index nCellBisected = mO1B->NumCellGlobal();
1768 CHECK(nCellOrig == 100);
1769 CHECK(nCellElevated == 100);
1770 CHECK(nCellBisected == 400);
1787 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1793 mO1->SetPeriodicGeometry(
1799 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1800 reader.Deduplicate1to1Periodic(1e-8);
1801 reader.BuildCell2Cell();
1808 reader.MeshPartitionCell2Cell(pOpt);
1809 reader.PartitionReorderToMeshCell2Cell();
1811 mO1->RecoverNode2CellAndNode2Bnd();
1812 mO1->RecoverCell2CellAndBnd2Cell();
1813 mO1->BuildGhostPrimary();
1814 mO1->AdjGlobal2LocalPrimary();
1817 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1818 mO2->BuildO2FromO1Elevation(*mO1);
1821 mO2->RecoverNode2CellAndNode2Bnd();
1822 mO2->RecoverCell2CellAndBnd2Cell();
1823 mO2->BuildGhostPrimary();
1824 mO2->AdjGlobal2LocalPrimary();
1825 mO2->AdjGlobal2LocalN2CB();
1827 mO2->InterpolateFace();
1828 mO2->AssertOnFaces();
1830 mO2->AdjLocal2GlobalN2CB();
1831 mO2->BuildGhostN2CB();
1832 mO2->AdjGlobal2LocalN2CB();
1834 mO2->BuildCell2CellFace();
1835 mO2->AdjGlobal2LocalC2CFace();
1841TEST_CASE(
"BoundarySmooth: NACA0012 moves boundary O2 nodes")
1844 const auto &cfg = g_configs[2];
1847 auto mO2 = buildO2MeshWithFaces(cfg, mO1);
1850 std::vector<tPoint> preSmooth(mO2->NumNode());
1852 preSmooth[
i] = mO2->coords[
i];
1855 mO2->ElevatedNodesGetBoundarySmooth(
1857 {
return Geom::FaceIDIsExternalBC(bndId); });
1860 CHECK(mO2->nTotalMoved > 0);
1863 CHECK(mO2->coordsElevDisp.father !=
nullptr);
1868 for (
int d = 0; d < 3; d++)
1869 CHECK(std::isfinite(mO2->coords[
i](d)));
1873TEST_CASE(
"BoundarySmooth: no NaN in displacement field")
1875 const auto &cfg = g_configs[2];
1878 auto mO2 = buildO2MeshWithFaces(cfg, mO1);
1880 mO2->ElevatedNodesGetBoundarySmooth(
1882 {
return Geom::FaceIDIsExternalBC(bndId); });
1886 for (
DNDS::index i = 0;
i < mO2->coordsElevDisp.father->Size();
i++)
1888 auto disp = mO2->coordsElevDisp[
i];
1890 if (disp(0) !=
UnInitReal && std::abs(disp(0)) < 1e30)
1892 for (
int d = 0; d < 3; d++)
1893 CHECK(std::isfinite(disp(d)));
1900 CHECK(nMovedGlobal > 0);
1903TEST_CASE(
"InternalSmooth V2: coordinates remain finite after solve")
1905 const auto &cfg = g_configs[2];
1908 auto mO2 = buildO2MeshWithFaces(cfg, mO1);
1911 mO2->elevationInfo.nIter = 10;
1912 mO2->elevationInfo.nSearch = 20;
1913 mO2->elevationInfo.RBFRadius = 5.0;
1916 mO2->ElevatedNodesGetBoundarySmooth(
1918 {
return Geom::FaceIDIsExternalBC(bndId); });
1920 if (mO2->nTotalMoved == 0)
1924 std::vector<tPoint> preSmoothCoords(mO2->NumNode());
1926 preSmoothCoords[
i] = mO2->coords[
i];
1929 mO2->ElevatedNodesSolveInternalSmoothV2();
1933 for (
int d = 0; d < 3; d++)
1934 CHECK(std::isfinite(mO2->coords[
i](d)));
1939 if ((mO2->coords[
i] - preSmoothCoords[
i]).norm() > 1e-15)
1943 MPI_Allreduce(&nChangedLocal, &nChangedGlobal, 1,
DNDS_MPI_INDEX, MPI_SUM, g_mpi.
comm);
1944 CHECK(nChangedGlobal > 0);
1947TEST_CASE(
"Elevation: O2 coordinates have no NaN")
1950 for (
int ci : {0, 2})
1952 const auto &cfg = g_configs[ci];
1955 auto mO1 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1958 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
1959 reader.Deduplicate1to1Periodic(1e-8);
1960 reader.BuildCell2Cell();
1967 reader.MeshPartitionCell2Cell(pOpt);
1968 reader.PartitionReorderToMeshCell2Cell();
1970 mO1->RecoverNode2CellAndNode2Bnd();
1971 mO1->RecoverCell2CellAndBnd2Cell();
1972 mO1->BuildGhostPrimary();
1973 mO1->AdjGlobal2LocalPrimary();
1975 auto mO2 = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
1976 mO2->BuildO2FromO1Elevation(*mO1);
1980 for (
int d = 0; d < 3; d++)
1981 CHECK(std::isfinite(mO2->coords[
i](d)));
1995 auto mesh = std::make_shared<UnstructuredMesh>(g_mpi, cfg.
dim);
2001 mesh->SetPeriodicGeometry(
2007 reader.ReadFromCGNSSerial(meshPath(cfg.
file));
2008 reader.Deduplicate1to1Periodic(1e-8);
2009 reader.BuildCell2Cell();
2016 reader.MeshPartitionCell2Cell(pOpt);
2017 reader.PartitionReorderToMeshCell2Cell();
2019 mesh->RecoverNode2CellAndNode2Bnd();
2020 mesh->RecoverCell2CellAndBnd2Cell();
2021 mesh->BuildGhostPrimary(nGhostLayers);
2022 mesh->AdjGlobal2LocalPrimary();
2023 mesh->AdjGlobal2LocalN2CB();
2025 mesh->InterpolateFace();
2026 mesh->AssertOnFaces();
2028 mesh->AdjLocal2GlobalN2CB();
2029 mesh->BuildGhostN2CB();
2030 mesh->AdjGlobal2LocalN2CB();
2035TEST_CASE(
"MultiLayerGhost: ghost counts increase with layer count")
2037 const auto &cfg = g_configs[0];
2039 DNDS::index cellGhost[3], nodeGhost[3], bndGhost[3];
2041 for (
int nL = 1; nL <= 3; nL++)
2043 auto m = buildMeshWithGhostLayers(cfg, nL);
2044 cellGhost[nL - 1] =
m->NumCellGhost();
2045 nodeGhost[nL - 1] =
m->NumNodeGhost();
2046 bndGhost[nL - 1] =
m->NumBndGhost();
2048 if (g_mpi.
rank == 0)
2049 fmt::print(
" nLayers={}: cellGhost={}, nodeGhost={}, bndGhost={}\n",
2050 nL, cellGhost[nL - 1], nodeGhost[nL - 1], bndGhost[nL - 1]);
2054 for (
int i = 0;
i < 2;
i++)
2056 CHECK(cellGhost[
i + 1] >= cellGhost[
i]);
2057 CHECK(nodeGhost[
i + 1] >= nodeGhost[
i]);
2058 CHECK(bndGhost[
i + 1] >= bndGhost[
i]);
2062 if (g_mpi.
size >= 2)
2064 CHECK(cellGhost[1] > cellGhost[0]);
2065 CHECK(cellGhost[2] > cellGhost[1]);
2066 CHECK(nodeGhost[1] > nodeGhost[0]);
2067 CHECK(nodeGhost[2] > nodeGhost[1]);
2071TEST_CASE(
"MultiLayerGhost: 2-layer cell2cell neighbors are all local")
2076 const auto &cfg = g_configs[0];
2077 auto m = buildMeshWithGhostLayers(cfg, 2);
2084 for (
DNDS::index iCell = 0; iCell < nCellOwned; iCell++)
2091 REQUIRE(iNeighbor < nCellProc);
2093 for (
DNDS::rowsize k = 0; k <
m->cell2cell.RowSize(iNeighbor); k++)
2098 CHECK(iNeighbor2 < nCellProc);
2104TEST_CASE(
"MultiLayerGhost: global cell count preserved across layers")
2106 const auto &cfg = g_configs[0];
2108 for (
int nL = 1; nL <= 3; nL++)
2111 auto m = buildMeshWithGhostLayers(cfg, nL);
#define DNDS_assert(expr)
Debug-only assertion (compiled out when DNDS_NDEBUG is defined). Prints the expression + file/line + ...
Eigen::Matrix< real, 3, 3 > m
the host side operators are provided as implemented
const MPI_Datatype DNDS_MPI_INDEX
MPI datatype matching index (= MPI_INT64_T).
DNDS_CONSTANT const index UnInitIndex
Sentinel "not initialised" index value (= INT64_MIN).
int32_t rowsize
Row-width / per-row element-count type (signed 32-bit).
DNDS_CONSTANT const real UnInitReal
Sentinel "not initialised" real value (NaN). Cheap to detect with std::isnan or IsUnInitReal; survive...
int64_t index
Global row / DOF index type (signed 64-bit; handles multi-billion-cell meshes).
std::shared_ptr< T > ssp
Shortened alias for std::shared_ptr used pervasively in DNDSR.
Convenience bundle of a father, son, and attached ArrayTransformer.
DNDS_DEVICE_CALLABLE constexpr Element ObtainFace(t_index iFace) const
DNDS_DEVICE_CALLABLE constexpr t_index GetNumVertices() const
Lightweight bundle of an MPI communicator and the calling rank's coordinates.
int size
Number of ranks in comm (-1 until initialised).
int rank
This rank's 0-based index within comm (-1 until initialised).
MPI_Comm comm
The underlying MPI communicator handle.
void setWorld()
Initialise the object to MPI_COMM_WORLD. Requires MPI_Init to have run.
DNDS::index expectedCells
Expected global cell count (-1 = skip check)
tPoint translation1
Periodic translation axis 1.
bool periodic
Has periodic boundaries.
int dim
Spatial dimension.
tPoint translation3
Periodic translation axis 3.
const char * file
Filename relative to data/mesh/.
tPoint translation2
Periodic translation axis 2.
const char * name
Human-readable label.
DNDS::index expectedBnds
Expected global bnd count (-1 = skip check)
REQUIRE(bool(result.parent2entityPbi.father))
#define FOR_EACH_MESH_CONFIG(body)
std::pair< DNDS::index, uint8_t > idx_pbi_t
input explicitMaps push_back(EntityReorderMap{EntityKind::Cell, cellPartition})
TEST_CASE("3D: VFV P2 HQM error < P1 on sinCos3D")