Salome HOME
Add test for .mesh file format
[tools/medcoupling.git] / src / MEDCoupling_Swig / MEDCouplingIntersectTest.py
index e02005e82d86659cdc3679c46899c4e9a24291e0..a67c8be966c75262dd78fad3053ab83e8bb2bddd 100644 (file)
@@ -1,5 +1,5 @@
 #  -*- coding: utf-8 -*-
-# Copyright (C) 2007-2019  CEA/DEN, EDF R&D
+# 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
 
 
 import sys
-if sys.platform == "win32":
-    from MEDCouplingCompat import *
-else:
-    from medcoupling import *
+from medcoupling import *
 import unittest
 from math import pi,e,sqrt,cos,sin
 from datetime import datetime
@@ -439,6 +436,32 @@ class MEDCouplingIntersectTest(unittest.TestCase):
         self.assertEqual(e2, res2Tool.getValues())
         pass
 
+    def testIntersect2DMeshes11(self):
+        """ Dealing properly with respective polygon orientation in QuadraticPolygon::haveIAChanceToBeCompletedBy()
+        The two polygons below have same orientation, but one edge of pol1 is colinear to pol2 in opposite directions.
+        """
+        eps = 1.0e-8
+        back = MEDCouplingUMesh('lback', 2)
+        coo = DataArrayDouble([(-2.5,-2.5),(-2.5,2.5),(2.5,2.5),(2.5,-2.5),(0,0),(0,1.66667),(1.66667,1.66667),(1.66667,0)])
+        back.setCoords(coo)
+        c = DataArrayInt([5, 6, 7, 4, 5, 2, 3, 7, 6, 5, 3, 0, 1, 2, 6, 4, 7])
+        cI = DataArrayInt([0, 4, 9, 17])
+        back.setConnectivity(c, cI)
+
+        tool = MEDCouplingUMesh('ltool', 2)
+        coo = DataArrayDouble([(0,0),(0,2.5),(2.5,2.5),(2.5,0)])
+        tool.setCoords(coo)
+        c = DataArrayInt([5, 0, 1, 2, 3])
+        cI = DataArrayInt([0, 5])
+        tool.setConnectivity(c, cI)
+
+        result, res2Back, res2Tool = MEDCouplingUMesh.Intersect2DMeshes(back, tool, eps)
+        self.assertEqual(result.getNodalConnectivity().getValues(), [5, 7, 8, 6, 5, 7, 6, 10, 11, 5, 11, 3, 7, 5, 9, 10, 6, 8, 5, 8, 7, 3, 0, 1, 9])
+        self.assertEqual(result.getNodalConnectivityIndex().getValues(), [0, 4, 9, 13, 18, 25])
+        self.assertEqual(res2Back.getValues(), [0, 1, 1, 2, 2])
+        self.assertEqual(res2Tool.getValues(), [0, 0, -1, 0, -1])
+        pass
+
     def testSwig2Intersect2DMeshesQuadra1(self):
         import cmath
         def createDiagCircle(lX, lY, R, cells=[0,1]):
@@ -671,6 +694,43 @@ class MEDCouplingIntersectTest(unittest.TestCase):
         self.assertEqual(res2Tool.getValues(), [0, -1, -1])
         pass
 
