feat: new crackAlong method

This new method is a more powerfull version of the
buildInnerBoundariesAlongM1Group method because it naturally handle non
connex and complex M1 groups

wip: fix cycling over nodal connectivity

fix shared_ptr<set<int>> allocation

wip: no more memleak !

TODO: for now this method only treat the mesh at level 0. The mesh at
level -1 and the groups lying on it should be taken care of.

wip: Duplicate faces inplace and replace old node number with new ones

wip: removing useless shared_ptr, adding doc

wip: moving duplicateFaces to CrackAlgo.cxx

wip: refactor CrackAlgo and treat elements families

missing node families

feat: functional new duplicateFaces method

fix: changing API and preparing methods for test

feat: separate tests on crackAlong method

wip: using M1 and not Mf

fix: case when cell has no neighbor

feat: test mesh at level -1

This test checks that the crack mesh has twice as many cells than the
starting crack mesh, and that they are both equivalent after merging the
nodes.

wip: fixing connectivity of faces which are not duplicated

fix: manage connectivity of m1 elements touching crack

wip: add python tests

Still figuring out why using M1 non connex group erease the families at lev -2

fix: crackAlong exceptions with boolean param

A boolean is added to change the method behavior in the case when a face
is duplicated but not modified (there is no real separation between
volumes)

feat: new tests and OpenCrack binding

fix: Copy family at lev0 and set them back

Because suprisingly they tend to disapear.

refactor: OpenCrack -> openCrack

fix: remove french comments, better use of stl

refactor: better documentation, less auto, more const