+    def testIntersect2DMeshesTmp13(self):
+        """ Bug fix: when part of mesh2 is fully included in mesh1 and some points of mesh2 are touching edges of mesh1, reconstruction of
+        residual cells was not properly done.
+        """
+        eps = 1.0e-6
+        # First case: only two points of mesh2 touching edges of mesh1
+        coo1 = DataArrayDouble([(0,0) , (1,1),  (1,0)])
+        mesh1 = MEDCouplingUMesh("mesh1", 2)
+        mesh1.setCoords(coo1)
+        c = DataArrayInt([NORM_TRI3, 0,1,2])
+        cI = DataArrayInt([0,4])
+        mesh1.setConnectivity(c, cI)
+
+        coo2 = DataArrayDouble([(0.5,0) , (0.5,0.5),  (0.7,0.25)])
+        mesh2 = MEDCouplingUMesh("mesh1", 2)
+        mesh2.setCoords(coo2)
+        c = DataArrayInt([NORM_TRI3, 0,1,2])
+        cI = DataArrayInt([0,4])
+        mesh2.setConnectivity(c, cI)
+
+        meshinter, res2Back, res2Tool = MEDCouplingUMesh.Intersect2DMeshes(mesh1,mesh2,eps)
+        self.assertEqual(meshinter.getNodalConnectivity().getValues(), [5, 3, 4, 5, 5, 4, 1, 2, 3, 5, 5, 3, 0, 4])
+        self.assertEqual(meshinter.getNodalConnectivityIndex().getValues(), [0, 4, 10, 14])
+        self.assertEqual(res2Back.getValues(), [0, 0, 0])
+        self.assertEqual(res2Tool.getValues(), [0, -1, -1])
+
+        # Second case: all three points of mesh2 touching edges of mesh1
+        coo2[2, 0] = 1.0
+        mesh2.setConnectivity(c, cI)
+
+        meshinter, res2Back, res2Tool = MEDCouplingUMesh.Intersect2DMeshes(mesh1,mesh2,eps)
+        self.assertEqual(meshinter.getNodalConnectivity().getValues(), [5, 3, 4, 5, 5, 4, 1, 5, 5, 5, 2, 3, 5, 3, 0, 4])
+        self.assertEqual(meshinter.getNodalConnectivityIndex().getValues(), [0, 4, 8, 12, 16])
+        self.assertEqual(res2Back.getValues(), [0,0,0,0])
+        self.assertEqual(res2Tool.getValues(), [0,-1,-1,-1])
+        pass
+
     def testSwig2Intersect2DMeshWith1DLine1(self):
         """A basic test with no colinearity between m1 and m2."""
         i=MEDCouplingIMesh("mesh",2,[5,5],[0.,0.],[1.,1.])
@@ -1121,30 +1181,36 @@ class MEDCouplingIntersectTest(unittest.TestCase):
 
     def testSwig2Intersect2DMeshWith1DLine18(self):
         """ Rare case: a *closed* line used as a tool, with the closing point inside a 2D cell ... """
-        tool = MEDCouplingUMesh('circle', 1)
-        coo = DataArrayDouble([(39.35,0),(27.8247,27.8247),(2.40949e-15,39.35),(-27.8247,27.8247),(-39.35,4.81899e-15),(-27.8247,-27.8247),(-7.22848e-15,-39.35),(27.8247,-27.8247),(39.35,7.39805e-15)])
-        tool.setCoords(coo)
-        c = DataArrayInt([2, 3, 5, 8,    2, 5, 3, 4])
-        cI = DataArrayInt([0, 4, 8])
-        tool.setConnectivity(c, cI)
-
-        meh = MEDCouplingUMesh('meh', 2)
-        coo = DataArrayDouble([(-26.4275,36.6199),(-23.5868,31.6996),(-34.1861,41.0993),(-30.3383,25.0214),(-40.1861,30.707),(-35.2622,27.8642),(-37.1861,35.9032),(-30.3068,38.8596),(-25.0071,34.1598),(-26.9625,28.3605),(-25.7138,32.5128),(-27.354,36.4726),(-36.9138,32.5128),(-27.354,28.553),(-26.8908,36.5462),(-28.8461,26.7872)])
-        meh.setCoords(coo)
-        c = DataArrayInt([32, 0, 1, 3, 13, 11, 8, 9, 15, 10, 14,
-                             32, 3, 4, 2, 0, 11, 13, 5, 6, 7, 14, 12, 15])
-        cI = DataArrayInt([0, 11, 24])
-        meh.setConnectivity(c, cI)
-
-        res2D, res1D, m1, m2 = MEDCouplingUMesh.Intersect2DMeshWith1DLine(meh, tool, 1e-12)
-        self.assertEqual(4, res2D.getNumberOfCells())
-        self.assertEqual(res2D.getNodalConnectivity().getValues(),[32, 13, 11, 0, 1, 25, 19, 26, 33, 34, 35, 36, 37, 38, 39, 32, 3, 26, 19, 25, 40, 41, 42, 43,
-                                                                   32, 4, 2, 0, 11, 13, 26, 27, 44, 45, 46, 47, 48, 49, 50, 32, 3, 27, 26, 51, 52, 53])
-        self.assertEqual(res2D.getNodalConnectivityIndex().getValues(),[0, 15, 24, 39, 46])
-        self.assertEqual(res1D.getNodalConnectivity().getValues(),[2, 19, 25, 28, 2, 25, 21, 29, 2, 21, 27, 30, 2, 27, 26, 31, 2, 26, 19, 32])
-        self.assertEqual(res1D.getNodalConnectivityIndex().getValues(),[0, 4, 8, 12, 16, 20])
-        self.assertEqual(m1.getValues(), [0,0,1,1])
-        self.assertEqual(m2.getValues(), [0,1,   -1,-1,  -1,-1,   2,3,  0,1])
+        eps = 1.0e-8
+        m_2d = MEDCouplingUMesh('Plan_IJ_IK_JK', 2)
+        coo = DataArrayDouble([(-2.55102,0.510204),(-1.53061,0.510204),(-0.510204,0.510204),(0.510204,0.510204),(1.53061,0.510204),(2.55102,0.510204),
+            (-2.55102,1.53061),(-1.53061,1.53061),(-0.510204,1.53061),(0.510204,1.53061),(1.53061,1.53061),(2.55102,1.53061),(-2.55102,2.55102),
+            (-1.53061,2.55102),(-0.510204,2.55102),(0.510204,2.55102),(1.53061,2.55102),(2.55102,2.55102)])
+        m_2d.setCoords(coo)
+        c = DataArrayInt([4, 1, 0, 6, 7, 4, 2, 1, 7, 8, 4, 3, 2, 8, 9, 4, 4, 3, 9, 10, 4, 5, 4, 10, 11, 4, 7, 6, 12, 13, 4, 8, 7, 13, 14, 4,
+            9, 8, 14, 15, 4, 10, 9, 15, 16, 4, 11, 10, 16, 17])
+        cI = DataArrayInt([0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50])
+        m_2d.setConnectivity(c, cI)
+        m_2d.checkConsistency()
+
+        m_1d = MEDCouplingUMesh('Slice3DSurf', 1)
+        coo = DataArrayDouble([(0.578179,2.05076),(-1.06484,1.77411),(1.16576,1.86369),(1.52253,1.65688),(-1.69839,1.14578),(1.7736,0.942747)])
+        m_1d.setCoords(coo)
+        c = DataArrayInt([1, 0, 1, 1, 1, 4, 1, 4, 5, 1, 5, 3, 1, 3, 2, 1, 2, 0])
+        cI = DataArrayInt([0, 3, 6, 9, 12, 15, 18])
+        m_1d.setConnectivity(c, cI)
+        m_1d.checkConsistency()
+
+        res2D, res1D, m1, m2 = MEDCouplingUMesh.Intersect2DMeshWith1DLine(m_2d, m_1d, eps)
+        res1D.writeVTK("/tmp/res1d.vtu")
+        self.assertEqual(20, res2D.getNumberOfCells())
+        self.assertEqual(res2D.getNodalConnectivity().getValues(),[4, 7, 6, 12, 13, 5, 10, 9, 32, 18, 20, 21, 33, 5, 15, 16, 33, 21, 20, 18, 32, 5, 9, 8, 31, 32, 5, 14, 15, 32, 31, 5, 8, 26, 19, 31, 5, 7, 13, 14, 31, 19, 26, 5, 24, 26, 8, 27,
+                                                                   25, 5, 2, 1, 25, 27, 5, 7, 26, 24, 5, 25, 22, 24, 5, 1, 0, 6, 7, 24, 22, 25, 5, 8, 9, 28, 27, 5, 3, 2, 27, 28, 5, 9, 10, 29, 28, 5, 4, 3, 28, 29, 5, 10, 30, 23, 29, 5, 11, 5, 4, 29, 23, 30, 5, 10, 33, 30, 5, 16, 17, 11, 30, 33])
+        self.assertEqual(res2D.getNodalConnectivityIndex().getValues(),[0, 5, 13, 21, 26, 31, 36, 43, 49, 54, 58, 62, 70, 75, 80, 85, 90, 95, 102, 106, 112])
+        self.assertEqual(res1D.getNodalConnectivity().getValues(),[1, 18, 32, 1, 32, 31, 1, 31, 19, 1, 19, 26, 1, 26, 24, 1, 24, 22, 1, 22, 25, 1, 25, 27, 1, 27, 28, 1, 28, 29, 1, 29, 23, 1, 23, 30, 1, 30, 33, 1, 33, 21, 1, 21, 20, 1, 20, 18])
+        self.assertEqual(res1D.getNodalConnectivityIndex().getValues(),[0, 3, 6, 9, 12, 15, 18, 21, 24, 27, 30, 33, 36, 39, 42, 45, 48])
+        self.assertEqual(m1.getValues(), [5, 8, 8, 7, 7, 6, 6, 1, 1, 1, 0, 0, 2, 2, 3, 3, 4, 4, 9, 9])
+        self.assertEqual(m2.getValues(), [1, 2, 3, 4, 5, 6, 5, 6, 7, 9, 10, 11, 10, 11, 7, 8, 12, 13, 14, 15, 16, 17, 16, 17, 18, 19, 1, 2, 1, 2, 1, 2])
 
     def testSwig2Intersect2DMeshWith1DLine19(self):
         """ Intersection arc of circle / segment was not properly detecting tangent cases """