diff --git a/src/MEDLoader/CMakeLists.txt b/src/MEDLoader/CMakeLists.txt
index aec75954d6e6560049c836a18447e0da81e5145f..6062d306203e1e2f9f33b8c2205b8dc6be9b308d 100644 (file)
--- a/src/MEDLoader/CMakeLists.txt
+++ b/src/MEDLoader/CMakeLists.txt
@@ -59,6 +59,7 @@ INCLUDE_DIRECTORIES(
   )
 
 SET(medloader_SOURCES
+  CrackAlgo.cxx
   MEDLoader.cxx
   MEDLoaderBase.cxx
   MEDLoaderTraits.cxx

diff --git a/src/MEDLoader/CrackAlgo.cxx b/src/MEDLoader/CrackAlgo.cxx
new file mode 100644 (file)
index 0000000..836c14f
--- /dev/null
+++ b/src/MEDLoader/CrackAlgo.cxx
@@ -0,0 +1,706 @@
+// Copyright (C) 2007-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email :
+// webmaster.salome@opencascade.com
+//
+// Author : Aymeric SONOLET (CEA/DES)
+
+#include "CrackAlgo.hxx"
+
+#include <algorithm>
+#include <unordered_map>
+#include <unordered_set>
+#include <string>
+#include <memory>
+#include <utility>
+
+#include "./MEDFileMesh.hxx"
+
+#include <MEDCouplingUMesh.hxx>
+#include <MEDCouplingMemArray.hxx>
+#include <InterpKernelException.hxx>
+#include <MCType.hxx>
+#include <MCAuto.hxx>
+
+using namespace MEDCoupling;
+using namespace std;
+using MCU = MCAuto<MEDCouplingUMesh>;
+using DAI = MCAuto<DataArrayIdType>;
+using DAD = MCAuto<DataArrayDouble>;
+
+
+CrackAlgo::Map2Map
+CrackAlgo::Compute(
+    MEDFileUMesh* mm,
+    const string& grp_name,
+    bool grpMustBeFullyDup
+) {
+
+    vector<int> levs = mm->getNonEmptyLevels();
+    if (find(levs.begin(), levs.end(), 0) == levs.end() ||
+            find(levs.begin(), levs.end(), -1) == levs.end())
+        throw INTERP_KERNEL::Exception(
+                "MEDFileUMesh::duplicateFaces : This method "
+                "works only for mesh defined on level 0 and -1 !");
+
+    MCU m0 = mm->getMeshAtLevel(0);
+
+    MCU crackingMesh = mm->getGroup(-1, grp_name);
+    MCU cleanCrackingMesh = CleanM1Mesh(*m0, *crackingMesh);
+
+    DAI c2f(DataArrayIdType::New()),
+        c2fIdx(DataArrayIdType::New()),
+        f2c(DataArrayIdType::New()),
+        f2cIdx(DataArrayIdType::New());
+    const MCU mf = m0->buildDescendingConnectivity(c2f, c2fIdx, f2c, f2cIdx);
+
+    DataArrayIdType * crackToMfP;
+    const bool crackCellsInMf = mf->areCellsIncludedIn(cleanCrackingMesh, 2, crackToMfP);
+    DAI f2dupIdInMf(crackToMfP);
+    if (!crackCellsInMf)
+        throw INTERP_KERNEL::Exception(
+            "crackAlong: All cells in crack are not part of Mf.");
+
+    DAI n2c(DataArrayIdType::New()),
+          n2cIdx(DataArrayIdType::New());
+    m0->getReverseNodalConnectivity(n2c, n2cIdx);
+
+    const Map2Set n2c_dup = GetNode2CellMap(mf, n2cIdx, n2c, f2dupIdInMf);
+
+    const Set cTouchingN_dup = GetCellsTouchingNodesToDup(
+        mf, n2cIdx, n2c, f2dupIdInMf);
+
+    const Graph c2c = BuildCutC2CGraph(
+        c2fIdx, c2f, f2cIdx, f2c, cTouchingN_dup, f2dupIdInMf);
+
+    Map2Map cellOld2NewNode = CreateNewNodesInTopLevelMesh(n2c_dup, c2c, m0);
+    // End of node creation, separation of M0
+
+    // Faces in M1 must be added/updated
+
+    MCU m1 = mm->getMeshAtLevel(-1);
+
+    // Before changing M1
+    const DAI f2dupIdInM1 = GetFacesToDupInM1(*cleanCrackingMesh, *m1);
+
+    const auto f2changeM1Mf = GetFacesInM1TouchingDuplicatedNodes(
+        n2c_dup, f2dupIdInM1, *mf, *m1);
+    const DAI f2changeIdInM1 = f2changeM1Mf.first;
+    const DAI f2changeIdInMf = f2changeM1Mf.second;
+
+    AddMissingElementsOnLevelM1AndChangeConnectivity(
+        f2cIdx, f2c, f2dupIdInM1, f2dupIdInMf, cellOld2NewNode, m1, grpMustBeFullyDup);
+
+    ChangeConnectivityOfM1Elements(
+            f2changeIdInM1, f2changeIdInMf, cellOld2NewNode, f2cIdx, f2c, m1);
+
+    // TODO(aymeric): If one wants to implement
+    // AddMissingElementsOnLevelM2AndChangeConnectivity and
+    // ChangeConnectivityOfM2Elements it should be there.
+
+    DAI famLev0 = CopyFamilyArrAtLev0(mm);
+    DAI newFam = CopyAndCompleteFamilyArrAtLevelM1(mm, m1, f2dupIdInM1);
+
+    mm->setMeshAtLevel(0, m0);
+    mm->setMeshAtLevel(-1, m1);
+    // mm->setMeshAtLevel(-2, m2);
+
+    mm->setFamilyFieldArr(0, famLev0);
+    mm->setFamilyFieldArr(-1, newFam);
+    // mm->setFamilyFieldArr(-2, newFam2);
+
+    const Map2Set addedNodes = BuildMap2Set(cellOld2NewNode);
+    CompleteFamilyArrAtNodeLevel(addedNodes, mm);
+
+    return cellOld2NewNode;
+}
+
+CrackAlgo::Map2Set
+CrackAlgo::BuildMap2Set(
+    const Map2Map& cellOld2NewNode
+) {
+    Map2Set res;
+    for (const auto & keyElem : cellOld2NewNode) {
+        const auto & old2NewMap = keyElem.second;
+        for (const auto & old2New : old2NewMap) {
+            res[old2New.first].insert(old2New.second);
+        }
+    }
+    return res;
+}
+
+DataArrayIdType *
+CrackAlgo::CopyAndCompleteFamilyArrAtLevelM1(
+    const MEDFileUMesh * mm,
+    const MEDCouplingUMesh * mf,
+    const DataArrayIdType * f2dup
+) {
+    DataArrayIdType* newFam = nullptr;
+    const DataArrayIdType *fam = mm->getFamilyFieldAtLevel(-1);
+    if (fam != nullptr) {
+        newFam = DataArrayIdType::New();
+        const auto fam_size = static_cast<size_t>(fam->getNumberOfTuples());
+        newFam->alloc(static_cast<size_t>(mf->getNumberOfCells()));
+        copy(fam->begin(), fam->end(), newFam->getPointer());
+        mcIdType * fam_ptr = newFam->getPointer();
+        size_t i_elem_added = 0;
+        for (const auto & face : *f2dup) {
+            // NOTE(aymeric): assumes all faces are duplicated
+            fam_ptr[fam_size + i_elem_added] = fam_ptr[face];
+            i_elem_added++;
+        }
+    }
+    return newFam;
+}
+
+DataArrayIdType *
+CrackAlgo::CopyFamilyArrAtLev0(
+    const MEDFileUMesh * mm
+) {
+    DataArrayIdType* famLev0 = nullptr;
+    const DataArrayIdType *fam = mm->getFamilyFieldAtLevel(0);
+    if (fam != nullptr) {
+        famLev0 = DataArrayIdType::New();
+        const auto fam_size = static_cast<size_t>(fam->getNumberOfTuples());
+        famLev0->alloc(fam_size);
+        copy(fam->begin(), fam->end(), famLev0->getPointer());
+    }
+    return famLev0;
+}
+
+void
+CrackAlgo::CompleteFamilyArrAtNodeLevel(
+    const Map2Set& addedNodes,
+    MEDFileUMesh * mm
+) {
+    DataArrayIdType *node_fam = mm->getFamilyFieldAtLevel(1);
+    if (node_fam != nullptr) {
+        node_fam->reAlloc(static_cast<size_t>(mm->getCoords()->getNumberOfTuples()));
+        mcIdType * node_fam_ptr = node_fam->getPointer();
+        for (const auto & old2NewPair : addedNodes) {
+            const auto & oldNode = old2NewPair.first;
+            const auto & newNodes = old2NewPair.second;
+            for (const auto & newNode : newNodes) {
+                node_fam_ptr[newNode] = node_fam_ptr[oldNode];
+            }
+        }
+    }
+}
+
+MEDCouplingUMesh *
+CrackAlgo::CleanM1Mesh(
+    const MEDCouplingUMesh & m,
+    const MEDCouplingUMesh & crackingMesh
+) {
+    MCU m0skin = m.computeSkin();
+    DataArrayIdType * idsToKeepP;
+    m0skin->areCellsIncludedIn(&crackingMesh, 2, idsToKeepP);
+    DAI idsToKeep(idsToKeepP);
+    // discard cells on the skin of M0
+    DAI ids2 = idsToKeep->findIdsNotInRange(0, m0skin->getNumberOfCells());
+    MCU otherDimM1OnSameCoords =
+        crackingMesh.buildPartOfMySelf(ids2->begin(), ids2->end(), true);
+    return otherDimM1OnSameCoords.retn();
+}
+
+void CrackAlgo::Dfs(
+    const Graph& graph,
+    const mcIdType &node,
+    const size_t &componentId,
+    map<mcIdType, bool>& visited,
+    map<mcIdType, size_t>& componentMap
+) {
+    visited[node] = true;
+    // Assign a unique component ID to the node
+    componentMap[node] = componentId;
+
+    for (const auto& neighbor : graph.at(node)) {
+        if (!visited[neighbor]) {
+            Dfs(graph, neighbor, componentId, visited, componentMap);
+        }
+    }
+}
+
+// Function to find connected components
+vector<shared_ptr<vector<mcIdType>>>
+CrackAlgo::FindConnectedComponents(
+    const Graph& graph
+) {
+    map<mcIdType, bool> visited;
+    map<mcIdType, size_t> componentMap;
+    size_t componentId = 0;
+
+    for (const auto &nodePair : graph)
+        visited[nodePair.first] = false;
+
+    for (const auto& nodePair : graph) {
+        if (!visited[nodePair.first]) {
+            Dfs(graph, nodePair.first, componentId, visited, componentMap);
+            componentId++;
+        }
+    }
+
+    using spv = shared_ptr<vector<mcIdType>>;
+    vector<spv> components;
+    for (size_t i = 0; i < componentId; i++)
+        components.push_back(make_shared<vector<mcIdType>>());
+    for (const auto& nodePair : componentMap) {
+        components[nodePair.second]->push_back(nodePair.first);
+    }
+
+    return components;
+}
+
+CrackAlgo::Set CrackAlgo::DataArrayToSet(
+    const DataArrayIdType& da
+) {
+    Set res;
+    for (const auto &elem : da) {
+        res.insert(elem);
+    }
+    return res;
+}
+
+DataArrayIdType * CrackAlgo::SetToDataArray(
+    const Set& s
+) {
+    DataArrayIdType * da(DataArrayIdType::New());
+    da->alloc(s.size());
+    mcIdType * da_p = da->rwBegin();
+    size_t i = 0;
+    for (const auto & elem : s) {
+        da_p[i] = elem;
+        i++;
+    }
+    return da;
+}
+
+CrackAlgo::Graph
+CrackAlgo::BuildCutC2CGraph(
+    const DataArrayIdType * c2fIdx,
+    const DataArrayIdType * c2f,
+    const DataArrayIdType * f2cIdx,
+    const DataArrayIdType * f2c,
+    const Set & cTouchingN_dup,
+    const DataArrayIdType * f2dup
+) {
+
+    Set f2dup_set = DataArrayToSet(*f2dup);
+    Graph c2c;
+
+    const mcIdType * c2fIdx_p {c2fIdx->begin()};
+    const mcIdType * c2f_p {c2f->begin()};
+    const mcIdType * f2cIdx_p {f2cIdx->begin()};
+    const mcIdType * f2c_p {f2c->begin()};
+
+    for (const auto &cell : cTouchingN_dup) {
+        c2c[cell] = std::unordered_set<mcIdType>();
+        for (auto faceIdx = c2fIdx_p[cell]; faceIdx < c2fIdx_p[cell + 1]; faceIdx++) {
+            const mcIdType &face = c2f_p[faceIdx];
+            // face is not in face to duplicate
+            if (f2dup_set.find(face) == f2dup_set.end()) {
+                for (mcIdType cellIdx {f2cIdx_p[face]}; cellIdx < f2cIdx_p[face +1]; cellIdx++) {
+                    // getting the id of the other cell (on face is
+                    // connected to two cells)
+                    const auto& cell2 { f2c_p[cellIdx] };
+                    if (
+                        (cell != cell2)
+                        & (cTouchingN_dup.find(cell2) != cTouchingN_dup.end())
+                    ) {
+                        c2c.at(cell).insert(cell2);
+                    }
+                }
+            }
+        }
+    }
+    return c2c;
+}
+
+CrackAlgo::Map2Map
+CrackAlgo::CreateNewNodesInTopLevelMesh(
+    const Map2Set & n2c_dup,
+    const Graph & c2c,
+    MEDCouplingUMesh* m0
+) {
+    DataArrayIdType* c2n = m0->getNodalConnectivity();
+    DataArrayIdType* c2nIdx = m0->getNodalConnectivityIndex();
+    mcIdType * c2n_ptr = c2n->rwBegin();
+    const mcIdType * c2nIdx_ptr = c2nIdx->begin();
+
+    DataArrayDouble* coords = m0->getCoords();
+    mcIdType i = coords->getNumberOfTuples();
+    const size_t coordsDim = coords->getNumberOfComponents();
+
+    Map2Map cellOld2NewNode;
+
+    for (const auto &nc_pair : n2c_dup) {
+        const auto &node = nc_pair.first;
+        const auto &cells = nc_pair.second;
+
+        // Building local connection graph
+        Graph c2c_local;
+        for (const auto &cell : cells) {
+            auto& val = c2c_local[cell];
+            if (c2c.find(cell) == c2c.end())
+                throw INTERP_KERNEL::Exception(
+                    "A cell touching a node is not part of c2c.");
+            for (const auto &neighbor : c2c.at(cell)) {
+                if (cells.find(neighbor) != cells.end())
+                    val.insert(neighbor);
+            }
+        }
+
+        const vector<shared_ptr<vector<mcIdType>>> compo_connex = FindConnectedComponents(c2c_local);
+
+        // Duplicate node for compo[1:]
+        const size_t i_nodes_to_add = compo_connex.size() - 1;
+        coords->reAlloc(static_cast<size_t>(i) + i_nodes_to_add);
+
+        int i_compo = 0;
+        for (const auto &compo : compo_connex) {
+            if (i_compo > 0) {
+                // Coords are copied at the end of the vector of coords, ie the
+                // node is duplicated
+                coords->getTuple(
+                        node,
+                        coords->getPointer() + static_cast<size_t>(i) * coordsDim);
+                for (const auto &cell : *compo) {
+                    // The node number is replaced in all corresponding cells
+                    std::replace(&c2n_ptr[c2nIdx_ptr[cell] + 1], &c2n_ptr[c2nIdx_ptr[cell + 1]], node, i);
+                    // This map is build in order to assign later new nodes to
+                    // the corresponding faces
+                    cellOld2NewNode[cell][node] = i;
+                }
+                i++;
+            }
+            i_compo++;
+        }
+    }
+    return cellOld2NewNode;
+}
+
+void
+CrackAlgo::AddMissingElementsOnLevelM1AndChangeConnectivity(
+    const DataArrayIdType * f2cIdx,
+    const DataArrayIdType * f2c,
+    const DataArrayIdType * f2dupIdInM1,
+    const DataArrayIdType * f2dupIdInMf,
+    const Map2Map & cellOld2NewNode,
+    MEDCouplingUMesh * m1,
+    bool grpMustBeFullyDup
+) {
+    DataArrayIdType* f2nIdx = m1->getNodalConnectivityIndex();
+    DataArrayIdType* f2n = m1->getNodalConnectivity();
+
+    // Change connectivity index size
+    const auto f2nIdx_size = static_cast<size_t>(f2nIdx->getNumberOfTuples());  // = nb_face + 1
+    const auto f2dup_size = static_cast<size_t>(f2dupIdInM1->getNumberOfTuples());
+    if (f2dup_size != static_cast<size_t>(f2dupIdInMf->getNumberOfTuples()))
+        throw INTERP_KERNEL::Exception(
+            "There is not as many faces to dup looking into Mf and into M1.");
+
+    f2nIdx->reAlloc(f2nIdx_size + f2dup_size);
+
+    // Change connectivity size
+    const auto f2n_size = static_cast<size_t>(f2n->getNumberOfTuples());
+    const mcIdType * const f2nIdx_p = f2nIdx->begin();
+    {
+        size_t new_connectivity_size = f2n_size;
+        for (const auto& face : *f2dupIdInM1)
+            new_connectivity_size += static_cast<size_t>(
+                f2nIdx_p[face + 1] - f2nIdx_p[face]);
+        f2n->reAlloc(new_connectivity_size);
+    }
+
+    mcIdType * const f2nIdx_pw = f2nIdx->rwBegin();
+    const mcIdType * const f2n_p = f2n->begin();
+    mcIdType * const f2n_pw = f2n->rwBegin();
+    size_t lastFace = f2nIdx_size - 1;
+
+    const mcIdType * const f2cIdx_p = f2cIdx->begin();
+    const mcIdType * const f2c_p = f2c->begin();
+
+    for (size_t iFace = 0; iFace < f2dup_size; iFace++) {
+        const auto & faceM1 = f2dupIdInM1->begin()[iFace];
+        const auto & faceMf = f2dupIdInMf->begin()[iFace];
+
+        // Copy the face appending it
+        //
+        //   1. Append last connectivity index number
+        const mcIdType face_start = f2nIdx_p[faceM1];
+        const mcIdType face_end = f2nIdx_p[faceM1 + 1];
+        const mcIdType lenFace = face_end - face_start;
+        f2nIdx_pw[lastFace + 1] = f2nIdx_p[lastFace] + lenFace;
+
+        //   2. Append last face connectivity
+        mcIdType *const new_f2n_p_start = f2n_pw + f2nIdx_p[lastFace];
+        copy(f2n_p + face_start, f2n_p + face_end, new_f2n_p_start);
+
+        // Change inplace the connectivity
+        // Check that this face is an inner face, ie it is connected to two
+        // cells
+        const mcIdType d = f2cIdx_p[faceMf + 1] - f2cIdx_p[faceMf];
+        if (d != 2)
+            throw INTERP_KERNEL::Exception(
+                "MEDFileMesh::duplicateFaces, the face to cell (or node to"
+                "segment) DataArray does not always adress two cells.");
+
+        bool noNodesAreChanged = true;
+
+        //   1. In original face, attaching it to cell0
+        const auto& cell0 = f2c_p[f2cIdx_p[faceMf]];
+        // If nodes where changed in cell0
+        if (cellOld2NewNode.find(cell0) != cellOld2NewNode.end()) {
+            const auto & mapO2N0 = cellOld2NewNode.at(cell0);
+            for (mcIdType i_node = f2nIdx_p[faceM1] + 1; i_node < f2nIdx_p[faceM1 + 1]; i_node++) {
+                const mcIdType node = f2n_p[i_node];
+                // If this node was modified
+                if (mapO2N0.find(node) != mapO2N0.end()) {
+                    const auto& new_node = mapO2N0.at(node);
+                    f2n_pw[i_node] = new_node;
+                    noNodesAreChanged = false;
+                }
+            }
+        }
+
+        //   2. In newly created face, attaching it to cell1
+        const auto& cell1 = f2c_p[f2cIdx_p[faceMf] + 1];
+        // If nodes where changed in cell1
+        if (cellOld2NewNode.find(cell1) != cellOld2NewNode.end()) {
+            const auto & mapO2N1 = cellOld2NewNode.at(cell1);
+            for (mcIdType i_node = f2nIdx_p[lastFace] + 1; i_node < f2nIdx_p[lastFace + 1]; i_node++) {
+                const mcIdType node = f2n_p[i_node];
+                // If this node was modified
+                if (mapO2N1.find(node) != mapO2N1.end()) {
+                    const auto& new_node = mapO2N1.at(node);
+                    f2n_pw[i_node] = new_node;
+                    noNodesAreChanged = false;
+                }
+            }
+        }
+
+        if (noNodesAreChanged && grpMustBeFullyDup)
+            throw INTERP_KERNEL::Exception(
+                "duplicateFaces: A face was supposed to be duplicated but could"
+                "not be. Please be aware that all M1 groups cannot be"
+                "duplicated, at least two adjecent inner faces are needed to"
+                "duplicate a node.");
+
+        lastFace++;
+    }
+}
+
+CrackAlgo::Set
+CrackAlgo::GetCellsTouchingNodesToDup(
+    const MEDCouplingUMesh * mf,
+    const DataArrayIdType * n2cIdx,
+    const DataArrayIdType * n2c,
+    const DataArrayIdType * f2dup
+) {
+    Set res{};
+
+    const DataArrayIdType* f2nIdx = mf->getNodalConnectivityIndex();
+    const DataArrayIdType* f2n = mf->getNodalConnectivity();
+
+    const mcIdType * f2nIdx_p {f2nIdx->begin()};
+    const mcIdType * f2n_p {f2n->begin()};
+    const mcIdType * n2cIdx_p {n2cIdx->begin()};
+    const mcIdType * n2c_p {n2c->begin()};
+
+    for (const auto face : *f2dup)
+        for (mcIdType nodeIdx {f2nIdx_p[face] + 1}; nodeIdx < f2nIdx_p[face + 1]; nodeIdx++) {
+            const auto &node = f2n_p[nodeIdx];
+            for (mcIdType cellIdx = n2cIdx_p[node]; cellIdx < n2cIdx_p[node + 1]; cellIdx++) {
+                const auto &cell = n2c_p[cellIdx];
+                res.insert(cell);
+            }
+        }
+    return res;
+}
+
+CrackAlgo::Map2Set CrackAlgo::GetNode2CellMap(const MEDCouplingUMesh * mf,
+        const DataArrayIdType * n2cIdx,
+        const DataArrayIdType * n2c,
+        const DataArrayIdType * f2dup
+) {
+    const DataArrayIdType* f2nIdx = mf->getNodalConnectivityIndex();
+    const DataArrayIdType* f2n = mf->getNodalConnectivity();
+
+    Map2Set n2c_dup{};
+    const mcIdType * f2nIdx_p {f2nIdx->begin()};
+    const mcIdType * f2n_p {f2n->begin()};
+    const mcIdType * n2cIdx_p {n2cIdx->begin()};
+    const mcIdType * n2c_p {n2c->begin()};
+
+    for (const auto face : *f2dup)
+        for (mcIdType nodeIdx {f2nIdx_p[face] + 1}; nodeIdx < f2nIdx_p[face + 1]; nodeIdx++) {
+            const auto &node = f2n_p[nodeIdx];
+            for (mcIdType cellIdx {n2cIdx_p[node]}; cellIdx < n2cIdx_p[node + 1]; cellIdx++) {
+                const auto &cell = n2c_p[cellIdx];
+                n2c_dup[node].insert(cell);
+            }
+        }
+    return n2c_dup;
+}
+
+void
+CrackAlgo::OpenCrack(
+    MEDFileUMesh * mf,
+    const Map2Map & cellOld2NewNode,
+    const double & factor
+) {
+    if ((factor <= 0.0) || (factor >= 1.0))
+        throw INTERP_KERNEL::Exception("factor should be between 0.0 and 1.0");
+    DataArrayDouble * coords = mf->getCoords();
+    const auto dim = static_cast<mcIdType>(coords->getNumberOfComponents());
+    const double * coords_p = coords->begin();
+    MCU m0 = mf->getMeshAtLevel(0);
+    DAD barys = m0->computeCellCenterOfMass();
+    const double * barys_p = barys->begin();
+    for (const auto & cellPair : cellOld2NewNode) {
+        const auto & cell = cellPair.first;
+        const auto & mapO2N = cellPair.second;
+        for (const auto & nodeO2N : mapO2N) {
+            const auto & oldN = nodeO2N.first;
+            const auto & newN = nodeO2N.second;
+            double * pos_new_p = coords->rwBegin() + dim * newN;
+            for (int i = 0; i < dim; i++) {
+                pos_new_p[i] += (1.0 - factor) * (barys_p[cell*dim+i] - coords_p[dim*oldN+i]);
+            }
+        }
+    }
+}
+
+DataArrayIdType *
+CrackAlgo::GetFacesToDupInM1(
+    const MEDCouplingUMesh & crackMesh,
+    const MEDCouplingUMesh & m1
+) {
+    DataArrayIdType * f2dupIdInM1P{};
+    const bool allIn = m1.areCellsIncludedIn(&crackMesh, 2, f2dupIdInM1P);
+    if (!allIn)
+        throw INTERP_KERNEL::Exception(
+            "crackAlong: cleanCrackMesh is not part of m1 mesh.");
+    return f2dupIdInM1P;
+}
+
+pair<DataArrayIdType*, DataArrayIdType*>
+CrackAlgo::GetFacesInM1TouchingDuplicatedNodes(
+    const Map2Set & n2c_dup,
+    const DataArrayIdType * f2dupIdInM1,
+    const MEDCouplingUMesh & mf,
+    const MEDCouplingUMesh & m1
+) {
+    Set f2dup_set = DataArrayToSet(*f2dupIdInM1);
+
+    DAI n2f(DataArrayIdType::New());
+    DAI n2fIdx(DataArrayIdType::New());
+    m1.getReverseNodalConnectivity(n2f, n2fIdx);
+
+    // 1. I retrieve all the faces of m1 which have a potentially duplicated
+    // node, I just use the reverse nodal connectivity of m1 for this.
+    Set f2changeIdInM1_set{};
+
+    for (const auto & nodeP : n2c_dup) {
+        const auto & node = nodeP.first;
+        for (mcIdType fId = n2fIdx->begin()[node]; fId < n2fIdx->begin()[node+1]; fId++) {
+            const auto & f = n2f->begin()[fId];
+            // f is not a face to duplicate
+            if (f2dup_set.find(f) == f2dup_set.end())
+                f2changeIdInM1_set.insert(f);
+        }
+    }
+    // 2. I retrieve M1 sub mesh
+    DataArrayIdType * f2changeIdInM1 = SetToDataArray(f2changeIdInM1_set);
+    MCU part = m1.buildPartOfMySelf(f2changeIdInM1->begin(), f2changeIdInM1->end());
+    // 3. I search for the ids in mf
+    DataArrayIdType * f2changeIdInMf{};
+    const bool ok = mf.areCellsIncludedIn(part, 2, f2changeIdInMf);
+    if (!ok)
+        throw INTERP_KERNEL::Exception(
+            "Some faces in -1 mesh next to the crack are not included in "
+            "the face mesh issued from m0.");
+
+    // 4. I return the 2 DAI
+    return pair<DataArrayIdType*, DataArrayIdType*>{f2changeIdInM1, f2changeIdInMf};
+}
+
+void
+CrackAlgo::ChangeConnectivityOfM1Elements(
+    const DataArrayIdType * f2changeIdInM1,
+    const DataArrayIdType * f2changeIdInMf,
+    const Map2Map & cellOld2NewNode,
+    const DataArrayIdType * f2cIdx,
+    const DataArrayIdType * f2c,
+    MEDCouplingUMesh * m1
+) {
+    // 1. I go through the faces in question. I go through all the nodes on the
+    // face. I will look in neighboring cells to see if there have been changes
+    // in connectivity on the nodes of the face.
+    // 2. if the changes are consistent, I change in place.
+    // 3. otherwise I return a consistency error.
+
+    DataArrayIdType* f2nIdx = m1->getNodalConnectivityIndex();
+    DataArrayIdType* f2n = m1->getNodalConnectivity();
+
+    // Change connectivity size
+    const mcIdType * const f2nIdx_p = f2nIdx->begin();
+    const mcIdType * const f2n_p = f2n->begin();
+    mcIdType * const f2n_pw = f2n->rwBegin();
+
+    const mcIdType * const f2cIdx_p {f2cIdx->begin()};
+    const mcIdType * const f2c_p {f2c->begin()};
+
+    const auto f2change_size = static_cast<size_t>(f2changeIdInM1->getNumberOfTuples());
+
+    for (size_t iFace = 0; iFace < f2change_size; iFace++) {
+        const auto & faceM1 = f2changeIdInM1->begin()[iFace];
+        const auto & faceMf = f2changeIdInMf->begin()[iFace];
+
+        // Change inplace the connectivity
+        // Check that this face is an inner face, ie it is connected to two
+        // cells
+        const mcIdType d = f2cIdx_p[faceMf + 1] - f2cIdx_p[faceMf];
+
+        //   1. In original face, attaching it to cell0
+        const auto& cell0 = f2c_p[f2cIdx_p[faceMf]];
+        // const auto& cell1 = f2c_p[f2cIdx_p[faceMf] + 1];
+
+        // If nodes where changed in cell0
+        if (cellOld2NewNode.find(cell0) != cellOld2NewNode.end()) {
+            const auto & mapO2N0 = cellOld2NewNode.at(cell0);
+            for (mcIdType i_node = f2nIdx_p[faceM1] + 1; i_node < f2nIdx_p[faceM1 + 1]; i_node++) {
+                const mcIdType node = f2n_p[i_node];
+                // If this node was modified
+                if (mapO2N0.find(node) != mapO2N0.end()) {
+                    const auto& new_node = mapO2N0.at(node);
+                    if ((d == 2)
+                        && ((cellOld2NewNode.find(f2c_p[f2cIdx_p[faceMf] + 1])
+                                == cellOld2NewNode.end())
+                            || (cellOld2NewNode.at(f2c_p[f2cIdx_p[faceMf] + 1]).find(node)
+                                == cellOld2NewNode.at(f2c_p[f2cIdx_p[faceMf] + 1]).end())
+                            || (cellOld2NewNode.at(f2c_p[f2cIdx_p[faceMf] + 1]).at(node) != new_node)))
+                            throw INTERP_KERNEL::Exception("There is an"
+                                "incoherent change of connectivity between both"
+                                "cell surrounding a face touching the crack but"
+                                "not duplicated.");
+                    f2n_pw[i_node] = new_node;
+                }
+            }
+        }
+    }
+}

diff --git a/src/MEDLoader/CrackAlgo.hxx b/src/MEDLoader/CrackAlgo.hxx
new file mode 100644 (file)
index 0000000..462b771
--- /dev/null
+++ b/src/MEDLoader/CrackAlgo.hxx
@@ -0,0 +1,193 @@
+// Copyright (C) 2007-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email :
+// webmaster.salome@opencascade.com
+//
+// Author : Aymeric SONOLET (CEA/DES)
+
+#ifndef SRC_MEDLOADER_CRACKALGO_HXX_
+#define SRC_MEDLOADER_CRACKALGO_HXX_
+
+#include <unordered_map>
+#include <vector>
+#include <memory>
+#include <string>
+#include <map>
+#include <unordered_set>
+#include <utility>
+
+#include <MCType.hxx>
+
+namespace MEDCoupling {
+
+class MEDFileUMesh;
+class DataArrayIdType;
+class MEDCouplingUMesh;
+
+class CrackAlgo {
+ public:
+     using Set = std::unordered_set<mcIdType>;
+     using Map = std::unordered_map<mcIdType, mcIdType>;
+     using Graph = std::unordered_map<mcIdType, Set>;
+     using Map2Set = std::unordered_map<mcIdType, Set>;
+     using Map2Map = std::unordered_map<mcIdType, Map>;
+
+     CrackAlgo() {}
+     ~CrackAlgo() {}
+     static Map2Map
+     Compute(MEDFileUMesh* mm, const std::string& grp_name, bool grpMustBeFullyDup = true);
+
+    static void
+    OpenCrack(
+        MEDFileUMesh * mf,
+        const Map2Map & cellOld2NewNode,
+        const double & factor = 0.9);
+
+ private:
+     /* Find connected components using a node to node graph.
+      */
+     static std::vector<std::shared_ptr<std::vector<mcIdType>>>
+     FindConnectedComponents(
+         const Graph& graph);
+
+     /* Depth First Search function to find connected components.
+      */
+     static void Dfs(
+         const Graph& graph,
+         const mcIdType &node,
+         const std::size_t &componentId,
+         std::map<mcIdType, bool>& visited,
+         std::map<mcIdType, std::size_t>& componentMap);
+
+     /* Converts DataArrayIdType to Set.
+      *
+      * Usefull for doing multiple search in the Set efficiently.
+      */
+     static Set
+     DataArrayToSet(
+         const DataArrayIdType &da);
+
+    static DataArrayIdType *
+    SetToDataArray(
+        const Set& s);
+
+     static Set
+     GetCellsTouchingNodesToDup(
+         const MEDCouplingUMesh * mf,
+         const DataArrayIdType * n2cIdx,
+         const DataArrayIdType * n2c,
+         const DataArrayIdType * f2dup);
+
+     static Map2Set
+     GetNode2CellMap(
+         const MEDCouplingUMesh * mf,
+         const DataArrayIdType * n2cIdx,
+         const DataArrayIdType * n2c,
+         const DataArrayIdType * f2dup);
+
+     /* Building the cell to cell graph.
+      *
+      * This graph concerns only cells touching nodes touching faces to duplicate.
+      * The connection between cells sharing a face to dup is cut in this graph.
+      */
+     static Graph
+     BuildCutC2CGraph(
+         const DataArrayIdType * c2fIdx,
+         const DataArrayIdType * c2f,
+         const DataArrayIdType * f2cIdx,
+         const DataArrayIdType * f2c,
+         const Set & cTouchingN_dup,
+         const DataArrayIdType * f2dup);
+
+     static Map2Map
+     CreateNewNodesInTopLevelMesh(
+         const Map2Set & n2c_dup,
+         const Graph & c2c,
+         MEDCouplingUMesh * m0);
+
+     static void
+     AddMissingElementsOnLevelM1AndChangeConnectivity(
+         const DataArrayIdType * f2cIdx,
+         const DataArrayIdType * f2c,
+         const DataArrayIdType * f2dupIdInM1,
+         const DataArrayIdType * f2dupIdInMf,
+         const Map2Map & cellOld2NewNode,
+         MEDCouplingUMesh * m1,
+         bool grpMustBeFullyDup = true);
+
+    /* Create new family array for level -1, which is the extended copy of the
+     * original as new elements are appended. Caller owns newFamily DAI.
+     * Supposes that the elements to duplicate are all duplicated and appended
+     * in the order at the end of mf.
+     */
+    static DataArrayIdType *
+    CopyAndCompleteFamilyArrAtLevelM1(
+        const MEDFileUMesh * mm,
+        const MEDCouplingUMesh * mf,
+        const DataArrayIdType * f2dup);
+
+    static DataArrayIdType *
+    CopyFamilyArrAtLev0(const MEDFileUMesh * mm);
+
+    /* Manage node families.
+     *
+     * Extend the family size inplace and set new nodes family to their
+     * predecessor.
+     */
+    static void
+    CompleteFamilyArrAtNodeLevel(
+        const Map2Set& addedNodes,
+        MEDFileUMesh * mm);
+
+    /* Aggregate CellOld2NewNodes into oldNode to newNodes.
+     */
+    static Map2Set
+    BuildMap2Set(
+        const Map2Map& cellOld2NewNode);
+
+    /* Remomves cells from crackingMesh which are part of m skin.
+     *
+     * Returns a new MEDCouplingUMesh. User is responsible to delete it.
+     */
+    static MEDCouplingUMesh *
+    CleanM1Mesh(
+        const MEDCouplingUMesh & m,
+        const MEDCouplingUMesh & crackingMesh);
+
+    static std::pair<DataArrayIdType*, DataArrayIdType*>
+    GetFacesInM1TouchingDuplicatedNodes(
+        const Map2Set & n2c_dup,
+        const DataArrayIdType * f2dupIdInM1,
+        const MEDCouplingUMesh & mf,
+        const MEDCouplingUMesh & m1);
+
+    static DataArrayIdType *
+    GetFacesToDupInM1(
+        const MEDCouplingUMesh & crackMesh,
+        const MEDCouplingUMesh & m1);
+
+    static void
+    ChangeConnectivityOfM1Elements(
+        const DataArrayIdType * f2changeIdInM1,
+        const DataArrayIdType * f2changeIdInMf,
+        const Map2Map & cellOld2New,
+        const DataArrayIdType * f2cIdx,
+        const DataArrayIdType * f2c,
+        MEDCouplingUMesh * m1);
+};
+}  // namespace MEDCoupling
+#endif  // SRC_MEDLOADER_CRACKALGO_HXX_

diff --git a/src/MEDLoader/MEDFileMesh.cxx b/src/MEDLoader/MEDFileMesh.cxx
index 7a0bd5f9728f72ec8abc9f7e71718b1fa7b102f7..59791f7cffdaf9811f8c46ef6255c36a27fcba53 100644 (file)
--- a/src/MEDLoader/MEDFileMesh.cxx
+++ b/src/MEDLoader/MEDFileMesh.cxx
@@ -25,6 +25,7 @@
 #include "MEDLoaderNS.hxx"
 #include "MEDFileSafeCaller.txx"
 #include "MEDLoaderBase.hxx"
+#include "CrackAlgo.hxx"
 
 #include "MEDCouplingUMesh.hxx"
 #include "MEDCouplingMappedExtrudedMesh.hxx"
@@ -32,8 +33,13 @@
 
 #include "InterpKernelAutoPtr.hxx"
 
+#include <cstddef>
 #include <limits>
 #include <cmath>
+#include <memory>
+#include <numeric>
+#include <unordered_map>
+#include <unordered_set>
 
 // From MEDLOader.cxx TU
 extern med_geometry_type                 typmai[MED_N_CELL_FIXED_GEO];
@@ -4201,6 +4207,51 @@ void MEDFileUMesh::optimizeFamilies()
     _groups.erase(*it);
 }
 