@@ -1201,6 +1267,44 @@ class MEDCouplingIntersectTest(unittest.TestCase):
         self.assertEqual(mapLeftRight.getValues(), [-1, -1, 0, 4, -1, -1, 1, 4, -1, -1, 1, 3, -1, -1, 2, 3, -1, -1])
         pass
 
+    def testSwig2Intersect2DMeshWith1DLine21(self):
+        """ A line intersecting a cell very close to one of its node (collinearity not detected) """
+        eps=1.0e-5  # was working at 1.0e-8, but should also really work with 1.0e-5
+        mesh = MEDCouplingUMesh('mesh', 2)
+        coo = DataArrayDouble([(110.65324,180.56968),(112.01128,182.78580),(113.36932,185.00192),(118.27200,181.90669),(118.27200,178.79852),(118.27200,175.67380)])
+        mesh.setCoords(coo)
+        c = DataArrayInt([NORM_QUAD4, 0, 1, 4, 5, NORM_QUAD4, 1, 2, 3, 4])
+        cI = DataArrayInt([0, 5, 10])
+        mesh.setConnectivity(c, cI)
+
+        tool = MEDCouplingUMesh('tool', 1)
+        coo = DataArrayDouble([(0.0, 182.78400),  (182.78400, 182.78400)])
+        tool.setCoords(coo)
+        c = DataArrayInt([NORM_SEG2,0,1])
+        cI = DataArrayInt([0, 3])
+        tool.setConnectivity(c, cI)
+
+        res2D, res1D, resToSelf, mapLeftRight = MEDCouplingUMesh.Intersect2DMeshWith1DLine(mesh, tool, eps)
+        self.assertEqual(res1D.getNodalConnectivity().getValues(), [1, 6, 1, 1, 1, 8, 1, 8, 9, 1, 9, 7])
+        self.assertEqual(res1D.getNodalConnectivityIndex().getValues(),[0, 3, 6, 9, 12])
+        self.assertEqual(res2D.getNodalConnectivity().getValues(), [5, 0, 1, 8, 4, 5, 5, 1, 2, 9, 8, 5, 3, 4, 8, 9])
+        self.assertEqual(res2D.getNodalConnectivityIndex().getValues(),[0, 6, 11, 16])
+
+        self.assertEqual(resToSelf.getValues(), [0, 1, 1])
+        self.assertEqual(mapLeftRight.getValues(), [-1, -1, 1, 0, 1, 2, -1, -1])
+
+    def testSwig2Intersect2DMeshWith1DLine22(self):
+        """ Simple check that an execption is thrown if 1D mesh is not a single line (one point
+        connected to more than 2 edges. """
+        m2 = MEDCouplingCMesh.New()
+        da = DataArrayDouble(3)
+        da.iota()
+        m2.setCoords(da, da)
+        m2 = m2.buildUnstructured()
+        m1,_,_,_,_ = m2.buildDescendingConnectivity() # m1 will have it central point connected to 4 edges
+        # MEDCouplingUMesh.Intersect2DMeshWith1DLine(m2, m1, 1.0e-8)
+        self.assertRaises(InterpKernelException, MEDCouplingUMesh.Intersect2DMeshWith1DLine, m2, m1, 1.0e-8)
+
     def testSwig2Conformize2D1(self):
         eps = 1.0e-8
         coo = [0.,-0.5,0.,0.,0.5,0.,0.5,-0.5,0.25,