+/**
+ * \b this must be filled at level 0 and -1, the -1 level being part
+ * of the descending connectivity of the top level. This method build a
+ * "crack", or an inner boundary, in \b this along the group of level -1 named
+ * grpNameM1. The boundary is built according to the following method:
+ *  - all nodes which are not connected anymore to its neighboring cells due to
+ *  the crack are duplicated/tripled/... as needed.
+ *  - new (-1)-level cells are built lying on those new nodes. So the
+ * edges/faces along the group are duplicated. The duplicated
+ * (-1)-level cells are added to the grpNameM1.
+ *  - connectivity modification at lev -1 according to the new connectivity at
+ *  lev 0.
+ * Note that some M1 group might lead to no some cells at level 0 being still
+ * connected despite sharing only a face in M1. The behavior of the method
+ * depends on the grpMustBeFullyDup boolean input parameter in this case.
+ * Groups at lev -1 are managed by adding the new duplicated faces to the same
+ * groups of their twin.
+ * Connectivity of mesh at lev -2 is not managed, so it could be not valid
+ * anymore.
+ * After this operation:
+ *  - a top-level cell bordering the group will loose some
+ * neighbors (typically the cell which is on the other side of the group is no
+ * more a neighbor)
+ *  - the connectivity of (part of) the top level-cells bordering the group is
+ * also modified so that some cells bordering the newly created boundary use
+ * the newly computed nodes.
+ *
+ *  \param[in] grpNameM1 name of the (-1)-level group defining the boundary
+ *  \param[in] grpMustBeFullyDup boolean indicating if the method should throw
+ *  an exeption in the case the M1 group leads to some cells at lev 0 not being
+ *  seperated. If set to false, duplicated faces are created but are identical.
+ *  They share the same connectivity.
+ *  \param[out] cell2old2newNode a dict of dict indicating for each modified
+ *  cells at lev 0 the new node connectivity.
+ */
+std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>>
+MEDFileUMesh::crackAlong(const std::string &grp_name, bool grpMustBeFullyDup) {
+    return CrackAlgo::Compute(this, grp_name, grpMustBeFullyDup);
+}
+
+void MEDFileUMesh::openCrack(const std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>> & c2o2nN, const double & factor) 
+{
+    CrackAlgo::OpenCrack(this, c2o2nN, factor);
+}
+
 /**
  * \b this must be filled at level 0 and -1, typically the -1 level being (part of) the descending connectivity
  * of the top level. This method build a "crack", or an inner boundary, in \b this along the group of level -1 named grpNameM1.

diff --git a/src/MEDLoader/MEDFileMesh.hxx b/src/MEDLoader/MEDFileMesh.hxx
index 91b377d8252588818e3fe514c6725f481cf2c8d1..4c2adbc4ca1b1b4113045b43b9bb1571b7f249ff 100644 (file)
--- a/src/MEDLoader/MEDFileMesh.hxx
+++ b/src/MEDLoader/MEDFileMesh.hxx
@@ -31,6 +31,9 @@
 
 #include <map>
 #include <list>
+#include <set>
+#include <memory>
+#include <unordered_map>
 
 namespace MEDCoupling
 {
@@ -363,6 +366,8 @@ MCAuto<DataArrayDouble>& coords, MCAuto<PartDefinition>& partCoords, MCAuto<Data
     MEDLOADER_EXPORT void setGroupsOnSetMesh(int meshDimRelToMax, const std::vector<const MEDCouplingUMesh *>& ms, bool renum=false);
     MEDLOADER_EXPORT void optimizeFamilies();
     // tools
+    MEDLOADER_EXPORT std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>> crackAlong(const std::string &grpNameM1, bool grpMustBeFullyDup=true);
+    MEDLOADER_EXPORT void openCrack(const std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>> & c2o2nN, const double & factor);
     MEDLOADER_EXPORT void buildInnerBoundaryAlongM1Group(const std::string& grpNameM1, DataArrayIdType *&nodesDuplicated, DataArrayIdType *&cellsModified, DataArrayIdType *&cellsNotModified);
     MEDLOADER_EXPORT bool unPolyze(std::vector<mcIdType>& oldCode, std::vector<mcIdType>& newCode, DataArrayIdType *& o2nRenumCell);
     MEDLOADER_EXPORT DataArrayIdType *zipCoords();
@@ -398,6 +403,17 @@ MCAuto<DataArrayDouble>& coords, MCAuto<PartDefinition>& partCoords, MCAuto<Data
     void changeFamilyIdArr(mcIdType oldId, mcIdType newId);
     std::list< MCAuto<DataArrayIdType> > getAllNonNullFamilyIds() const;
     MCAuto<MEDFileUMeshSplitL1>& checkAndGiveEntryInSplitL1(int meshDimRelToMax, MEDCouplingPointSet *m);
+    static std::vector<std::shared_ptr<std::vector<mcIdType>>>
+    findConnectedComponents(
+        const std::map<mcIdType, std::set<mcIdType>> &graph
+    );
+    static void dfs(
+        const std::map<mcIdType, std::set<mcIdType>> &graph,
+        const mcIdType &node,
+        const std::size_t& componentId,
+        std::map<mcIdType, bool>& visited,
+        std::map<mcIdType, std::size_t>& componentMap
+    );
 
   private:
     static const char SPE_FAM_STR_EXTRUDED_MESH[];

diff --git a/src/MEDLoader/Swig/CMakeLists.txt b/src/MEDLoader/Swig/CMakeLists.txt
index 2cc09ced6c1f8a28480fcc6677066af3b9c76d96..479afcf421c314de23543ae5883ceeeb83a92020 100644 (file)
--- a/src/MEDLoader/Swig/CMakeLists.txt
+++ b/src/MEDLoader/Swig/CMakeLists.txt
@@ -99,7 +99,7 @@ SALOME_INSTALL_SCRIPTS(${CMAKE_CURRENT_BINARY_DIR}/MEDLoader.py ${MEDCOUPLING_IN
 INSTALL(FILES MEDLoaderSplitter.py MEDLoaderFinalize.py DESTINATION ${MEDCOUPLING_INSTALL_PYTHON})
 
 
-INSTALL(FILES ${ALL_TESTS} MEDLoaderDataForTest.py MEDLoaderTest1.py MEDLoaderTest2.py MEDLoaderTest3.py DESTINATION ${MEDCOUPLING_INSTALL_SCRIPT_SCRIPTS})
+INSTALL(FILES ${ALL_TESTS} MEDLoaderDataForTest.py MEDLoaderTest1.py MEDLoaderTest2.py MEDLoaderTest3.py CrackAlongTest.py DESTINATION ${MEDCOUPLING_INSTALL_SCRIPT_SCRIPTS})
 
 INSTALL(FILES CaseIO.py CaseReader.py CaseWriter.py VTKReader.py DESTINATION ${MEDCOUPLING_INSTALL_PYTHON})
 
@@ -138,6 +138,7 @@ SET(MEDLOADER_TEST_FILES
     MEDLoaderTest1.py
     MEDLoaderTest2.py
     MEDLoaderTest3.py
+    CrackAlongTest.py
     CaseIO.py
     CaseReader.py
     CaseWriter.py

diff --git a/src/MEDLoader/Swig/CrackAlongTest.py b/src/MEDLoader/Swig/CrackAlongTest.py
new file mode 100644 (file)
index 0000000..d93ffae
--- /dev/null
+++ b/src/MEDLoader/Swig/CrackAlongTest.py
@@ -0,0 +1,648 @@
+#  -*- coding: iso-8859-1 -*-
+# Copyright (C) 2007-2024  CEA, EDF
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License, or (at your option) any later version.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author : Aymeric SONOLET (CEA)
+
+import unittest
+
+from medcoupling import (
+    DataArrayDouble,
+    DataArrayInt,
+    MEDCouplingCMesh,
+    MEDCouplingUMesh,
+)
+from MEDLoader import MEDFileUMesh
+
+from MEDLoaderDataForTest import WriteInTmpDir
+
+
+class CrackAlongTest(unittest.TestCase):
+    @WriteInTmpDir
+    def testBuildInnerBoundaryAlongM1Group1(self):
+        fname = "Pyfile44.med"
+        m = MEDCouplingCMesh.New()
+        m.setCoordsAt(0, DataArrayDouble.New([0.0, 1.1, 2.3, 3.6, 5.0, 6.5]))
+        m.setCoordsAt(1, DataArrayDouble.New([0.0, 1.1, 2.3, 3.6, 5.0]))
+        m = m.buildUnstructured()
+        m.setName("AnthonyDuplicate")
+        m.getCoords().setInfoOnComponents(["X [km]", "Z [mm]"])
+        m2 = m.buildDescendingConnectivity()[0][
+            [8, 11, 14, 20, 21, 22, 23, 24, 25, 26, 31, 32, 33, 34, 35, 36, 37]
+        ]
+        m2.setName(m.getName())
+        grp = DataArrayInt.New([4, 6, 8])
+        grp.setName("Grp")
+        grp2 = DataArrayInt.New([9, 16])
+        grp2.setName("Grp2")
+        mm = MEDFileUMesh.New()
+        mm.setMeshAtLevel(0, m)
+        mm.setMeshAtLevel(-1, m2)
+        mm.setGroupsAtLevel(-1, [grp, grp2])
+        grpNode = DataArrayInt.New([4, 21, 23])
+        grpNode.setName("GrpNode")
+        mm.setGroupsAtLevel(1, [grpNode])
+        ref0 = [4, 15, 14, 20, 21, 4, 16, 15, 21, 22, 4, 17, 16, 22, 23]
+        ref1 = [4, 9, 8, 14, 15, 4, 10, 9, 15, 16, 4, 11, 10, 16, 17]
+        ref2 = [4, 32, 14, 20, 21, 4, 31, 32, 21, 22, 4, 30, 31, 22, 23]
+        #
+        self.assertEqual(30, mm.getNumberOfNodes())
+        self.assertEqual(
+            ref0, mm.getMeshAtLevel(0)[[12, 13, 14]].getNodalConnectivity().getValues()
+        )
+        self.assertEqual(
+            ref1, mm.getMeshAtLevel(0)[[7, 8, 9]].getNodalConnectivity().getValues()
+        )
+        #
+        c2o2nN = mm.crackAlong("Grp")
+        self.assertEqual(
+            {12: {15: 32}, 13: {15: 32, 16: 31}, 14: {16: 31, 17: 30}}, c2o2nN
+        )
+        self.assertEqual(33, mm.getNumberOfNodes())
+        self.assertEqual([4, 6, 8, 17, 18, 19], mm.getGroupArr(-1, "Grp").getValues())
+        self.assertEqual([9, 16], mm.getGroupArr(-1, "Grp2").getValues())
+        self.assertEqual([4, 21, 23], mm.getGroupArr(1, "GrpNode").getValues())
+        self.assertEqual(
+            ref2, mm.getMeshAtLevel(0)[[12, 13, 14]].getNodalConnectivity().getValues()
+        )  # cells 7,8,9 and 12,13,14 are lying on "Grp" but only 12,13,14 are renumbered
+        self.assertEqual(
+            ref1, mm.getMeshAtLevel(0)[[7, 8, 9]].getNodalConnectivity().getValues()
+        )  #
+        # fmt: off
+        refValues = DataArrayDouble.New(
+            [1.21, 1.32, 1.43, 1.54, 1.65, 1.32, 1.44, 1.56, 1.68, 1.80,
+             1.43, 1.56, 1.69, 1.82, 1.95, 1.54, 1.68, 1.82, 1.96, 2.10]
+        )
+        # fmt: on
+        valsToTest = mm.getMeshAtLevel(0).getMeasureField(True).getArray()
+        delta = valsToTest - refValues
+        delta.abs()
+        self.assertTrue(delta.getMaxValue()[0] < 1e-12)
+
+    @WriteInTmpDir
+    def testBuildInnerBoundaryAlongM1Group2(self):
+        fname = "Pyfile45.med"
+        m = MEDCouplingCMesh.New()
+        m.setCoordsAt(0, DataArrayDouble.New([0.0, 1.1, 2.3, 3.6, 5.0, 6.5]))
+        m.setCoordsAt(1, DataArrayDouble.New([0.0, 1.1, 2.3, 3.6, 5.0]))
+        m = m.buildUnstructured()
+        m.setName("AnthonyDuplicate")
+        m.getCoords().setInfoOnComponents(["X [km]", "Z [mm]"])
+        m2 = m.buildDescendingConnectivity()[0][
+            [8, 11, 14, 20, 21, 22, 23, 24, 25, 26, 31, 32, 33, 34, 35, 36, 37]
+        ]
+        m2.setName(m.getName())
+        grp = DataArrayInt.New([4, 6])
+        grp.setName("Grp")
+        grp2 = DataArrayInt.New([9, 16])
+        grp2.setName("Grp2")
+        mm = MEDFileUMesh.New()
+        mm.setMeshAtLevel(0, m)
+        mm.setMeshAtLevel(-1, m2)
+        mm.setGroupsAtLevel(-1, [grp, grp2])
+        grpNode = DataArrayInt.New([4, 21, 23])
+        grpNode.setName("GrpNode")
+        mm.setGroupsAtLevel(1, [grpNode])
+        ref0 = [4, 15, 14, 20, 21, 4, 16, 15, 21, 22, 4, 17, 16, 22, 23]
+        ref1 = [4, 9, 8, 14, 15, 4, 10, 9, 15, 16]
+        ref2 = [4, 30, 14, 20, 21, 4, 16, 30, 21, 22, 4, 17, 16, 22, 23]
+        #
+        self.assertEqual(30, mm.getNumberOfNodes())
+        self.assertEqual(
+            ref0, mm.getMeshAtLevel(0)[[12, 13, 14]].getNodalConnectivity().getValues()
+        )
+        self.assertEqual(
+            ref1, mm.getMeshAtLevel(0)[[7, 8]].getNodalConnectivity().getValues()
+        )
+        #
+        c2o2nN = mm.crackAlong("Grp")
+        self.assertEqual({13: {15: 30}, 12: {15: 30}}, c2o2nN)
+        self.assertEqual(31, mm.getNumberOfNodes())
+        self.assertEqual([4, 6, 17, 18], mm.getGroupArr(-1, "Grp").getValues())
+        self.assertEqual([9, 16], mm.getGroupArr(-1, "Grp2").getValues())
+        self.assertEqual([4, 21, 23], mm.getGroupArr(1, "GrpNode").getValues())
+        self.assertEqual(
+            ref2, mm.getMeshAtLevel(0)[[12, 13, 14]].getNodalConnectivity().getValues()
+        )  # cells 7,8,9 and 12,13,14 are lying on "Grp" but only 12 and 13 are renumbered
+        self.assertEqual(
+            ref1, mm.getMeshAtLevel(0)[[7, 8]].getNodalConnectivity().getValues()
+        )
+        # fmt: off
+        refValues = DataArrayDouble.New(
+            [1.21, 1.32, 1.43, 1.54, 1.65, 1.32, 1.44, 1.56, 1.68, 1.80,
+             1.43, 1.56, 1.69, 1.82, 1.95, 1.54, 1.68, 1.82, 1.96, 2.10]
+        )
+        # fmt: on
+        valsToTest = mm.getMeshAtLevel(0).getMeasureField(True).getArray()
+        delta = valsToTest - refValues
+        delta.abs()
+        self.assertTrue(delta.getMaxValue()[0] < 1e-12)
+
+    @WriteInTmpDir
+    def testBuildInnerBoundaryAlongM1Group3(self):
+        """Test buildInnerBoundaryAlongM1Group() with *non-connex* cracks"""
+        fname = "Pyfile73.med"
+        m = MEDCouplingCMesh.New()
+        m.setCoordsAt(0, DataArrayDouble([0.0, 1.1, 2.3, 3.6, 5.0]))
+        m.setCoordsAt(1, DataArrayDouble([0.0, 1.0, 2.0]))
+        m = m.buildUnstructured()
+        m.setName("simple")
+        m2 = m.buildDescendingConnectivity()[0]
+        m2.setName(m.getName())
+
+        # A crack in two non connected parts of the mesh:
+        grpSeg = DataArrayInt([3, 19])
+        grpSeg.setName("Grp")
+
+        mm = MEDFileUMesh.New()
+        mm.setMeshAtLevel(0, m)
+        mm.setMeshAtLevel(-1, m2)
+        mm.setGroupsAtLevel(-1, [grpSeg])
+        c2o2nN = mm.crackAlong("Grp")
+        self.assertEqual({1: {1: 16}, 7: {13: 15}}, c2o2nN)
+        # self.assertEqual([1, 13], nodes.getValues())
+        # self.assertEqual([0, 6], cellsMod.getValues())
+        # self.assertEqual([1, 7], cellsNotMod.getValues())
+        self.assertEqual(17, mm.getNumberOfNodes())
+        self.assertEqual({3, 19, 22, 23}, set(mm.getGroupArr(-1, "Grp").getValues()))
+
+        refValues = DataArrayDouble([1.1, 1.2, 1.3, 1.4, 1.1, 1.2, 1.3, 1.4])
+        valsToTest = mm.getMeshAtLevel(0).getMeasureField(True).getArray()
+        delta = valsToTest - refValues
+        delta.abs()
+        self.assertTrue(delta.getMaxValue()[0] < 1e-10)
+
+    @WriteInTmpDir
+    def testBuildInnerBoundaryAlongM1Group4(self):
+        """Test case where cells touch the M1 group on some nodes only and not
+        on full edges (triangle mesh for ex)
+        """
+        # fmt: off
+        coo = DataArrayDouble(
+            [0.0, 0.0, 1.0, 0.0, 2.0, 0.0, 3.0, 0.0, 0.0, 1.0, 1.0, 1.0,
+             2.0, 1.0, 3.0, 1.0, 0.0, 2.0, 1.0, 2.0, 2.0, 2.0, 3.0, 2.0],
+            12,
+            2,
+        )
+        conn = [3, 0, 4, 1,
+                3, 1, 4, 5,
+                3, 5, 9, 10,
+                3, 5, 10, 6,
+                3, 2, 6, 7,
+                3, 2, 7, 3,
+                3, 4, 8, 9,
+                3, 4, 9, 5,
+                3, 1, 5, 6,
+                3, 1, 6, 2,
+                3, 6, 10, 11,
+                3, 6, 11, 7]
+        # fmt: on
+        # Only TRI3:
+        connI = DataArrayInt()
+        connI.alloc(13, 1)
+        connI.iota()
+        connI *= 4
+        m2 = MEDCouplingUMesh("2D", 2)
+        m2.setCoords(coo)
+        m2.setConnectivity(DataArrayInt(conn), connI)
+        m2.checkConsistency()
+        m1, _, _, _, _ = m2.buildDescendingConnectivity()
+        grpIds = DataArrayInt([9, 11])
+        grpIds.setName("group")
+        grpIds2 = DataArrayInt([0, 1])
+        grpIds2.setName("group2")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m2)
+        mfu.setMeshAtLevel(-1, m1)
+        mfu.setGroupsAtLevel(-1, [grpIds, grpIds2])
+        nNod = m2.getNumberOfNodes()
+
+        ref0 = [3, 5, 10, 6, 3, 2, 7, 3, 3, 6, 10, 11]
+        ref1 = [3, 2, 6, 7, 3, 1, 5, 6, 3, 1, 6, 2, 3, 6, 11, 7]
+        ref2 = [3, 2, 13, 7, 3, 1, 5, 13, 3, 1, 13, 2, 3, 6, 11, 12]
+
+        self.assertEqual(ref0, m2[[3, 5, 10]].getNodalConnectivity().getValues())
+        self.assertEqual(ref1, m2[[4, 8, 9, 11]].getNodalConnectivity().getValues())
+
+        c2o2nN = mfu.crackAlong("group")
+        self.assertEqual({9: {6: 13}, 8: {6: 13}, 4: {6: 13}, 11: {7: 12}}, c2o2nN)
+
+        self.assertEqual(ref0, m2[[3, 5, 10]].getNodalConnectivity().getValues())
+        self.assertEqual(ref2, m2[[4, 8, 9, 11]].getNodalConnectivity().getValues())
+
+        m2_bis = mfu.getMeshAtLevel(0)
+        m2_bis.checkConsistency()
+        m1_bis = mfu.getMeshAtLevel(-1)
+        m1_bis.checkConsistency()
+        self.assertEqual(nNod + 2, mfu.getNumberOfNodes())
+        self.assertEqual(nNod + 2, m2_bis.getNumberOfNodes())
+        self.assertEqual(nNod + 2, m1_bis.getNumberOfNodes())
+        self.assertEqual([9, 11, 23, 24], mfu.getGroupArr(-1, "group").getValues())
+        self.assertEqual([0, 1], mfu.getGroupArr(-1, "group2").getValues())
+
+        m_bis0 = mfu.getMeshAtLevel(-1)
+        m_desc, _, _, _, _ = m_bis0.buildDescendingConnectivity()
+        m_bis0.checkDeepEquivalOnSameNodesWith(mfu.getMeshAtLevel(-1), 2, 9.9999999)
+
+    @WriteInTmpDir
+    def testBuildInnerBoundary5(self):
+        """Full 3D test with tetras only. In this case a tri from the group is not duplicated because it is made only
+        of non duplicated nodes. The tri in question is hence not part of the final new "dup" group."""
+        # fmt: off
+        coo = DataArrayDouble(
+            [200.0, 200.0, 0.0, 200.0, 200.0, 200.0, 200.0, 0.0, 200.0, 200.0, 0.0, 0.0, 0.0, 200.0, 0.0, 0.0, 200.0, 200.0, 0.0, 0.0, 0.0, 0.0, 0.0, 200.0, 400.0, 200.0, 0.0, 400.0, 200.0, 200.0, 400.0, 0.0, 0.0, 400.0, 0.0, 200.0, 0.0, 100.00000000000016, 200.0, 63.15203310314546, 200.0, 200.0, 134.45205700643342, 200.0, 200.0, 200.0, 100.00000000000016, 200.0, 63.15203310314546, 0.0, 200.0, 134.45205700643342, 0.0, 200.0, 0.0, 100.00000000000016, 0.0, 63.15203310314546, 200.0, 0.0, 134.45205700643342, 200.0, 0.0, 200.0, 100.00000000000016, 0.0, 63.15203310314546, 0.0, 0.0, 134.45205700643342, 0.0, 0.0, 200.0, 200.0, 100.02130053568538, 0.0, 200.0, 100.00938163175135, 200.0, 0.0, 100.02130053568538, 0.0, 0.0, 100.00938163175135, 299.3058739933347, 200.0, 200.0, 400.0, 98.68100542924483, 200.0, 302.8923433403344, 0.0, 200.0, 302.8923433403344, 200.0, 0.0, 400.0, 100.00000000000016, 0.0, 302.8923433403344, 0.0, 0.0, 400.0, 200.0, 98.55126825835082, 400.0, 0.0, 100.02162286181577, 99.31624553977466, 99.99999998882231, 200.0, 99.31624576683302, 100.00000010178034, 0.0, 99.31624560596512, 200.0, 100.0050761312483, 99.31624560612883, 0.0, 100.00507613125338, 200.0, 99.99999995813045, 100.00950673487786, 0.0, 99.99999989928207, 100.0041870621175, 301.29063354383015, 100.0000000093269, 0.0, 301.29063360689975, 0.0, 100.00957769061164, 140.52853868782435, 99.99999963972768, 100.00509135751312, 297.87779091770784, 97.16750463405486, 97.18018457127863],
+            46,
+            3,
+        )
+        c0 = [14, 45, 31, 21, 42, 14, 37, 38, 20, 44, 14, 39, 36, 41, 44, 14, 5, 25, 12, 13, 14, 38, 36, 44, 41, 14, 21, 20, 24, 44, 14, 38, 25, 41, 19, 14, 37, 38, 44, 41, 14, 16, 27, 39, 41, 14, 21, 45, 26, 40, 14, 39, 37, 44, 41, 14, 14, 15, 24, 44, 14, 25, 38, 41, 13, 14, 27, 18, 6, 22, 14, 38, 36, 41, 13, 14, 44, 14, 15, 36, 14, 44, 23, 39, 26, 14, 21, 26, 23, 44, 14, 38, 44, 14, 24, 14, 39, 37, 41, 22, 14, 21, 33, 45, 42, 14, 27, 22, 39, 41, 14, 23, 26, 21, 3, 14, 27, 18, 22, 41, 14, 39, 36, 44, 17, 14, 21, 26, 44, 40, 14, 39, 37, 22, 23, 14, 37, 38, 41, 19, 14, 25, 12, 13, 41, 14, 30, 26, 43, 45, 14, 38, 36, 13, 14, 14, 12, 36, 13, 41, 14, 20, 44, 21, 37, 14, 16, 36, 12, 41, 14, 39, 36, 17, 16, 14, 44, 20, 24, 38, 14, 27, 16, 12, 41, 14, 26, 15, 17, 44, 14, 19, 18, 41, 37, 14, 40, 45, 26, 15, 14, 37, 38, 19, 20, 14, 17, 15, 26, 2, 14, 39, 36, 16, 41, 14, 24, 21, 44, 40, 14, 16, 7, 27, 12, 14, 22, 18, 37, 41, 14, 21, 31, 45, 24, 14, 44, 40, 15, 24, 14, 24, 45, 15, 28, 14, 44, 40, 26, 15, 14, 24, 20, 21, 0, 14, 38, 36, 14, 44, 14, 39, 37, 23, 44, 14, 45, 31, 42, 32, 14, 25, 18, 19, 4, 14, 36, 44, 17, 15, 14, 25, 19, 18, 41, 14, 24, 15, 14, 1, 14, 45, 24, 34, 28, 14, 35, 45, 30, 43, 14, 17, 44, 39, 26, 14, 44, 23, 21, 37, 14, 30, 45, 29, 15, 14, 45, 35, 33, 43, 14, 30, 15, 26, 45, 14, 31, 21, 0, 24, 14, 33, 35, 32, 10, 14, 29, 45, 34, 28, 14, 32, 45, 34, 29, 14, 45, 31, 32, 34, 14, 33, 26, 45, 43, 14, 45, 31, 34, 24, 14, 33, 26, 21, 45, 14, 11, 30, 35, 29, 14, 33, 35, 45, 32, 14, 33, 45, 42, 32, 14, 32, 8, 34, 31, 14, 21, 26, 33, 3, 14, 35, 45, 32, 29, 14, 29, 34, 9, 28, 14, 15, 45, 24, 40, 14, 29, 45, 28, 15, 14, 21, 24, 45, 40, 14, 24, 15, 1, 28, 14, 35, 45, 29, 30, 14, 26, 15, 30, 2]
+        cI0 = [0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75, 80, 85, 90, 95, 100, 105, 110, 115, 120, 125, 130, 135, 140, 145, 150, 155, 160, 165, 170, 175, 180, 185, 190, 195, 200, 205, 210, 215, 220, 225, 230, 235, 240, 245, 250, 255, 260, 265, 270, 275, 280, 285, 290, 295, 300, 305, 310, 315, 320, 325, 330, 335, 340, 345, 350, 355, 360, 365, 370, 375, 380, 385, 390, 395, 400, 405, 410, 415, 420, 425, 430]
+        # fmt: on
+        m3 = MEDCouplingUMesh("3D", 3)
+        m3.setCoords(coo)
+        m3.setConnectivity(DataArrayInt(c0), DataArrayInt(cI0))
+        m3.checkConsistency()
+        m2, _, _, _, _ = m3.buildDescendingConnectivity()
+        grpIds = DataArrayInt([36, 74])
+        grpIds.setName("group")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m3)
+        mfu.setMeshAtLevel(-1, m2)
+        grpIds3D = DataArrayInt([0, 1])
+        grpIds3D.setName("group_3d")
+        mfu.setGroupsAtLevel(0, [grpIds3D])  # just to check preservation of 3D group
+        mfu.setGroupsAtLevel(-1, [grpIds])
+        nNod = m3.getNumberOfNodes()
+
+        c2o2nN = mfu.crackAlong("group", grpMustBeFullyDup=False)
+        self.assertEqual({77: {3: 46}}, c2o2nN)
+
+        m3_bis = mfu.getMeshAtLevel(0)
+        m3_bis.checkConsistency()
+        m2_bis = mfu.getMeshAtLevel(-1)
+        m2_bis.checkConsistency()
+        self.assertEqual(nNod + 1, mfu.getNumberOfNodes())
+        self.assertEqual(nNod + 1, m3_bis.getNumberOfNodes())
+        self.assertEqual(nNod + 1, m2_bis.getNumberOfNodes())
+        self.assertEqual(
+            m3_bis.getCoords()[3].getValues(), m3_bis.getCoords()[nNod:].getValues()
+        )
+        self.assertEqual([36, 74, 213, 214], mfu.getGroupArr(-1, "group").getValues())
+        self.assertEqual([0, 1], mfu.getGroupArr(0, "group_3d").getValues())
+        m_bis0 = mfu.getMeshAtLevel(-1)
+        m_desc, _, _, _, _ = m_bis0.buildDescendingConnectivity()
+        m_bis0.checkDeepEquivalOnSameNodesWith(mfu.getMeshAtLevel(-1), 2, 9.9999999)
+
+    @WriteInTmpDir
+    def testBuildInnerBoundary6(self):
+        """3D test where the crack has a funny shape with a singular point (i.e. two faces of the M1 group are only connected by one point, not a full segment)
+        The singular point was wrongly duplicated.
+        """
+        # fmt: off
+        coo = DataArrayDouble(
+            [(-1.38778e-17, 0.226, 0), (-1.38778e-17, -1.38778e-17, 0), (0.226, 0.226, 0), (0.226, -1.38778e-17, 0), (0.452, 0.226, 0), (0.452, -1.38778e-17, 0), (-1.38778e-17, 0.452, 0), (0.226, 0.452, 0), (0.452, 0.452, 0), (-1.38778e-17, 0.226, 0.25), (0.226, 0.226, 0.25), (0.226, -1.38778e-17, 0.25), (-1.38778e-17, -1.38778e-17, 0.25), (-1.38778e-17, 0.226, 0.779375), (0.226, 0.226, 0.779375), (0.226, -1.38778e-17, 0.779375), (-1.38778e-17, -1.38778e-17, 0.779375), (-1.38778e-17, 0.226, 1.30875), (0.226, 0.226, 1.30875), (0.226, -1.38778e-17, 1.30875), (-1.38778e-17, -1.38778e-17, 1.30875), (0.452, 0.226, 0.25), (0.452, -1.38778e-17, 0.25), (0.452, 0.226, 0.779375), (0.452, -1.38778e-17, 0.779375), (0.452, 0.226, 1.30875), (0.452, -1.38778e-17, 1.30875), (-1.38778e-17, 0.452, 0.25), (0.226, 0.452, 0.25), (-1.38778e-17, 0.452, 0.779375), (0.226, 0.452, 0.779375), (-1.38778e-17, 0.452, 1.30875), (0.226, 0.452, 1.30875), (0.452, 0.452, 0.25), (0.452, 0.452, 0.779375), (0.452, 0.452, 1.30875), (0.146, 0.226, 0.779375), (0.146, -1.38778e-17, 0.779375), (0.146, 0.226, 1.30875), (0.146, -1.38778e-17, 1.30875), (0.146, 0.452, 0.779375), (0.146, 0.452, 1.30875)]
+        )
+        c0 = [18, 0, 2, 3, 1, 9, 10, 11, 12, 18, 9, 10, 11, 12, 13, 36, 37, 16, 18, 13, 36, 37, 16, 17, 38, 39, 20, 18, 2, 4, 5, 3, 10, 21, 22, 11, 18, 10, 21, 22, 11, 14, 23, 24, 15, 18, 14, 23, 24, 15, 18, 25, 26, 19, 18, 6, 7, 2, 0, 27, 28, 10, 9, 18, 27, 28, 10, 9, 29, 40, 36, 13, 18, 29, 40, 36, 13, 31, 41, 38, 17, 18, 7, 8, 4, 2, 28, 33, 21, 10, 18, 28, 33, 21, 10, 30, 34, 23, 14, 18, 30, 34, 23, 14, 32, 35, 25, 18]
+        # fmt: on
+        cI0 = [0, 9, 18, 27, 36, 45, 54, 63, 72, 81, 90, 99, 108]
+        m3 = MEDCouplingUMesh("3D", 3)
+        m3.setCoords(coo)
+        m3.setConnectivity(DataArrayInt(c0), DataArrayInt(cI0))
+        m3.checkConsistency()
+        m2, _, _, _, _ = m3.buildDescendingConnectivity()
+        grpIds = DataArrayInt([7, 12, 22, 27])
+        grpIds.setName("group")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m3)
+        mfu.setMeshAtLevel(-1, m2)
+        mfu.setGroupsAtLevel(-1, [grpIds])
+        nNod = m3.getNumberOfNodes()
+        c2o2nN = mfu.crackAlong("group")
+        self.assertEqual(
+            {
+                7: {13: 49, 36: 48},
+                8: {13: 49, 36: 48, 17: 46, 38: 45},
+                10: {23: 47, 14: 43},
+                11: {23: 47, 25: 44, 14: 43, 18: 42},
+            },
+            c2o2nN,
+        )
+        m3_bis = mfu.getMeshAtLevel(0)
+        m3_bis.checkConsistency()
+        m2_bis = mfu.getMeshAtLevel(-1)
+        m2_bis.checkConsistency()
+        self.assertEqual(nNod + 8, mfu.getNumberOfNodes())
+        self.assertEqual(nNod + 8, m3_bis.getNumberOfNodes())
+        self.assertEqual(nNod + 8, m2_bis.getNumberOfNodes())
+        # self.assertEqual([13, 14, 17, 18, 23, 25, 36, 38], nodesDup.getValues())
+        # self.assertEqual(
+        #     m3_bis.getCoords()[nodesDup].getValues(),
+        #     m3_bis.getCoords()[nNod:].getValues(),
+        # )
+        # self.assertEqual(set([1, 2, 4, 5]), set(cells1.getValues()))
+        # self.assertEqual(set([7, 8, 10, 11]), set(cells2.getValues()))
+        self.assertEqual(
+            [7, 12, 22, 27, 56, 57, 58, 59], mfu.getGroupArr(-1, "group").getValues()
+        )
+        m_desc, _, _, _, _ = m3_bis.buildDescendingConnectivity()
+        m_desc.checkDeepEquivalOnSameNodesWith(m2_bis, 2, 9.9999)
+
+    @WriteInTmpDir
+    def testBuildInnerBoundary7(self):
+        """3D test where the crack has another funny shape with another singular point (i.e. two faces of the M1 group are only connected by one point, not a full segment)
+        Once the crack is inserted, the cells on either side of the crack do not necessarily form a connex spread zone. This was not properly handled either.
+        """
+        # fmt: off
+        coo = DataArrayDouble(
+            [(5, 17, 0), (0, 17, 0), (0, 12, 0), (5, 12, 0), (15, 17, 0), (15, 12, 0), (20, 12, 0), (20, 17, 0), (20, 2, 0), (15, 2, 0), (15, -3, 0), (20, -3, 0), (5, -3, 0), (5, 2, 0), (0, -3, 0), (0, 2, 0), (5, 17, 10), (5, 17, 20), (5, 17, 30), (5, 17, 40), (0, 17, 10), (0, 17, 20), (0, 17, 30), (0, 17, 40), (0, 12, 10), (0, 12, 20), (0, 12, 30), (0, 12, 40), (5, 12, 10), (5, 12, 20), (5, 12, 30), (5, 12, 40), (15, 17, 10), (15, 17, 20), (15, 17, 30), (15, 17, 40), (15, 12, 10), (15, 12, 20), (15, 12, 30), (15, 12, 40), (20, 12, 10), (20, 12, 20), (20, 12, 30), (20, 12, 40), (20, 17, 10), (20, 17, 20), (20, 17, 30), (20, 17, 40), (20, 2, 10), (20, 2, 20), (20, 2, 30), (20, 2, 40), (15, 2, 10), (15, 2, 20), (15, 2, 30), (15, 2, 40), (15, -3, 10), (15, -3, 20), (15, -3, 30), (15, -3, 40), (20, -3, 10), (20, -3, 20), (20, -3, 30), (20, -3, 40), (5, -3, 10), (5, -3, 20), (5, -3, 30), (5, -3, 40), (5, 2, 10), (5, 2, 20), (5, 2, 30), (5, 2, 40), (0, -3, 10), (0, -3, 20), (0, -3, 30), (0, -3, 40), (0, 2, 10), (0, 2, 20), (0, 2, 30), (0, 2, 40), (20, 8, 0), (0, 8, 0), (20, 8, 10), (20, 8, 20), (20, 8, 30), (20, 8, 40), (15, 8, 30), (15, 8, 40), (5, 8, 30), (5, 8, 40), (0, 8, 10), (0, 8, 20), (0, 8, 30), (0, 8, 40)]
+        )
+        c = DataArrayInt(
+            [31, 0, 3, 2, 1, -1, 16, 20, 24, 28, -1, 0, 16, 28, 3, -1, 3, 28, 24, 2, -1, 2, 24, 20, 1, -1, 1, 20, 16, 0, 31, 16, 28, 24, 20, -1, 17, 21, 25, 29, -1, 16, 17, 29, 28, -1, 28, 29, 25, 24, -1, 24, 25, 21, 20, -1, 20, 21, 17, 16, 31, 17, 29, 25, 21, -1, 18, 22, 26, 30, -1, 17, 18, 30, 29, -1, 29, 30, 26, 25, -1, 25, 26, 22, 21, -1, 21, 22, 18, 17, 31, 18, 30, 26, 22, -1, 19, 23, 27, 31, -1, 18, 19, 31, 30, -1, 30, 31, 27, 26, -1, 26, 27, 23, 22, -1, 22, 23, 19, 18, 31, 4, 5, 3, 0, -1, 32, 16, 28, 36, -1, 4, 32, 36, 5, -1, 5, 36, 28, 3, -1, 3, 28, 16, 0, -1, 0, 16, 32, 4, 31, 32, 36, 28, 16, -1, 33, 17, 29, 37, -1, 32, 33, 37, 36, -1, 36, 37, 29, 28, -1, 28, 29, 17, 16, -1, 16, 17, 33, 32, 31, 33, 37, 29, 17, -1, 34, 18, 30, 38, -1, 33, 34, 38, 37, -1, 37, 38, 30, 29, -1, 29, 30, 18, 17, -1, 17, 18, 34, 33, 31, 34, 38, 30, 18, -1, 35, 19, 31, 39, -1, 34, 35, 39, 38, -1, 38, 39, 31, 30, -1, 30, 31, 19, 18, -1, 18, 19, 35, 34, 31, 6, 5, 4, 7, -1, 40, 44, 32, 36, -1, 6, 40, 36, 5, -1, 5, 36, 32, 4, -1, 4, 32, 44, 7, -1, 7, 44, 40, 6, 31, 40, 36, 32, 44, -1, 41, 45, 33, 37, -1, 40, 41, 37, 36, -1, 36, 37, 33, 32, -1, 32, 33, 45, 44, -1, 44, 45, 41, 40, 31, 41, 37, 33, 45, -1, 42, 46, 34, 38, -1, 41, 42, 38, 37, -1, 37, 38, 34, 33, -1, 33, 34, 46, 45, -1, 45, 46, 42, 41, 31, 42, 38, 34, 46, -1, 43, 47, 35, 39, -1, 42, 43, 39, 38, -1, 38, 39, 35, 34, -1, 34, 35, 47, 46, -1, 46, 47, 43, 42, 31, 80, 9, 5, 6, -1, 82, 40, 36, 52, -1, 80, 82, 52, 9, -1, 9, 52, 36, 5, -1, 5, 36, 40, 6, -1, 6, 40, 82, 80, 31, 82, 52, 36, 40, -1, 83, 41, 37, 53, -1, 82, 83, 53, 52, -1, 52, 53, 37, 36, -1, 36, 37, 41, 40, -1, 40, 41, 83, 82, 31, 83, 53, 37, 41, -1, 84, 42, 38, 86, -1, 83, 84, 86, 53, -1, 53, 86, 38, 37, -1, 37, 38, 42, 41, -1, 41, 42, 84, 83, 31, 84, 86, 38, 42, -1, 85, 43, 39, 87, -1, 84, 85, 87, 86, -1, 86, 87, 39, 38, -1, 38, 39, 43, 42, -1, 42, 43, 85, 84, 31, 10, 9, 8, 11, -1, 56, 60, 48, 52, -1, 10, 56, 52, 9, -1, 9, 52, 48, 8, -1, 8, 48, 60, 11, -1, 11, 60, 56, 10, 31, 56, 52, 48, 60, -1, 57, 61, 49, 53, -1, 56, 57, 53, 52, -1, 52, 53, 49, 48, -1, 48, 49, 61, 60, -1, 60, 61, 57, 56, 31, 57, 53, 49, 61, -1, 58, 62, 50, 54, -1, 57, 58, 54, 53, -1, 53, 54, 50, 49, -1, 49, 50, 62, 61, -1, 61, 62, 58, 57, 31, 58, 54, 50, 62, -1, 59, 63, 51, 55, -1, 58, 59, 55, 54, -1, 54, 55, 51, 50, -1, 50, 51, 63, 62, -1, 62, 63, 59, 58, 31, 12, 13, 9, 10, -1, 64, 56, 52, 68, -1, 12, 64, 68, 13, -1, 13, 68, 52, 9, -1, 9, 52, 56, 10, -1, 10, 56, 64, 12, 31, 64, 68, 52, 56, -1, 65, 57, 53, 69, -1, 64, 65, 69, 68, -1, 68, 69, 53, 52, -1, 52, 53, 57, 56, -1, 56, 57, 65, 64, 31, 65, 69, 53, 57, -1, 66, 58, 54, 70, -1, 65, 66, 70, 69, -1, 69, 70, 54, 53, -1, 53, 54, 58, 57, -1, 57, 58, 66, 65, 31, 66, 70, 54, 58, -1, 67, 59, 55, 71, -1, 66, 67, 71, 70, -1, 70, 71, 55, 54, -1, 54, 55, 59, 58, -1, 58, 59, 67, 66, 31, 14, 15, 13, 12, -1, 72, 64, 68, 76, -1, 14, 72, 76, 15, -1, 15, 76, 68, 13, -1, 13, 68, 64, 12, -1, 12, 64, 72, 14, 31, 72, 76, 68, 64, -1, 73, 65, 69, 77, -1, 72, 73, 77, 76, -1, 76, 77, 69, 68, -1, 68, 69, 65, 64, -1, 64, 65, 73, 72, 31, 73, 77, 69, 65, -1, 74, 66, 70, 78, -1, 73, 74, 78, 77, -1, 77, 78, 70, 69, -1, 69, 70, 66, 65, -1, 65, 66, 74, 73, 31, 74, 78, 70, 66, -1, 75, 67, 71, 79, -1, 74, 75, 79, 78, -1, 78, 79, 71, 70, -1, 70, 71, 67, 66, -1, 66, 67, 75, 74, 31, 2, 3, 13, 81, -1, 24, 90, 68, 28, -1, 2, 24, 28, 3, -1, 3, 28, 68, 13, -1, 13, 68, 90, 81, -1, 81, 90, 24, 2, 31, 24, 28, 68, 90, -1, 25, 91, 69, 29, -1, 24, 25, 29, 28, -1, 28, 29, 69, 68, -1, 68, 69, 91, 90, -1, 90, 91, 25, 24, 31, 25, 29, 69, 91, -1, 26, 92, 88, 30, -1, 25, 26, 30, 29, -1, 29, 30, 88, 69, -1, 69, 88, 92, 91, -1, 91, 92, 26, 25, 31, 26, 30, 88, 92, -1, 27, 93, 89, 31, -1, 26, 27, 31, 30, -1, 30, 31, 89, 88, -1, 88, 89, 93, 92, -1, 92, 93, 27, 26, 31, 13, 3, 5, 9, -1, 68, 52, 36, 28, -1, 13, 68, 28, 3, -1, 3, 28, 36, 5, -1, 5, 36, 52, 9, -1, 9, 52, 68, 13, 31, 68, 28, 36, 52, -1, 69, 53, 37, 29, -1, 68, 69, 29, 28, -1, 28, 29, 37, 36, -1, 36, 37, 53, 52, -1, 52, 53, 69, 68, 31, 69, 29, 37, 53, -1, 88, 86, 38, 30, -1, 69, 88, 30, 29, -1, 29, 30, 38, 37, -1, 37, 38, 86, 53, -1, 53, 86, 88, 69, 31, 88, 30, 38, 86, -1, 89, 87, 39, 31, -1, 88, 89, 31, 30, -1, 30, 31, 39, 38, -1, 38, 39, 87, 86, -1, 86, 87, 89, 88]
+        )
+        cI = DataArrayInt(
+            [0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330, 360, 390, 420, 450, 480, 510, 540, 570, 600, 630, 660, 690, 720, 750, 780, 810, 840, 870, 900, 930, 960, 990, 1020, 1050, 1080]
+        )
+        # fmt: on
+        m3 = MEDCouplingUMesh("box", 3)
+        m3.setCoords(coo)
+        m3.setConnectivity(c, cI)
+        m3.checkConsistency()
+        m2, _, _, _, _ = m3.buildDescendingConnectivity()
+        grpIds = DataArrayInt([2, 7, 12, 17, 95, 99, 103, 107, 129, 133, 137, 141])
+        grpIds.setName("group")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m3)
+        mfu.setMeshAtLevel(-1, m2)
+        mfu.setGroupsAtLevel(-1, [grpIds])
+        nNod = m3.getNumberOfNodes()
+
+        c2o2nN = mfu.crackAlong("group")
+
+        self.assertEqual(
+            {
+                7: {19: 115, 31: 114, 18: 112, 30: 111},
+                27: {70: 116, 71: 99, 66: 95, 67: 94},
+                26: {70: 116, 65: 98, 69: 96, 66: 95},
+                25: {64: 118, 68: 106, 65: 98, 69: 96},
+                29: {68: 107, 69: 97},
+                30: {69: 97},
+                24: {64: 118, 12: 117, 68: 106, 13: 100},
+                28: {68: 107, 13: 101},
+                32: {28: 108, 3: 103},
+                35: {31: 114, 89: 113, 30: 111, 88: 110},
+                6: {18: 112, 30: 111, 17: 109, 29: 105},
+                4: {28: 108, 0: 104, 3: 103, 16: 102},
+                33: {28: 108, 29: 105},
+                5: {17: 109, 28: 108, 29: 105, 16: 102},
+                34: {30: 111, 88: 110, 29: 105},
+            },
+            c2o2nN,
+        )
+        m3_bis = mfu.getMeshAtLevel(0)
+        m3_bis.checkConsistency()
+        m2_bis = mfu.getMeshAtLevel(-1)
+        m2_bis.checkConsistency()
+        self.assertEqual(nNod + 25, mfu.getNumberOfNodes())
+        self.assertEqual(nNod + 25, m3_bis.getNumberOfNodes())
+        self.assertEqual(nNod + 25, m2_bis.getNumberOfNodes())
+        # self.assertEqual(
+        #     [0, 3, 12, 13, 16, 17, 18, 19, 28, 29, 30, 31, 64, 65, 66, 67, 68, 69, 70, 71, 88, 89],
+        #     nodesDup.getValues(),
+        # )  # fmt: skip
+        # self.assertEqual(
+        #     m3_bis.getCoords()[nodesDup].getValues(),
+        #     m3_bis.getCoords()[nNod:].getValues(),
+        # )
+        # self.assertEqual(
+        #     set([0, 1, 2, 3, 24, 25, 26, 27, 28, 29, 30, 31]), set(cells1.getValues())
+        # )
+        # self.assertEqual(
+        #     set([4, 5, 6, 7, 20, 21, 22, 23, 32, 33, 34, 35]), set(cells2.getValues())
+        # )
+        self.assertEqual(
+            [2  , 7  , 12 , 17 , 95 , 99 , 103, 107, 129, 133, 137, 141,
+             151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162],
+            mfu.getGroupArr(-1, "group").getValues(),
+        )  # fmt: skip
+        m_desc, _, _, _, _ = m3_bis.buildDescendingConnectivity()
+        m_desc.checkDeepEquivalOnSameNodesWith(m2_bis, 2, 9.9999)
+
+    # def testBuildInnerBoundary8(self):
+    #     """3D test where the crack leaves 'naked' cells. If we call a 'close-to-crack cell' a cell which shares a face with the M1 group,
+    #     a 'naked cell' is a cell that has some node duplicated, but which do not share any face with a 'close-to-crack cell'. In this case
+    #     it is tricky to decide whether this cell should be renumbered or not ...
+    #     Warning: on the mesh below some points have already been doubled by a previous cut.
+    #     """
+    #     # fmt: off
+    #     coo = DataArrayDouble(
+    #         [(0, 15, 0), (0, 5, 0), (3, 5, 0), (5, 5, 0), (5, 15, 0), (5, 20, 0), (0, 20, 0), (15, 20, 0), (15, 15, 0), (20, 15, 0), (20, 20, 0), (20, 5, 0), (15, 5, 0), (15, 0, 0), (20, 0, 0), (5, -1.60551e-25, 0), (5, 3, 0), (3, 0, 0), (3, 3, 0), (0, 0, 0), (0, 3, 0), (0, 15, 10), (0, 15, 20), (0, 15, 30), (0, 15, 40), (0, 5, 10), (0, 5, 20), (0, 5, 30), (0, 5, 40), (3, 5, 10), (3, 5, 20), (3, 5, 30), (3, 5, 40), (5, 5, 10), (5, 5, 20), (5, 5, 30), (5, 5, 40), (5, 15, 10), (5, 15, 20), (5, 15, 30), (5, 15, 40), (5, 20, 10), (5, 20, 20), (5, 20, 30), (5, 20, 40), (0, 20, 10), (0, 20, 20), (0, 20, 30), (0, 20, 40), (15, 20, 10), (15, 20, 20), (15, 20, 30), (15, 20, 40), (15, 15, 10), (15, 15, 20), (15, 15, 30), (15, 15, 40), (20, 15, 10), (20, 15, 20), (20, 15, 30), (20, 15, 40), (20, 20, 10), (20, 20, 20), (20, 20, 30), (20, 20, 40), (20, 5, 10), (20, 5, 20), (20, 5, 30), (20, 5, 40), (15, 5, 10), (15, 5, 20), (15, 5, 30), (15, 5, 40), (15, 0, 10), (15, 0, 20), (15, 0, 30), (15, 0, 40), (20, 0, 10), (20, 0, 20), (20, 0, 30), (20, 0, 40), (5, -1.60551e-25, 10), (5, -1.60551e-25, 20), (5, -1.60551e-25, 30), (5, -1.60551e-25, 40), (5, 3, 10), (5, 3, 20), (5, 3, 30), (5, 3, 40), (3, 0, 10), (3, 0, 20), (3, 0, 30), (3, 0, 40), (3, 3, 10), (3, 3, 20), (3, 3, 30), (3, 3, 40), (0, 0, 10), (0, 0, 20), (0, 0, 30), (0, 0, 40), (0, 3, 10), (0, 3, 20), (0, 3, 30), (0, 3, 40), (0, 9, 0), (3, 9, 0), (20, 9, 0), (0, 9, 10), (0, 9, 20), (0, 9, 30), (0, 9, 40), (3, 9, 10), (3, 9, 20), (3, 9, 30), (3, 9, 40), (5, 9, 30), (5, 9, 40), (20, 9, 10), (20, 9, 20), (20, 9, 30), (20, 9, 40), (15, 9, 30), (15, 9, 40), (0, 15, 0), (20, 15, 0), (0, 15, 10), (0, 15, 20), (0, 15, 30), (0, 15, 40), (5, 15, 30), (5, 15, 40), (15, 15, 30), (15, 15, 40), (20, 15, 10), (20, 15, 20), (20, 15, 30), (20, 15, 40)]
+    #     )
+    #     c = DataArrayInt(
+    #         [31, 5, 4, 124, 6, -1, 41, 45, 126, 37, -1, 5, 41, 37, 4, -1, 4, 37, 126, 124, -1, 124, 126, 45, 6, -1, 6, 45, 41, 5, 31, 41, 37, 126, 45, -1, 42, 46, 127, 38, -1, 41, 42, 38, 37, -1, 37, 38, 127, 126, -1, 126, 127, 46, 45, -1, 45, 46, 42, 41, 31, 42, 38, 127, 46, -1, 43, 47, 128, 130, -1, 42, 43, 130, 38, -1, 38, 130, 128, 127, -1, 127, 128, 47, 46, -1, 46, 47, 43, 42, 31, 43, 130, 128, 47, -1, 44, 48, 129, 131, -1, 43, 44, 131, 130, -1, 130, 131, 129, 128, -1, 128, 129, 48, 47, -1, 47, 48, 44, 43, 31, 7, 8, 4, 5, -1, 49, 41, 37, 53, -1, 7, 49, 53, 8, -1, 8, 53, 37, 4, -1, 4, 37, 41, 5, -1, 5, 41, 49, 7, 31, 49, 53, 37, 41, -1, 50, 42, 38, 54, -1, 49, 50, 54, 53, -1, 53, 54, 38, 37, -1, 37, 38, 42, 41, -1, 41, 42, 50, 49, 31, 50, 54, 38, 42, -1, 51, 43, 130, 132, -1, 50, 51, 132, 54, -1, 54, 132, 130, 38, -1, 38, 130, 43, 42, -1, 42, 43, 51, 50, 31, 51, 132, 130, 43, -1, 52, 44, 131, 133, -1, 51, 52, 133, 132, -1, 132, 133, 131, 130, -1, 130, 131, 44, 43, -1, 43, 44, 52, 51, 31, 125, 8, 7, 10, -1, 134, 61, 49, 53, -1, 125, 134, 53, 8, -1, 8, 53, 49, 7, -1, 7, 49, 61, 10, -1, 10, 61, 134, 125, 31, 134, 53, 49, 61, -1, 135, 62, 50, 54, -1, 134, 135, 54, 53, -1, 53, 54, 50, 49, -1, 49, 50, 62, 61, -1, 61, 62, 135, 134, 31, 135, 54, 50, 62, -1, 136, 63, 51, 132, -1, 135, 136, 132, 54, -1, 54, 132, 51, 50, -1, 50, 51, 63, 62, -1, 62, 63, 136, 135, 31, 136, 132, 51, 63, -1, 137, 64, 52, 133, -1, 136, 137, 133, 132, -1, 132, 133, 52, 51, -1, 51, 52, 64, 63, -1, 63, 64, 137, 136, 31, 107, 12, 8, 9, -1, 118, 57, 53, 69, -1, 107, 118, 69, 12, -1, 12, 69, 53, 8, -1, 8, 53, 57, 9, -1, 9, 57, 118, 107, 31, 118, 69, 53, 57, -1, 119, 58, 54, 70, -1, 118, 119, 70, 69, -1, 69, 70, 54, 53, -1, 53, 54, 58, 57, -1, 57, 58, 119, 118, 31, 119, 70, 54, 58, -1, 120, 59, 55, 122, -1, 119, 120, 122, 70, -1, 70, 122, 55, 54, -1, 54, 55, 59, 58, -1, 58, 59, 120, 119, 31, 120, 122, 55, 59, -1, 121, 60, 56, 123, -1, 120, 121, 123, 122, -1, 122, 123, 56, 55, -1, 55, 56, 60, 59, -1, 59, 60, 121, 120, 31, 13, 12, 11, 14, -1, 73, 77, 65, 69, -1, 13, 73, 69, 12, -1, 12, 69, 65, 11, -1, 11, 65, 77, 14, -1, 14, 77, 73, 13, 31, 73, 69, 65, 77, -1, 74, 78, 66, 70, -1, 73, 74, 70, 69, -1, 69, 70, 66, 65, -1, 65, 66, 78, 77, -1, 77, 78, 74, 73, 31, 74, 70, 66, 78, -1, 75, 79, 67, 71, -1, 74, 75, 71, 70, -1, 70, 71, 67, 66, -1, 66, 67, 79, 78, -1, 78, 79, 75, 74, 31, 75, 71, 67, 79, -1, 76, 80, 68, 72, -1, 75, 76, 72, 71, -1, 71, 72, 68, 67, -1, 67, 68, 80, 79, -1, 79, 80, 76, 75, 31, 17, 18, 16, 15, -1, 89, 81, 85, 93, -1, 17, 89, 93, 18, -1, 18, 93, 85, 16, -1, 16, 85, 81, 15, -1, 15, 81, 89, 17, 31, 89, 93, 85, 81, -1, 90, 82, 86, 94, -1, 89, 90, 94, 93, -1, 93, 94, 86, 85, -1, 85, 86, 82, 81, -1, 81, 82, 90, 89, 31, 90, 94, 86, 82, -1, 91, 83, 87, 95, -1, 90, 91, 95, 94, -1, 94, 95, 87, 86, -1, 86, 87, 83, 82, -1, 82, 83, 91, 90, 31, 91, 95, 87, 83, -1, 92, 84, 88, 96, -1, 91, 92, 96, 95, -1, 95, 96, 88, 87, -1, 87, 88, 84, 83, -1, 83, 84, 92, 91, 31, 19, 20, 18, 17, -1, 97, 89, 93, 101, -1, 19, 97, 101, 20, -1, 20, 101, 93, 18, -1, 18, 93, 89, 17, -1, 17, 89, 97, 19, 31, 97, 101, 93, 89, -1, 98, 90, 94, 102, -1, 97, 98, 102, 101, -1, 101, 102, 94, 93, -1, 93, 94, 90, 89, -1, 89, 90, 98, 97, 31, 98, 102, 94, 90, -1, 99, 91, 95, 103, -1, 98, 99, 103, 102, -1, 102, 103, 95, 94, -1, 94, 95, 91, 90, -1, 90, 91, 99, 98, 31, 99, 103, 95, 91, -1, 100, 92, 96, 104, -1, 99, 100, 104, 103, -1, 103, 104, 96, 95, -1, 95, 96, 92, 91, -1, 91, 92, 100, 99, 31, 1, 2, 18, 20, -1, 25, 101, 93, 29, -1, 1, 25, 29, 2, -1, 2, 29, 93, 18, -1, 18, 93, 101, 20, -1, 20, 101, 25, 1, 31, 25, 29, 93, 101, -1, 26, 102, 94, 30, -1, 25, 26, 30, 29, -1, 29, 30, 94, 93, -1, 93, 94, 102, 101, -1, 101, 102, 26, 25, 31, 26, 30, 94, 102, -1, 27, 103, 95, 31, -1, 26, 27, 31, 30, -1, 30, 31, 95, 94, -1, 94, 95, 103, 102, -1, 102, 103, 27, 26, 31, 27, 31, 95, 103, -1, 28, 104, 96, 32, -1, 27, 28, 32, 31, -1, 31, 32, 96, 95, -1, 95, 96, 104, 103, -1, 103, 104, 28, 27, 31, 3, 4, 8, 12, -1, 33, 69, 53, 37, -1, 3, 33, 37, 4, -1, 4, 37, 53, 8, -1, 8, 53, 69, 12, -1, 12, 69, 33, 3, 31, 33, 37, 53, 69, -1, 34, 70, 54, 38, -1, 33, 34, 38, 37, -1, 37, 38, 54, 53, -1, 53, 54, 70, 69, -1, 69, 70, 34, 33, 31, 34, 38, 54, 70, -1, 116, 122, 55, 39, -1, 34, 116, 39, 38, -1, 38, 39, 55, 54, -1, 54, 55, 122, 70, -1, 70, 122, 116, 34, 31, 116, 39, 55, 122, -1, 117, 123, 56, 40, -1, 116, 117, 40, 39, -1, 39, 40, 56, 55, -1, 55, 56, 123, 122, -1, 122, 123, 117, 116, 31, 16, 18, 2, 3, -1, 85, 33, 29, 93, -1, 16, 85, 93, 18, -1, 18, 93, 29, 2, -1, 2, 29, 33, 3, -1, 3, 33, 85, 16, 31, 85, 93, 29, 33, -1, 86, 34, 30, 94, -1, 85, 86, 94, 93, -1, 93, 94, 30, 29, -1, 29, 30, 34, 33, -1, 33, 34, 86, 85, 31, 86, 94, 30, 34, -1, 87, 35, 31, 95, -1, 86, 87, 95, 94, -1, 94, 95, 31, 30, -1, 30, 31, 35, 34, -1, 34, 35, 87, 86, 31, 87, 95, 31, 35, -1, 88, 36, 32, 96, -1, 87, 88, 96, 95, -1, 95, 96, 32, 31, -1, 31, 32, 36, 35, -1, 35, 36, 88, 87, 31, 4, 3, 106, 105, 0, -1, 37, 21, 108, 112, 33, -1, 3, 4, 37, 33, -1, 106, 3, 33, 112, -1, 105, 106, 112, 108, -1, 0, 105, 108, 21, -1, 4, 0, 21, 37, 31, 37, 33, 112, 108, 21, -1, 38, 22, 109, 113, 34, -1, 33, 37, 38, 34, -1, 112, 33, 34, 113, -1, 108, 112, 113, 109, -1, 21, 108, 109, 22, -1, 37, 21, 22, 38, 31, 38, 34, 113, 109, 22, -1, 39, 23, 110, 114, 116, -1, 34, 38, 39, 116, -1, 113, 34, 116, 114, -1, 109, 113, 114, 110, -1, 22, 109, 110, 23, -1, 38, 22, 23, 39, 31, 39, 116, 114, 110, 23, -1, 40, 24, 111, 115, 117, -1, 116, 39, 40, 117, -1, 114, 116, 117, 115, -1, 110, 114, 115, 111, -1, 23, 110, 111, 24, -1, 39, 23, 24, 40, 31, 16, 3, 12, 13, 15, -1, 85, 81, 73, 69, 33, -1, 3, 16, 85, 33, -1, 12, 3, 33, 69, -1, 13, 12, 69, 73, -1, 15, 13, 73, 81, -1, 16, 15, 81, 85, 31, 85, 33, 69, 73, 81, -1, 86, 82, 74, 70, 34, -1, 33, 85, 86, 34, -1, 69, 33, 34, 70, -1, 73, 69, 70, 74, -1, 81, 73, 74, 82, -1, 85, 81, 82, 86, 31, 86, 34, 70, 74, 82, -1, 87, 83, 75, 71, 35, -1, 34, 86, 87, 35, -1, 70, 34, 35, 71, -1, 74, 70, 71, 75, -1, 82, 74, 75, 83, -1, 86, 82, 83, 87, 31, 87, 35, 71, 75, 83, -1, 88, 84, 76, 72, 36, -1, 35, 87, 88, 36, -1, 71, 35, 36, 72, -1, 75, 71, 72, 76, -1, 83, 75, 76, 84, -1, 87, 83, 84, 88]
+    #     )
+    #     cI = DataArrayInt(
+    #         [0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330, 360, 390, 420, 450, 480, 510, 540, 570, 600, 630, 660, 690, 720, 750, 780, 810, 840, 870, 900, 930, 960, 990, 1020, 1050, 1080, 1110, 1140, 1170, 1200, 1237, 1274, 1311, 1348, 1385, 1422, 1459, 1496]
+    #     )
+    #     # fmt: on
+    #     m3 = MEDCouplingUMesh("box", 3)
+    #     m3.setCoords(coo)
+    #     m3.setConnectivity(c, cI)
+    #     m3.mergeNodes(1e-12)
+    #     m3.checkConsistency()
+    #     m2, _, _, _, _ = m3.buildDescendingConnectivity()
+    #     # grpIds = DataArrayInt(
+    #     #     [2, 7, 12, 17, 101, 106, 111, 116, 160, 164, 170, 173, 176, 179]
+    #     # )
+    #     grpIds = DataArrayInt(
+    #         [16, 149, 12, 146, 8, 143, 137, 108, 134, 103, 131, 98, 128, 93]
+    #     )
+    #     grpIds.setName("group")
+    #     mfu = MEDFileUMesh()
+    #     mfu.setMeshAtLevel(0, m3)
+    #     mfu.setMeshAtLevel(-1, m2)
+    #     mfu.setGroupsAtLevel(-1, [grpIds])
+    #     nNod = m3.getNumberOfNodes()
+
+    #     mfu.write("case8dupnode_in.med", 2)
+
+    #     c2o2nN = mfu.crackAlong("group", False)
+
+    #     mfu.write("case8dupnode_out.med", 2)
+
+    #     mfu.OpenCrack(c2o2nN)
+
+    #     mfu.write("case8dupnode_cracked.med", 2)
+
+    #     self.assertEqual(
+    #         {
+    #             44: {85: 166, 15: 163, 81: 159, 16: 158},
+    #             6: {130: 162, 43: 161, 38: 160, 42: 157},
+    #             41: {38: 160, 37: 154},
+    #             47: {83: 146, 84: 143, 87: 148, 35: 139, 88: 145, 36: 138},
+    #             46: {86: 150, 82: 149, 83: 146, 87: 148, 35: 139},
+    #             7: {131: 165, 44: 164, 130: 162, 43: 161},
+    #             36: {33: 155, 3: 140},
+    #             43: {116: 147, 117: 144, 39: 142, 40: 141},
+    #             33: {38: 160, 37: 154},
+    #             4: {4: 156, 37: 154, 41: 153, 5: 152},
+    #             42: {38: 160, 116: 147, 39: 142},
+    #             45: {85: 166, 81: 159, 86: 150, 82: 149},
+    #             32: {4: 156, 37: 154},
+    #             37: {33: 155, 34: 151},
+    #             38: {34: 151},
+    #             34: {38: 160},
+    #             5: {38: 160, 42: 157, 37: 154, 41: 153},
+    #             40: {4: 156, 37: 154},
+    #         },
+    #         c2o2nN,
+    #     )
+    #     m3_bis = mfu.getMeshAtLevel(0)
+    #     m3_bis.checkConsistency()
+    #     m2_bis = mfu.getMeshAtLevel(-1)
+    #     m2_bis.checkConsistency()
+    #     self.assertEqual(nNod + 23, mfu.getNumberOfNodes())
+    #     self.assertEqual(nNod + 23, m3_bis.getNumberOfNodes())
+    #     self.assertEqual(nNod + 23, m2_bis.getNumberOfNodes())
+    #     self.assertEqual(
+    #         [5, 15, 16, 35, 36, 39, 40, 41, 42, 43, 44, 81, 82, 83, 84, 85, 86, 87, 88, 116, 117, 130, 131],
+    #         nodesDup.getValues(),
+    #     )  # fmt: skip
+    #     self.assertEqual(
+    #         m3_bis.getCoords()[nodesDup].getValues(),
+    #         m3_bis.getCoords()[nNod:].getValues(),
+    #     )
+    #     self.assertEqual(
+    #         set([0, 1, 2, 3, 20, 21, 22, 23, 34, 35, 36, 37, 38, 39]),
+    #         set(cells1.getValues()),
+    #     )
+    #     self.assertEqual(
+    #         set([4, 5, 6, 7, 42, 43, 44, 45, 46, 47]), set(cells2.getValues())
+    #     )
+    #     self.assertEqual(
+    #         [2, 7, 12, 17, 101, 106, 111, 116, 160, 164, 170, 173, 176, 179],
+    #         mfu.getGroupArr(-1, "group").getValues(),
+    #     )
+    #     self.assertEqual(
+    #         [212, 213, 214, 215, 216, 217, 218, 219, 220, 221, 222, 223, 224, 225],
+    #         mfu.getGroupArr(-1, "group_dup").getValues(),
+    #     )  # fmt: skip
+    #     m_desc, _, _, _, _ = m3_bis.buildDescendingConnectivity()
+    #     m_desc.checkDeepEquivalOnSameNodesWith(m2_bis, 2, 9.9999)
+
+    def testBuildInnerBoundary9(self):
+        """3D test where the crack is performed so that two non-connex parts are found facing one single connex part on the other side
+        of the crack.
+        """
+        # fmt: off
+        coo = DataArrayDouble(
+            [(0, 4.6, 0), (3, 4.6, 0), (5, 4.6, 0), (15, 4.6, 0), (15, 0, 0), (5, -1.60551e-25, 0), (5, 3, 0), (3, 0, 0), (3, 3.8, 0), (0, 0, 0), (0, 3.8, 0), (0, 4.6, 10), (0, 4.6, 20), (3, 4.6, 10), (3, 4.6, 20), (5, 4.6, 10), (5, 4.6, 20), (15, 4.6, 10), (15, 4.6, 20), (15, 0, 10), (15, 0, 20), (5, -1.60551e-25, 10), (5, -1.60551e-25, 20), (5, 3, 10), (5, 3, 20), (3, 0, 10), (3, 0, 20), (3, 3.8, 10), (3, 3.8, 20), (0, 0, 10), (0, 0, 20), (0, 3.8, 10), (0, 3.8, 20), (3, 3, 0), (0, 3, 0), (3, 3, 10), (3, 3, 20), (0, 3, 10), (0, 3, 20), ]
+        )
+        c = DataArrayInt(
+            [31, 7, 33, 6, 5, -1, 25, 21, 23, 35, -1, 7, 25, 35, 33, -1, 33, 35, 23, 6, -1, 6, 23, 21, 5, -1, 5, 21, 25, 7, 31, 25, 35, 23, 21, -1, 26, 22, 24, 36, -1, 25, 26, 36, 35, -1, 35, 36, 24, 23, -1, 23, 24, 22, 21, -1, 21, 22, 26, 25, 31, 9, 34, 33, 7, -1, 29, 25, 35, 37, -1, 9, 29, 37, 34, -1, 34, 37, 35, 33, -1, 33, 35, 25, 7, -1, 7, 25, 29, 9, 31, 29, 37, 35, 25, -1, 30, 26, 36, 38, -1, 29, 30, 38, 37, -1, 37, 38, 36, 35, -1, 35, 36, 26, 25, -1, 25, 26, 30, 29, 31, 0, 1, 8, 10, -1, 11, 31, 27, 13, -1, 0, 11, 13, 1, -1, 1, 13, 27, 8, -1, 8, 27, 31, 10, -1, 10, 31, 11, 0, 31, 11, 13, 27, 31, -1, 12, 32, 28, 14, -1, 11, 12, 14, 13, -1, 13, 14, 28, 27, -1, 27, 28, 32, 31, -1, 31, 32, 12, 11, 31, 6, 8, 1, 2, -1, 23, 15, 13, 27, -1, 6, 23, 27, 8, -1, 8, 27, 13, 1, -1, 1, 13, 15, 2, -1, 2, 15, 23, 6, 31, 23, 27, 13, 15, -1, 24, 16, 14, 28, -1, 23, 24, 28, 27, -1, 27, 28, 14, 13, -1, 13, 14, 16, 15, -1, 15, 16, 24, 23, 31, 6, 2, 3, 4, 5, -1, 23, 21, 19, 17, 15, -1, 2, 6, 23, 15, -1, 3, 2, 15, 17, -1, 4, 3, 17, 19, -1, 5, 4, 19, 21, -1, 6, 5, 21, 23, 31, 23, 15, 17, 19, 21, -1, 24, 22, 20, 18, 16, -1, 15, 23, 24, 16, -1, 17, 15, 16, 18, -1, 19, 17, 18, 20, -1, 21, 19, 20, 22, -1, 23, 21, 22, 24]
+        )
+        # fmt: on
+        cI = DataArrayInt([0, 30, 60, 90, 120, 150, 180, 210, 240, 277, 314])
+        m3 = MEDCouplingUMesh("box", 3)
+        m3.setCoords(coo)
+        m3.setConnectivity(c, cI)
+        m3.checkConsistency()
+        m2, _, _, _, _ = m3.buildDescendingConnectivity()
+        grpIds = DataArrayInt([4, 9, 35, 39])
+        grpIds.setName("group")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m3)
+        mfu.setMeshAtLevel(-1, m2)
+        mfu.setGroupsAtLevel(-1, [grpIds])
+        m2, _, _, _, _ = m3.buildDescendingConnectivity()
+        grpIds = DataArrayInt([4, 9, 35, 39])
+        grpIds.setName("group")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m3)
+        mfu.setMeshAtLevel(-1, m2)
+        mfu.setGroupsAtLevel(-1, [grpIds])
+        nNod = m3.getNumberOfNodes()
+
+        c2o2nN = mfu.crackAlong("group")
+
+        self.assertEqual({6: {6: 49, 23: 47}, 7: {23: 47, 24: 43}, 8: {6: 50, 23: 48, 21: 46, 5: 45, 2: 41, 15: 40}, 9: {23: 48, 21: 46, 24: 44, 22: 42, 15: 40, 16: 39}}, c2o2nN)
+
+        m3_bis = mfu.getMeshAtLevel(0)
+        m3_bis.checkConsistency()
+        m2_bis = mfu.getMeshAtLevel(-1)
+        m2_bis.checkConsistency()
+        self.assertEqual(nNod + 12, mfu.getNumberOfNodes())
+        self.assertEqual(nNod + 12, m3_bis.getNumberOfNodes())
+        self.assertEqual(nNod + 12, m2_bis.getNumberOfNodes())
+        # self.assertEqual([2, 5, 6, 15, 16, 21, 22, 23, 24], nodesDup.getValues())
+        # self.assertEqual(
+        #     m3_bis.getCoords()[nodesDup].getValues(),
+        #     m3_bis.getCoords()[nNod:].getValues(),
+        # )
+        # self.assertEqual(set([0, 1, 6, 7]), set(cells1.getValues()))
+        # self.assertEqual(set([8, 9]), set(cells2.getValues()))
+        self.assertEqual([4, 9, 35, 39, 49, 50, 51, 52], mfu.getGroupArr(-1, "group").getValues())
+        m_desc, _, _, _, _ = m3_bis.buildDescendingConnectivity()
+        m_desc.checkDeepEquivalOnSameNodesWith(m2_bis, 2, 9.9999)
+
+    def testBuildInnerBoundary10(self):
+        """2D tests where some cells are touching the M1 group with just a node, and are **not** neighbor
+        with any cells touching the M1 group by a face.
+        """
+        # fmt: off
+        coo = DataArrayDouble(
+            [0.0, 0.0, 1.0, 0.0, 2.0, 0.0, 3.0, 0.0, 0.0, 1.0, 1.0, 1.0, 2.0, 1.0, 3.0, 1.0, 0.0, 2.0, 1.0, 2.0, 2.0, 2.0, 3.0, 2.0, 2.5, 0.0, 3.0, 0.5, ],
+            14,
+            2,
+        )
+        c = DataArrayInt(
+            [3, 12, 2, 6, 3, 3, 12, 6, 3, 13, 3, 6, 3, 7, 13, 6, 4, 1, 0, 4, 5, 4, 2, 1, 5, 6, 4, 5, 4, 8, 9, 4, 6, 5, 9, 10, 4, 7, 6, 10, 11]
+        )
+        # fmt: on
+        cI = DataArrayInt([0, 4, 8, 12, 16, 21, 26, 31, 36, 41])
+        m2 = MEDCouplingUMesh("mesh", 2)
+        m2.setCoords(coo)
+        m2.setConnectivity(c, cI)
+        m2.checkConsistency()
+        m1, _, _, _, _ = m2.buildDescendingConnectivity()
+        grpIds = DataArrayInt([8, 14])
+        grpIds.setName("group")
+        mfu = MEDFileUMesh()
+        mfu.setMeshAtLevel(0, m2)
+        mfu.setMeshAtLevel(-1, m1)
+        mfu.setGroupsAtLevel(-1, [grpIds])
+        nNod = m2.getNumberOfNodes()
+
+        c2o2nN = mfu.crackAlong("group")
+
+        self.assertEqual({7: {6: 15}, 8: {6: 15, 7: 14}}, c2o2nN)
+
+        m2_bis = mfu.getMeshAtLevel(0)
+        m2_bis.checkConsistency()
+        m1_bis = mfu.getMeshAtLevel(-1)
+        m1_bis.checkConsistency()
+        self.assertEqual(nNod + 2, mfu.getNumberOfNodes())
+        self.assertEqual(nNod + 2, m2_bis.getNumberOfNodes())
+        self.assertEqual(nNod + 2, m1_bis.getNumberOfNodes())
+        # self.assertEqual([6, 7], nodesDup.getValues())
+        # self.assertEqual([2.0, 1.0, 3.0, 1.0], m2_bis.getCoords()[nNod:].getValues())
+        # self.assertEqual(set([0, 1, 2, 3, 5]), set(cells1.getValues()))
+        # self.assertEqual(set([7, 8]), set(cells2.getValues()))
+        self.assertEqual([8, 14, 22, 23], mfu.getGroupArr(-1, "group").getValues())
+
+
+if __name__ == "__main__":
+    unittest.main()

diff --git a/src/MEDLoader/Swig/MEDLoaderCommon.i b/src/MEDLoader/Swig/MEDLoaderCommon.i
index aebbe75c3522d93c3eb8cbc4a085cef3bb8d2de0..6f97b005fcffaa2485e4e50620e5fc5b650e92ee 100644 (file)
--- a/src/MEDLoader/Swig/MEDLoaderCommon.i
+++ b/src/MEDLoader/Swig/MEDLoaderCommon.i
@@ -1729,6 +1729,60 @@ namespace MEDCoupling
            return ret;
          }
 
+         PyObject * crackAlong(const std::string& grpNameM1, bool grpMustBeFullyDup=true)
+         {
+           const std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>> cellOld2NewNode = self->crackAlong(grpNameM1, grpMustBeFullyDup);
+           PyObject * cellOld2NewNodePy = PyDict_New();
+           for (const auto & cellIt: cellOld2NewNode) {
+               PyObject * old2NewPy = PyDict_New();
+               for (const auto & keyVal: cellIt.second) {
+                   PyDict_SetItem(old2NewPy, PyInt_FromLong(keyVal.first), PyInt_FromLong(keyVal.second));
+               }
+               PyDict_SetItem(cellOld2NewNodePy, PyInt_FromLong(cellIt.first), old2NewPy);
+           }
+           return cellOld2NewNodePy;
+         }
+
+         // void OpenCrack(const std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>>& c2o2nN, double factor=0.9)
+         // {
+         //   self->OpenCrack(c2o2nN, factor);
+         // }
+
+         void openCrack(PyObject* c2o2nNPy, double factor=0.9)
+         {
+             if (!PyDict_Check(c2o2nNPy))
+                 throw INTERP_KERNEL::Exception("c2o2nNPy is not a python dict");
+
+             std::unordered_map<mcIdType, std::unordered_map<mcIdType, mcIdType>> c2o2nN;
+             Py_ssize_t i_cells = 0;
+             PyObject *cell, *val;
+
+             while (PyDict_Next(c2o2nNPy, &i_cells, &cell, &val)) {
+                 if (!PyInt_Check(cell))
+                     throw INTERP_KERNEL::Exception("One of c2o2nNPy cell id is not an int.");
+
+                 if (!PyDict_Check(val))
+                     throw INTERP_KERNEL::Exception("One of c2o2nNPy val is not a python dict");
+                 std::unordered_map<mcIdType, mcIdType> o2nN {};
+                 Py_ssize_t i_node = 0;
+                 PyObject *oldN, *newN;
+
+                 while (PyDict_Next(val, &i_node, &oldN, &newN)) {
+                     if (!PyInt_Check(oldN))
+                         throw INTERP_KERNEL::Exception(
+                             "One of c2o2nNPy old node id is not an int.");
+                     if (!PyInt_Check(newN))
+                         throw INTERP_KERNEL::Exception(
+                             "One of c2o2nNPy new node id is not an int.");
+                     o2nN[PyInt_AsLong(oldN)] = PyInt_AsLong(newN);
+                 }
+
+                 c2o2nN[PyInt_AsLong(cell)] = o2nN;
+             }
+
+             self->openCrack(c2o2nN, factor);
+         }
+
          MEDCoupling1GTUMesh *getDirectUndergroundSingleGeoTypeMesh(INTERP_KERNEL::NormalizedCellType gt) const
          {
            MEDCoupling1GTUMesh *ret(self->getDirectUndergroundSingleGeoTypeMesh(gt));

diff --git a/src/MEDLoader/Swig/tests.set b/src/MEDLoader/Swig/tests.set
index 2e58e373954f9a2d2ec2ca56c608c2296e81ad5b..f57a3e06b4ff78c35b57f612607fbde73d88947a 100644 (file)
--- a/src/MEDLoader/Swig/tests.set
+++ b/src/MEDLoader/Swig/tests.set
@@ -23,6 +23,7 @@ SET(BASE_TESTS
   MEDLoaderExamplesTest.py
   UsersGuideExamplesTest_ML.py
   MeshFormatWriterTest1.py
+  CrackAlongTest.py
 )
 
 # if numpy is used

diff --git a/src/MEDLoader/Test/CMakeLists.txt b/src/MEDLoader/Test/CMakeLists.txt
index b1b1b907e85b88b400dee006e1388bc568d204e6..274ce2fec82a4f71240d35cb278e7376ab135da4 100644 (file)
--- a/src/MEDLoader/Test/CMakeLists.txt
+++ b/src/MEDLoader/Test/CMakeLists.txt
@@ -31,6 +31,7 @@ INCLUDE_DIRECTORIES(
 
 SET(TestMEDLoader_SOURCES
   TestMEDLoader.cxx
+  CrackAlgoTest.cxx
   MEDLoaderTest.cxx
   )
 

diff --git a/src/MEDLoader/Test/CrackAlgoTest.cxx b/src/MEDLoader/Test/CrackAlgoTest.cxx
new file mode 100644 (file)
index 0000000..7fdb906
--- /dev/null
+++ b/src/MEDLoader/Test/CrackAlgoTest.cxx
@@ -0,0 +1,657 @@
+// Copyright (C) 2007-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author : Aymeric SONOLET (CEA/DES)
+
+#include "CrackAlgoTest.hxx"
+
+#include <string>
+#include <algorithm>
+#include <vector>
+#include <utility>
+
+#include <InterpKernelException.hxx>
+#include <CrackAlgo.hxx>
+#include <MEDFileMesh.hxx>
+#include <MEDCouplingMemArray.hxx>
+#include <MEDCouplingCMesh.hxx>
+#include <MEDCouplingUMesh.hxx>
+#include <MCAuto.hxx>
+#include <MCType.hxx>
+
+using std::vector;
+using std::string;
+using std::pair;
+using std::copy;
+
+using MEDCoupling::CrackAlgoTest;
+using MEDCoupling::MEDFileUMesh;
+using MEDCoupling::MEDCouplingCMesh;
+using MEDCoupling::MEDCouplingUMesh;
+using MEDCoupling::DataArrayIdType;
+using MEDCoupling::MCAuto;
+
+using MFU = MCAuto<MEDFileUMesh>;
+using MCC = MCAuto<MEDCouplingCMesh>;
+using MCU = MCAuto<MEDCouplingUMesh>;
+using DAI = MCAuto<DataArrayIdType>;
+
+
+/*
+ *            TOP z=2              MIDDLE z=1           BOTTOM z=0
+ *   2 +--------+--------+   +--------+--------+   +--------+--------+
+ *     |        |        |   |        |        |   |        o        |
+ *     |        |        |   |        |        |   |        o        |
+ *     |        |        |   |        |        |   |       +-+       |
+ *   1 +--------+--------+   +--------+--------+   +ooooooo|o|ooooooo+
+ *     |        |        |   |        |        |   |       +-+       |
+ *     |        |        |   |        |        |   |        o        |
+ *     |        |        |   |        |        |   |        o        |
+ *   0 +--------+--------+   +--------+--------+   +--------+--------+
+ *
+ *     0        1        2
+ */
+void CrackAlgoTest::test3DHalfCrossCut(
+) {
+    MFU mf = make2x2Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box[6] {1.0, 1.0, 1.0, 1.0, 0.1, 0.9};
+    MCAuto<DataArrayIdType> g1 = m1->getCellsInBoundingBox(box, 0.01);
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 4);
+    const auto res = TestCrack(mf, "M1", "halfCrossCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+/*
+ *            TOP z=2              MIDDLE z=1           BOTTOM z=0
+ *   2 +--------+--------+   ooooooooooooooooooo   +--------+--------+
+ *     |        o        |   ooooooooooooooooooo   |        o        |
+ *     |        o        |   ooooooooooooooooooo   |        o        |
+ *     |       +-+       |   oooooooo+-+oooooooo   |       +-+       |
+ *   1 +ooooooo|o|ooooooo+   oooooooo|o|oooooooo   +ooooooo|o|ooooooo+
+ *     |       +-+       |   oooooooo+-+oooooooo   |       +-+       |
+ *     |        o        |   ooooooooooooooooooo   |        o        |
+ *     |        o        |   ooooooooooooooooooo   |        o        |
+ *   0 +--------+--------+   ooooooooooooooooooo   +--------+--------+
+ *
+ *     0        1        2
+ */
+void CrackAlgoTest::test3DFullFullCut(
+) {
+    MFU mf = make2x2Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box[6] {1.0, 1.0, 1.0, 1.0, 0.1, 1.9};
+    DAI g1 = m1->getCellsInBoundingBox(box, 0.01);
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    const auto res = TestCrack(mf, "M1", "fullFullCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 12);
+}
+
+/*
+ *            TOP z=2              MIDDLE z=1           BOTTOM z=0
+ *   2 +--------+--------+   +--------+--------+   +--------+--------+
+ *     |        o        |   |        o        |   |        o        |
+ *     |        o        |   |        o        |   |        o        |
+ *     |       +-+       |   |        o        |   |       +-+       |
+ *   1 +ooooooo|o|ooooooo+   +ooooooooooooooooo+   +ooooooo|o|ooooooo+
+ *     |       +-+       |   |        o        |   |       +-+       |
+ *     |        o        |   |        o        |   |        o        |
+ *     |        o        |   |        o        |   |        o        |
+ *   0 +--------+--------+   +--------+--------+   +--------+--------+
+ *
+ *     0        1        2
+ */
+void CrackAlgoTest::test3DFullCrossCut(
+) {
+    MFU mf = make2x2Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box[6] {1.0, 1.0, 1.0, 1.0, 0.1, 0.9};
+    MCAuto<DataArrayIdType> g1 = m1->getCellsInBoundingBox(box, 0.01);
+    double box2[6] {1.0, 1.0, 1.0, 1.0, 1.1, 1.9};
+    MCAuto<DataArrayIdType> g2 = m1->getCellsInBoundingBox(box, 0.01);
+    g1->aggregate(g2);
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 8);
+
+    const auto res = TestCrack(mf, "M1", "fullCrossCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+/*
+ *            TOP z=2              MIDDLE z=1           BOTTOM z=0
+ *   2 +--------+--------+   +--------+--------+   +--------+--------+
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |  +-+   |   +-+  |
+ *   1 +--------+--------+   +--------+--------+   +oo|o|ooooooo|o|oo+
+ *     |        |        |   |        |        |   |  +-+   |   +-+  |
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *   0 +--------+--------+   +--------+--------+   +--------+--------+
+ *
+ *     0        1        2
+ */
+void CrackAlgoTest::test3DHalfCut(
+) {
+    MFU mf = make2x2Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box1[6] {0.5, 0.5, 1.0, 1.0, 0.1, 0.9};
+    DAI g1 = m1->getCellsInBoundingBox(box1, 0.01);
+    double box2[6] {1.5, 1.5, 1.0, 1.0, 0.1, 0.9};
+    DAI g2 = m1->getCellsInBoundingBox(box2, 0.01);
+
+    g1->aggregate(g2);
+
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 2);
+    const auto res = TestCrack(mf, "M1", "halfCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+/*
+ *            TOP z=2              MIDDLE z=1           BOTTOM z=0
+ *   2 +--------+--------+   +--------+--------+   +--------+--------+
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *     |  +-+   |   +-+  |   |        |        |   |  +-+   |   +-+  |
+ *   1 +oo|o|ooooooo|o|oo+   +ooooooooooooooooo+   +oo|o|ooooooo|o|oo+
+ *     |  +-+   |   +-+  |   |        |        |   |  +-+   |   +-+  |
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *   0 +--------+--------+   +--------+--------+   +--------+--------+
+ *
+ *     0        1        2
+ */
+void CrackAlgoTest::test3DFullCut(
+) {
+    MFU mf = make2x2Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box1[6] {0.5, 0.5, 1.0, 1.0, 0.1, 0.9};
+    DAI g1 = m1->getCellsInBoundingBox(box1, 0.01);
+    double box2[6] {1.5, 1.5, 1.0, 1.0, 0.1, 0.9};
+    DAI g2 = m1->getCellsInBoundingBox(box2, 0.01);
+    double box3[6] {0.5, 0.5, 1.0, 1.0, 1.1, 1.9};
+    DAI g3 = m1->getCellsInBoundingBox(box1, 0.01);
+    double box4[6] {1.5, 1.5, 1.0, 1.0, 1.1, 1.9};
+    DAI g4 = m1->getCellsInBoundingBox(box2, 0.01);
+
+    g1->aggregate(g2);
+    g1->aggregate(g3);
+    g1->aggregate(g4);
+
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 4);
+    const auto res = TestCrack(mf, "M1", "fullCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+/*
+ *            TOP z=2              MIDDLE z=1           BOTTOM z=0
+ *   2 +--------+--------+   +--------+--------+   +--------+--------+
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *     | +---+  |        |   |        |        |   |        |        |
+ *   1 oo|ooo|oo+--------+   +--------+--------+   +--------+--------+
+ *     | +---+  |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *     |        |        |   |        |        |   |        |        |
+ *   0 +--------+--------+   +--------+--------+   +--------+--------+
+ *
+ *     0        1        2
+ */
+void CrackAlgoTest::test3DAngleCut(
+) {
+    MFU mf = make2x2Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box[6] {0.1, 0.9, 1.0, 1.0, 1.1, 1.9};
+    MCAuto<DataArrayIdType> g1 = m1->getCellsInBoundingBox(box, 0.01);
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 1);
+    const auto res = TestCrack(mf, "M1", "angleCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+void CrackAlgoTest::testInnerCrossCut(
+) {
+    MFU mf = make4x4Voxel();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box[6] {2.0, 2.0, 2.0, 2.0, 2.0, 2.0};
+    MCAuto<DataArrayIdType> g1 = m1->getCellsInBoundingBox(box, 0.01);
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 12);
+    const auto res = TestCrack(mf, "M1", "innerCrossCut");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+void CrackAlgoTest::test2DGrp(
+) {
+    MCAuto<MEDFileUMesh> mf = make2DMesh();
+    const auto res = TestCrack(mf, "Grp", "2dGrp");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+void CrackAlgoTest::test1DMesh(
+) {
+    MFU mf = make1DMesh();
+
+    const MCU m1 = mf->getMeshAtLevel(-1);
+    double box[2] {1.9, 2.1};
+    MCAuto<DataArrayIdType> g1 = m1->getCellsInBoundingBox(box, 0.01);
+    g1->setName("M1");
+    mf->setGroupsAtLevel(-1, vector<const DataArrayIdType*>{g1});
+
+    CPPUNIT_ASSERT(g1->getNumberOfTuples() == 1);
+    const auto res = TestCrack(mf, "M1", "1DMesh");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+void CrackAlgoTest::test2DMeshNonConnexCut(
+) {
+    MCAuto<MEDFileUMesh> mf = make2DMesh2();
+    const auto res = TestCrack(mf, "Grp", "2dGrpNonConnex");
+    CPPUNIT_ASSERT(res.first);
+    CPPUNIT_ASSERT(res.second);
+}
+
+MEDFileUMesh * CrackAlgoTest::make2x2Voxel(
+) {
+    double coords[3]  {0., 1., 2.};
+    MCAuto<DataArrayDouble> x_coords = DataArrayDouble::New();
+    x_coords->alloc(3);
+    copy(coords, coords + 3, x_coords->getPointer());
+
+    MCC m_cmesh = MEDCouplingCMesh::New();
+    m_cmesh->setCoords(x_coords, x_coords, x_coords);
+    MCU m0 = m_cmesh->buildUnstructured();
+    m0->setName("UMesh");
+
+    DAI desc(DataArrayIdType::New()),
+        descIdx(DataArrayIdType::New()),
+        revDesc(DataArrayIdType::New()),
+        revDescIdx(DataArrayIdType::New());
+    MCU m1 = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+
+    MFU mf = MEDFileUMesh::New();
+    mf->setCoords(m0->getCoords());
+    mf->setMeshAtLevel(0, m0);
+    mf->setMeshAtLevel(-1, m1);
+
+    return mf.retn();
+}
+
+MEDFileUMesh * CrackAlgoTest::make1DMesh(
+) {
+    double coords[5]  {0., 1., 2., 3., 4.};
+    MCAuto<DataArrayDouble> x_coords = DataArrayDouble::New();
+    x_coords->alloc(5);
+    copy(coords, coords + 5, x_coords->getPointer());
+
+    MCC m_cmesh = MEDCouplingCMesh::New();
+    m_cmesh->setCoords(x_coords);
+    MCU m0 = m_cmesh->buildUnstructured();
+    m0->setName("UMesh");
+
+    DAI desc(DataArrayIdType::New()),
+        descIdx(DataArrayIdType::New()),
+        revDesc(DataArrayIdType::New()),
+        revDescIdx(DataArrayIdType::New());
+    MCU m1 = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+
+    MFU mf = MEDFileUMesh::New();
+    mf->setCoords(m0->getCoords());
+    mf->setMeshAtLevel(0, m0);
+    mf->setMeshAtLevel(-1, m1);
+
+    return mf.retn();
+}
+
+MEDFileUMesh * CrackAlgoTest::make4x4Voxel(
+) {
+    double coords[5]  {0., 1., 2., 3., 4.};
+    MCAuto<DataArrayDouble> x_coords = DataArrayDouble::New();
+    x_coords->alloc(5);
+    copy(coords, coords + 5, x_coords->getPointer());
+
+    MCAuto<MEDCouplingCMesh> m_cmesh = MEDCouplingCMesh::New();
+    m_cmesh->setCoords(x_coords, x_coords, x_coords);
+    MCAuto<MEDCouplingUMesh> m0 = m_cmesh->buildUnstructured();
+    m0->setName("UMesh");
+
+    DAI desc(DataArrayIdType::New()), descIdx(DataArrayIdType::New()), revDesc(DataArrayIdType::New()), revDescIdx(DataArrayIdType::New());
+    MCU m1 = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+
+    MEDFileUMesh * mf = MEDFileUMesh::New();
+    mf->setCoords(m0->getCoords());
+    mf->setMeshAtLevel(0, m0);
+    mf->setMeshAtLevel(-1, m1);
+
+    return mf;
+}
+
+MEDFileUMesh * CrackAlgoTest::make2DMesh(
+) {
+    using DAD = MCAuto<DataArrayDouble>;
+
+    double c0[6] {0.0, 1.1, 2.3, 3.6, 5.0, 6.5};
+    DAD coords0 = DataArrayDouble::New();
+    coords0->alloc(6);
+    copy(c0, c0 + 6, coords0->rwBegin());
+
+    double c1[5] {0.0, 1.1, 2.3, 3.6, 5.0};
+    DAD coords1 = DataArrayDouble::New();
+    coords1->alloc(5);
+    copy(c1, c1+ 5, coords1->rwBegin());
+
+    MCC m = MEDCouplingCMesh::New();
+    m->setCoordsAt(0, coords0);
+    m->setCoordsAt(1, coords1);
+
+    MCU m0 = m->buildUnstructured();
+    m0->setName("duplicate");
+
+    DAI desc(DataArrayIdType::New()),
+        descIdx(DataArrayIdType::New()),
+        revDesc(DataArrayIdType::New()),
+        revDescIdx(DataArrayIdType::New());
+    MCU m2 = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+    m2->setName(m0->getName());
+
+    mcIdType ids_[17] {8,11,14,20,21,22,23,24,25,26,31,32,33,34,35,36,37};
+    MCU m3 = m2->buildPartOfMySelf(ids_, ids_+17);
+
+    mcIdType grp_[3] {4, 6, 8};
+    DAI grp = DataArrayIdType::New();
+    grp->alloc(3);
+    copy(grp_, grp_+3, grp->rwBegin());
+    grp->setName("Grp");
+
+    mcIdType grp2_[2] {9, 16};
+    DAI grp2 = DataArrayIdType::New();
+    grp2->alloc(2);
+    copy(grp2_, grp2_+2, grp2->rwBegin());
+    grp2->setName("Grp2");
+
+    MEDFileUMesh * mm = MEDFileUMesh::New();
+    mm->setMeshAtLevel(0, m0);
+    mm->setMeshAtLevel(-1, m3);
+    mm->setGroupsAtLevel(-1, {grp, grp2});
+
+    mcIdType grpNode_[3] {4, 21, 23};
+    DAI grpNode = DataArrayIdType::New();
+    grpNode->alloc(3);
+    copy(grpNode_, grpNode_+3, grpNode->rwBegin());
+    grpNode->setName("GrpNode");
+    mm->setGroupsAtLevel(1, {grpNode});
+
+    return mm;
+}
+
+MEDFileUMesh * CrackAlgoTest::make2DMesh2(
+) {
+    using DAD = MCAuto<DataArrayDouble>;
+
+    double c0[5] {0.0, 1.1, 2.3, 3.6, 5.0};
+    DAD coords0 = DataArrayDouble::New();
+    coords0->alloc(5);
+    copy(c0, c0 + 5, coords0->rwBegin());
+
+    double c1[3] {0.0, 1.0, 2.0};
+    DAD coords1 = DataArrayDouble::New();
+    coords1->alloc(3);
+    copy(c1, c1 + 3, coords1->rwBegin());
+
+    MCC m = MEDCouplingCMesh::New();
+    m->setCoordsAt(0, coords0);
+    m->setCoordsAt(1, coords1);
+
+    MCU m0 = m->buildUnstructured();
+    m0->setName("simple");
+
+    DAI desc(DataArrayIdType::New()),
+        descIdx(DataArrayIdType::New()),
+        revDesc(DataArrayIdType::New()),
+        revDescIdx(DataArrayIdType::New());
+    MCU m2 = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+    m2->setName(m0->getName());
+
+    mcIdType grp_[2] {3, 19};
+    DAI grp = DataArrayIdType::New();
+    grp->alloc(2);
+    copy(grp_, grp_+2, grp->rwBegin());
+    grp->setName("Grp");
+
+    MEDFileUMesh * mm = MEDFileUMesh::New();
+    mm->setMeshAtLevel(0, m0);
+    mm->setMeshAtLevel(-1, m2);
+    mm->setGroupsAtLevel(-1, {grp});
+
+    return mm;
+}
+
+
+pair<bool, bool> CrackAlgoTest::TestCrack(
+    MEDFileUMesh * mm_init,
+    const string & grp_name,
+    const string & test_name
+) {
+    const auto c2cBroken = GetC2CBroken(mm_init, grp_name);
+    const auto c2cPreserved = GetC2CPreserved(mm_init, grp_name);
+
+    const MCU f2dup = mm_init->getGroup(-1, grp_name);
+
+    mm_init->write(test_name + "_in.med", 2);
+    const auto cellOld2NewNode = CrackAlgo::Compute(mm_init, grp_name);
+    mm_init->write(test_name + "_out.med", 2);
+
+    const MCU f2dup_b = mm_init->getGroup(-1, grp_name);
+
+    const auto res = pair<bool, bool> {
+        CheckM0Mesh(mm_init, c2cBroken, c2cPreserved),
+        CheckM1Mesh(f2dup, f2dup_b)
+        };
+
+    CrackAlgo::OpenCrack(mm_init, cellOld2NewNode);
+    mm_init->write(test_name + "_cracked.med", 2);
+
+    return res;
+}
+
+bool CrackAlgoTest::CheckM1Mesh(
+    const MEDCouplingUMesh * f2dup_before,
+    const MEDCouplingUMesh * f2dup_after
+) {
+    const mcIdType nbFaces_0 = f2dup_before->getNumberOfCells();
+    const mcIdType nbFaces_1 = f2dup_after->getNumberOfCells();
+
+    MCU f2dup_before_copy = f2dup_after->deepCopy();
+    f2dup_before_copy->zipCoords();
+
+    MCU f2dup_after_copy = f2dup_after->deepCopy();
+    bool nodesAreMerged = false;
+    mcIdType newNbOfNodes = 0;
+    DAI o2n = f2dup_after_copy->mergeNodes(1e-12, nodesAreMerged, newNbOfNodes);
+    f2dup_after_copy->zipCoords();
+
+    DataArrayIdType * cellCor;
+    DataArrayIdType * nodeCor;
+    f2dup_after_copy->checkGeoEquivalWith(f2dup_before_copy, 12, 1e-12, cellCor, nodeCor);
+    DAI cellCorrI(cellCor);
+    DAI nodeCorrI(nodeCor);
+
+    return (2 * nbFaces_0 == nbFaces_1);
+}
+
+vector<pair<mcIdType, mcIdType>>
+CrackAlgoTest::GetC2CBroken(
+    const MEDFileUMesh * mm,
+    const string & grp_name
+) {
+    vector<pair<mcIdType, mcIdType>> res;
+    MCU m0 = mm->getMeshAtLevel(0);
+    DAI desc(DataArrayIdType::New()),
+        descIdx(DataArrayIdType::New()),
+        revDesc(DataArrayIdType::New()),
+        revDescIdx(DataArrayIdType::New());
+    MCU mf = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+    MCU f2dupMesh = mm->getGroup(-1, grp_name);
+    DataArrayIdType * f2dup_p;
+    mf->areCellsIncludedIn(f2dupMesh, 2, f2dup_p);
+    DAI f2dup(f2dup_p);
+
+    for (const auto & f : *f2dup) {
+        const auto & cellsId_start = revDescIdx->begin()[f];
+        const auto & cellsId_end = revDescIdx->begin()[f+1];
+        if (cellsId_end - cellsId_start != 2)
+            throw INTERP_KERNEL::Exception(
+                "crackAlong: A face of group M1 does not have two neighbors.");
+        const auto & cell0 = revDesc->begin()[cellsId_start];
+        const auto & cell1 = revDesc->begin()[cellsId_start + 1];
+        res.push_back(pair<mcIdType, mcIdType>(cell0, cell1));
+    }
+    return res;
+}
+
+vector<pair<mcIdType, mcIdType>>
+CrackAlgoTest::GetC2CPreserved(
+    const MEDFileUMesh * mm,
+    const string & grp_name
+) {
+    vector<pair<mcIdType, mcIdType>> res;
+    MCU m0 = mm->getMeshAtLevel(0);
+    DAI desc(DataArrayIdType::New()),
+        descIdx(DataArrayIdType::New()),
+        revDesc(DataArrayIdType::New()),
+        revDescIdx(DataArrayIdType::New());
+    MCU mf = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+    MCU f2dupMesh = mm->getGroup(-1, grp_name);
+
+    DataArrayIdType * f2dup_p;
+    mf->areCellsIncludedIn(f2dupMesh, 2, f2dup_p);
+    DAI f2dup(f2dup_p);
+
+    mcIdType nbOfFaces = mf->getNumberOfCells();
+    for (mcIdType f = 0; f < nbOfFaces; f++) {
+        if (std::find(f2dup->begin(), f2dup->end(), f) != f2dup->end())
+            continue;
+
+        const auto & cellsId_start = revDescIdx->begin()[f];
+        const auto & cellsId_end = revDescIdx->begin()[f+1];
+        if (cellsId_end - cellsId_start == 1)
+            continue;
+        if (cellsId_end - cellsId_start > 2)
+            throw INTERP_KERNEL::Exception(
+                "crackAlong: A face of the mf has more than two neighbors.");
+
+        const auto & cell0 = revDesc->begin()[cellsId_start];
+        const auto & cell1 = revDesc->begin()[cellsId_start + 1];
+        res.push_back(pair<mcIdType, mcIdType>(cell0, cell1));
+    }
+    return res;
+}
+
+bool
+CrackAlgoTest::CheckM0Mesh(
+    const MEDFileUMesh * mm,
+    const Connections& c2cBrokenConnection,
+    const Connections& c2cPreservedConnection
+) {
+    MCU m0 = mm->getMeshAtLevel(0);
+    DAI desc(DataArrayIdType::New()), descIdx(DataArrayIdType::New()), revDesc(DataArrayIdType::New()), revDescIdx(DataArrayIdType::New());
+    MCU mf = m0->buildDescendingConnectivity(desc, descIdx, revDesc, revDescIdx);
+
+    CrackAlgo::Map2Set c2c;
+
+    const mcIdType * descIdx_ptr = descIdx->begin();
+    const mcIdType * desc_ptr = desc->begin();
+    const mcIdType * revDescIdx_ptr = revDescIdx->begin();
+    const mcIdType * revDesc_ptr = revDesc->begin();
+    for (mcIdType cell = 0; cell < m0->getNumberOfCells(); cell++) {
+        auto & neighbors = c2c[cell];
+        for (
+            mcIdType descId = descIdx_ptr[cell];
+            descId < descIdx_ptr[cell+1];
+            descId++
+        ) {
+            const auto & face = desc_ptr[descId];
+            for (
+                mcIdType revDescId = revDescIdx_ptr[face];
+                revDescId < revDescIdx_ptr[face+1];
+                revDescId++
+            ) {
+                const auto & otherCell = revDesc_ptr[revDescId];
+                if (otherCell != cell)
+                    neighbors.insert(otherCell);
+            }
+        }
+    }
+
+    bool res = true;
+    for (const auto & conn : c2cBrokenConnection) {
+        const auto & cell0 = conn.first;
+        const auto & cell1 = conn.second;
+        // NOTE: all cells are in c2c so cell0 and cell1 are in c2c
+        if (c2c.at(cell0).find(cell1) != c2c.at(cell0).end())
+            res = false;
+        if (c2c.at(cell1).find(cell0) != c2c.at(cell1).end())
+            res = false;
+    }
+
+    for (const auto & conn : c2cPreservedConnection) {
+        const auto & cell0 = conn.first;
+        const auto & cell1 = conn.second;
+        // NOTE: all cells are in c2c so cell0 and cell1 are in c2c
+        if (c2c.at(cell0).find(cell1) == c2c.at(cell0).end())
+            res = false;
+        if (c2c.at(cell1).find(cell0) == c2c.at(cell1).end())
+            res = false;
+    }
+    return res;
+}

diff --git a/src/MEDLoader/Test/CrackAlgoTest.hxx b/src/MEDLoader/Test/CrackAlgoTest.hxx
new file mode 100644 (file)
index 0000000..549942d
--- /dev/null
+++ b/src/MEDLoader/Test/CrackAlgoTest.hxx
@@ -0,0 +1,108 @@
+// Copyright (C) 2007-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author : Aymeric SONOLET (CEA/DES)
+
+#ifndef __CRACKALGOTEST_HXX__
+#define __CRACKALGOTEST_HXX__
+
+#include <cppunit/extensions/HelperMacros.h>
+#include <string>
+#include <vector>
+#include <utility>
+
+#include <MCIdType.hxx>
+
+namespace MEDCoupling {
+
+class MEDFileUMesh;
+class MEDCouplingUMesh;
+
+class CrackAlgoTest : public CppUnit::TestFixture {
+    CPPUNIT_TEST_SUITE(CrackAlgoTest);
+    CPPUNIT_TEST(test3DHalfCrossCut);
+    CPPUNIT_TEST(test3DFullCrossCut);
+    CPPUNIT_TEST(test3DHalfCut);
+    CPPUNIT_TEST(test3DFullCut);
+    CPPUNIT_TEST(test3DFullFullCut);
+    CPPUNIT_TEST(test3DAngleCut);
+    CPPUNIT_TEST(test2DGrp);
+    CPPUNIT_TEST(test2DMeshNonConnexCut);
+    CPPUNIT_TEST(test1DMesh);
+    CPPUNIT_TEST(testInnerCrossCut);
+    CPPUNIT_TEST_SUITE_END();
+
+ public:
+
+    void test3DHalfCrossCut();
+    void test3DFullCrossCut();
+    void test3DAngleCut();
+    void test3DFullCut();
+    void test3DFullFullCut();
+    void test3DHalfCut();
+    void test2DGrp();
+    void test2DMeshNonConnexCut();
+    void test1DMesh();
+    void testInnerCrossCut();
+
+ private:
+    using Connections = std::vector<std::pair<mcIdType, mcIdType>>;
+    static std::pair<bool, bool>
+    TestCrack(
+        MEDFileUMesh * mm_init,
+        const std::string & grp_name,
+        const std::string & test_name);
+
+    static bool
+    CheckM0Mesh(
+        const MEDFileUMesh * mm,
+        const Connections& c2cBrokenConnection,
+        const Connections& c2cPreservedConnection);
+
+    static bool
+    CheckM1Mesh(
+        const MEDCouplingUMesh * f2dup_before,
+        const MEDCouplingUMesh * f2dup_after);
+
+    static Connections
+    GetC2CBroken(
+        const MEDFileUMesh * mm,
+        const std::string & grp_name);
+
+    static std::vector<std::pair<mcIdType, mcIdType>>
+    GetC2CPreserved(
+        const MEDFileUMesh * mm,
+        const std::string & grp_name);
+
+    static MEDFileUMesh *
+    make2x2Voxel();
+
+    static MEDFileUMesh *
+    make4x4Voxel();
+
+    static MEDFileUMesh *
+    make2DMesh();
+
+    static MEDFileUMesh *
+    make2DMesh2();
+
+    static MEDFileUMesh *
+    make1DMesh();
+};
+}  // namespace MEDCoupling
+#endif  // __CRACKALGOTEST_HXX__

diff --git a/src/MEDLoader/Test/MEDLoaderTest.cxx b/src/MEDLoader/Test/MEDLoaderTest.cxx
index 0d1add2c1b986b70f3c8fa7e2677a419947e9fdc..5d1c77ddb147dac736f2b9ae846534e3da945200 100644 (file)
--- a/src/MEDLoader/Test/MEDLoaderTest.cxx
+++ b/src/MEDLoader/Test/MEDLoaderTest.cxx
@@ -19,6 +19,7 @@
 // Author : Anthony Geay (CEA/DEN)
 
 #include "MEDLoaderTest.hxx"
+#include "MEDCouplingCMesh.hxx"
 #include "MEDLoader.hxx"
 #include "MEDLoaderBase.hxx"
 #include "MEDCouplingUMesh.hxx"
@@ -27,8 +28,9 @@
 #include "MEDCouplingFieldInt64.hxx"
 #include "MEDCouplingMemArray.hxx"
 #include "TestInterpKernelUtils.hxx"  // getResourceFile()
+#include "MEDFileMesh.hxx"
 
-#include <cmath>
+#include <algorithm>
 #include <numeric>
 
 using namespace MEDCoupling;
@@ -1569,6 +1571,3 @@ MEDCouplingFieldDouble *MEDLoaderTest::buildVecFieldOnGaussNE_1()
   m->decrRef();
   return f;
 }
-
-
-

diff --git a/src/MEDLoader/Test/MEDLoaderTest.hxx b/src/MEDLoader/Test/MEDLoaderTest.hxx
index baccfc28257c4d6a1f2ae32f6ad7d83c4f5deb08..dc3332af2410d5e33c2d8bf144a97821ab4c210c 100644 (file)
--- a/src/MEDLoader/Test/MEDLoaderTest.hxx
+++ b/src/MEDLoader/Test/MEDLoaderTest.hxx
@@ -61,6 +61,7 @@ namespace MEDCoupling
     CPPUNIT_TEST(testMEDLoaderPolygonRead);
     CPPUNIT_TEST(testMEDLoaderPolyhedronRead);
 
+
     CPPUNIT_TEST_SUITE_END();
   public:
     void testMesh1DRW();
@@ -89,6 +90,7 @@ namespace MEDCoupling
     void testMEDLoaderRead1();
     void testMEDLoaderPolygonRead();
     void testMEDLoaderPolyhedronRead();
+
   private:
     MEDCouplingUMesh *build1DMesh_1();
     MEDCouplingUMesh *build2DCurveMesh_1();

diff --git a/src/MEDLoader/Test/TestMEDLoader.cxx b/src/MEDLoader/Test/TestMEDLoader.cxx
index 68387f1f22be6f5a2e54bdb7c4f6dbeb2b52e7b1..9e73982d18314e7adda4c61fc72f3bd8533d8cc2 100644 (file)
--- a/src/MEDLoader/Test/TestMEDLoader.cxx
+++ b/src/MEDLoader/Test/TestMEDLoader.cxx
@@ -20,7 +20,9 @@
 
 #include "CppUnitTest.hxx"
 #include "MEDLoaderTest.hxx"
+#include "CrackAlgoTest.hxx"
 
 CPPUNIT_TEST_SUITE_REGISTRATION( MEDCoupling::MEDLoaderTest );
+CPPUNIT_TEST_SUITE_REGISTRATION( MEDCoupling::CrackAlgoTest );
 
 #include "BasicMainTest.hxx"