# -*- coding: utf-8 -*-
-# Copyright (C) 2007-2019 CEA/DEN, EDF R&D
+# Copyright (C) 2007-2023 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
from MEDCouplingDataForTest import MEDCouplingDataForTest
-import rlcompleter,readline # this line has to be here, to ensure a usability of MEDCoupling/MEDLoader. B4 removing it please notify to anthony.geay@cea.fr
+import rlcompleter,readline # this line has to be here,to ensure a usability of MEDCoupling/MEDLoader. B4 removing it please notify to anthony.geay@cea.fr
class MEDCouplingBasicsTest7(unittest.TestCase):
def testDAIAggregateMulti1(self):
a=DataArrayInt64.New()
- a.setValues(list(range(4)), 2, 2)
+ a.setValues(list(range(4)),2, 2)
a.setName("aa")
b=DataArrayInt64.New()
b.setValues(list(range(6)), 3, 2)
# non regression test in python wrapping
rg=DataArrayInt64([0,10,29,56,75,102,121,148,167,194,213,240,259,286,305,332,351,378,397,424,443,470,489,516])
a,b,c=DataArrayInt64([75]).splitByValueRange(rg)
- assert(a.isEqual(DataArrayInt64([4])))
- assert(b.isEqual(DataArrayInt64([0])))
- assert(c.isEqual(DataArrayInt64([4])))
+ self.assertTrue(a.isEqual(DataArrayInt64([4])))
+ self.assertTrue(b.isEqual(DataArrayInt64([0])))
+ self.assertTrue(c.isEqual(DataArrayInt64([4])))
pass
def testDAIBuildExplicitArrByRanges1(self):
self.assertRaises(InterpKernelException,a.indicesOfSubPart,f) # 12 in f does not exist in a
pass
+ def testDAIsortToHaveConsecutivePairs(self):
+ dref=DataArrayInt64([(6, 216), (216, 218), (218, 220), (220, 222), (222, 224), (224, 226)])
+ dtest=DataArrayInt64([(6, 216), (218, 216), (224, 226), (222, 220), (218, 220), (222, 224)])
+ dtest.sortToHaveConsecutivePairs()
+ self.assertTrue(dtest.isEqual(dref))
+
def testDAIFromLinkedListOfPairToList1(self):
d=DataArrayInt64([(5,7),(7,3),(3,12),(12,17)])
zeRes=DataArrayInt64([5,7,3,12,17])
self.assertTrue(f_voro.getArray().isEqual(f.getArray(),1e-8))
pass
+ def testDAIOccurenceRankInThis(self):
+ arr=DataArrayInt([5,3,2,1,4,5,2,1,0,11,5,4])
+ self.assertTrue(arr.occurenceRankInThis().isEqual(DataArrayInt([0,0,0,0,0,1,1,1,0,0,2,1])))
+
+ def testDAIFindPermutationFromFirstToSecondDuplicate(self):
+ arr0 = DataArrayInt([5,3,2,1,4,5,2,1,0,11,5,4])
+ arr1 = DataArrayInt([0,1,1,2,2,3,4,4,5,5,5,11])
+ self.assertTrue(DataArrayInt.FindPermutationFromFirstToSecondDuplicate(arr0,arr1).isEqual(DataArrayInt([8,5,3,1,6,9,4,2,0,11,10,7])))
+ self.assertTrue(DataArrayInt.FindPermutationFromFirstToSecondDuplicate(arr1,arr0).isEqual(DataArrayInt([8,3,7,2,6,1,4,11,0,5,10,9])))
+
+ def testDAIIndexOfSameConsecutiveValueGroups(self):
+ arr = DataArrayInt([0,1,1,2,2,3,4,4,5,5,5,11])
+ self.assertTrue(arr.indexOfSameConsecutiveValueGroups().isEqual(DataArrayInt([0,1,3,5,6,8,11,12])))
+
+ def testSkyLineGroupPacks(self):
+ arr = DataArrayInt([1,4,5,0,2,4,5,6,1,3,5,6,7,2,6,7,0,1,5,8,9,0,1,2,4,6,8,9,10,1,2,3,5,7,9,10,11,2,3,6,10,11,4,5,9,12,13,4,5,6,8,10,12,13,14,5,6,7,9,11,13,14,15,6,7,10,14,15,8,9,13,8,9,10,12,14,9,10,11,13,15,10,11,14])
+ arrI = DataArrayInt([0,3,8,13,16,21,29,37,42,47,55,63,68,71,76,81,84])
+ sk = MEDCouplingSkyLineArray(arrI,arr)
+ part = DataArrayInt([0,3,4,7,16])
+ sk2 = sk.groupPacks(part)
+ self.assertTrue(sk2.getValuesArray().isEqual(arr))
+ self.assertTrue(sk2.getIndexArray().isEqual(DataArrayInt([0,13,16,37,84])))
+
+ def testSkyLineUniqueNotSortedByPack(self):
+ arrI = DataArrayInt([0,3,9,15,18,24,36,48,54])
+ arr = DataArrayInt([1,4,5,0,4,5,2,5,6,3,6,7,1,5,6,2,6,7,0,1,5,5,8,9,0,1,4,6,9,10,1,2,4,6,8,9,2,3,5,7,9,10,1,2,5,7,10,11,2,3,6,6,10,11])
+ sk = MEDCouplingSkyLineArray(arrI,arr)
+ sk2 = sk.uniqueNotSortedByPack()
+ self.assertTrue(sk2.getIndexArray().isEqual(DataArrayInt([0,3,8,13,16,21,29,37,42])))
+ self.assertTrue(sk2.getValuesArray().isEqual(DataArrayInt([1,4,5,0,2,4,5,6,1,3,5,6,7,2,6,7,0,1,5,8,9,0,1,2,4,6,8,9,10,1,2,3,5,7,9,10,11,2,3,6,10,11])))
+
+ def testSkyLineAggregatePacks1(self):
+ arr = DataArrayDouble(3) ; arr.iota()
+ m = MEDCouplingCMesh() ; m.setCoords(arr,arr) ; m = m.buildUnstructured()
+ a,b = m.computeEnlargedNeighborsOfNodes()
+ sk = MEDCouplingSkyLineArray(b,a)
+ sk1 = sk.deepCopy()
+ sk1.getValuesArray()[:] *= 2
+ sk2 = sk.deepCopy()
+ sk2.getValuesArray()[:] *= 3
+ skOut = MEDCouplingSkyLineArray.AggregatePacks([sk,sk1,sk2])
+ self.assertTrue(skOut.getIndexArray().isEqual(DataArrayInt([0,9,24,33,48,72,87,96,111,120])))
+ self.assertTrue(skOut.getValuesArray().isEqual(DataArrayInt([1,3,4,2,6,8,3,9,12,0,2,3,4,5,0,4,6,8,10,0,6,9,12,15,1,4,5,2,8,10,3,12,15,0,1,4,6,7,0,2,8,12,14,0,3,12,18,21,0,1,2,3,5,6,7,8,0,2,4,6,10,12,14,16,0,3,6,9,15,18,21,24,1,2,4,7,8,2,4,8,14,16,3,6,12,21,24,3,4,7,6,8,14,9,12,21,3,4,5,6,8,6,8,10,12,16,9,12,15,18,24,4,5,7,8,10,14,12,15,21])))
+
+ def testDACopySorted1(self):
+ d = DataArrayInt32([5,1,100,20])
+ self.assertTrue(d.copySorted().isEqual(DataArrayInt32([1,5,20,100])))
+ d = DataArrayInt64([5,1,100,20])
+ self.assertTrue(d.copySorted().isEqual(DataArrayInt64([1,5,20,100])))
+ d = DataArrayInt([5,1,100,20])
+ self.assertTrue(d.copySorted().isEqual(DataArrayInt([1,5,20,100])))
+ d = DataArrayDouble([5,1,100,20])
+ self.assertTrue(d.copySorted().isEqual(DataArrayDouble([1,5,20,100]),1e-10))
+
+ def testFieldAreStrictlyCompatible(self):
+ arr=DataArrayDouble(10) ; arr.iota()
+ m=MEDCouplingCMesh() ; m.setCoords(arr,arr) ; m=m.buildUnstructured()
+ f=MEDCouplingFieldDouble(ON_CELLS) ; f.setMesh(m)
+ f2=MEDCouplingFieldDouble(ON_CELLS) ; f2.setMesh(m)
+ self.assertTrue(f.areStrictlyCompatible(f2))
+ self.assertTrue(f.areStrictlyCompatibleForMulDiv(f2))
+ f2.setMesh(f2.getMesh().deepCopy())
+ self.assertTrue(not f.areStrictlyCompatible(f2))
+ self.assertTrue(not f.areStrictlyCompatibleForMulDiv(f2))
+ f3=MEDCouplingFieldDouble(ON_NODES) ; f3.setMesh(m)
+ self.assertTrue(not f.areStrictlyCompatible(f3))
+ self.assertTrue(not f.areStrictlyCompatibleForMulDiv(f3))
+
+ def testBugZipConnectivityTraducer(self):
+ """
+ Non regression test : here cell #1 and cell #2 are nearly the same but not the same. zipConnectivityTraducer called by areCellsIncludedIn
+ failed to capture that.
+ """
+ coo = DataArrayDouble([0,1,2,3],4,1)
+ m = MEDCouplingUMesh("",1)
+ m.setCoords(coo)
+ m.allocateCells()
+ m.insertNextCell(NORM_SEG2,[0,1])
+ m.insertNextCell(NORM_SEG2,[1,2])
+ m.insertNextCell(NORM_SEG2,[2,1])
+ #
+ a,b = m.areCellsIncludedIn(m,0)
+ self.assertTrue(a)
+ self.assertTrue(b.isIota(3))
+ #
+ self.assertTrue(m.deepCopy().zipConnectivityTraducer(0).isIota(3))
+ self.assertTrue(m.deepCopy().zipConnectivityTraducer(1).isIota(3))
+ self.assertTrue(m.deepCopy().zipConnectivityTraducer(2).isEqual(DataArrayInt([0,1,1])))
+
+ def testBugAreCellsIncludedIn1(self):
+ """
+ Non regression test: a.areCellsIncludedIn(b) was buggy when some cells in b were duplicated into a following specified policy.
+ """
+ coo = DataArrayDouble([0,1,2,3,4],5,1)
+ m = MEDCouplingUMesh("",1)
+ m.setCoords(coo)
+ m.allocateCells()
+ # m contains several duplicated cells - some of those duplicated cells will be in m1
+ for i in range(3):
+ m.insertNextCell(NORM_SEG2,[0,1])
+ for i in range(4):
+ m.insertNextCell(NORM_SEG2,[1,2])
+ for i in range(2):
+ m.insertNextCell(NORM_SEG2,[2,3])
+ for i in range(2):
+ m.insertNextCell(NORM_SEG2,[3,4])
+ #
+ bexp = DataArrayInt([0,1,2, 3,4,5,6, 9,10])
+ m1 = m[bexp]
+ #
+ a,b = m.areCellsIncludedIn(m1,0)
+ self.assertTrue(a)
+ self.assertTrue(b.isEqual(DataArrayInt([2,2,2, 6,6,6,6, 10,10])))
+ #
+ bexp2 = DataArrayInt([0,1,2, 3,4,0, 6, 9,10])
+ m2 = m[bexp2]
+ a,b = m.areCellsIncludedIn(m2,0)
+ self.assertTrue(a)
+ self.assertTrue(b.isEqual(DataArrayInt([2,2,2,6,6,2,6,10,10])))
+
+ def testSkyLineArrayThreshold(self):
+ x = DataArrayInt([0, 1, 2, 11, 12, 13, 3, 4, 5, 6, 14, 15, 16, 17, 9, 10, 18, 19])
+ xi = DataArrayInt([0, 6, 14, 18])
+ sk = MEDCouplingSkyLineArray(xi,x)
+ lsk,rsk = sk.thresholdPerPack(11)
+ self.assertTrue(lsk.getValuesArray().isEqual(DataArrayInt([0, 1, 2, 3, 4, 5, 6, 9, 10])))
+ self.assertTrue(lsk.getIndexArray().isEqual(DataArrayInt([0, 3, 7, 9])))
+ self.assertTrue(rsk.getValuesArray().isEqual(DataArrayInt([11, 12, 13, 14, 15, 16, 17, 18, 19])))
+ self.assertTrue(rsk.getIndexArray().isEqual(DataArrayInt([0, 3, 7, 9])))
+
+ def testPenta18GaussNE(self):
+ conn = [1,0,2,4,3,5,6,7,8,9,13,14,11,10,15,12,17,16]
+ coo = DataArrayDouble([(27.237499999999997, 9.8, 0.0), (26.974999999999994, 9.8, 0.0), (27.111517409545634, 9.532083869948877, 0.0), (27.237499999999997, 9.8, 0.5000000000000001), (26.974999999999994, 9.8, 0.5000000000000002), (27.111517409545634, 9.532083869948877, 0.5), (27.106249999999996, 9.8, 0.0), (27.17450870477282, 9.666041934974439, 0.0), (27.04325870477281, 9.666041934974439, 0.0), (27.106249999999996, 9.8, 0.5000000000000001), (27.237499999999997, 9.8, 0.25), (26.974999999999994, 9.8, 0.2500000000000001), (27.106249999999996, 9.8, 0.2500000000000001), (27.174508704772816, 9.666041934974439, 0.5), (27.043258704772814, 9.666041934974439, 0.5000000000000001), (27.111517409545634, 9.532083869948877, 0.25), (27.043258704772818, 9.666041934974436, 0.25000000000000006), (27.174508704772816, 9.666041934974436, 0.25)])
+ m = MEDCouplingUMesh("mesh",3)
+ m.setCoords(coo)
+ m.allocateCells()
+ m.insertNextCell(NORM_PENTA18,conn)
+ f = MEDCouplingFieldDouble(ON_GAUSS_NE)
+ f.setMesh(m)
+ f.setArray(DataArrayDouble(18*[0.]))
+ self.assertTrue(f.getLocalizationOfDiscr().isEqual(coo[conn],1e-10))
+
+ def testDADEigenValuesPb(self):
+ """EDF22126 : eigen values with Identity matrix returned nan. Now it returns correct eigen values 1.0 """
+ valuesExp = DataArrayDouble([(1.,1.,1.),(2.,-1.,0.),(2.,0.,1.),(3.,0.,0.)])
+ d = DataArrayDouble(4, 6)
+ for i,(v0, v1, v2, v3, v4, v5,) in enumerate([
+ (1, 1, 1, 0, 0, 0),
+ (1, 0, 0, 1, 0, 1),
+ (1, 1, 1, 0, 1, 0),
+ (1, 1, 1, 1, 1, 1)]):
+ d[i] = [v0, v1, v2, v3, v4, v5]
+ self.assertTrue(d.eigenValues().isEqual(valuesExp,1e-12))
+ pass
+
+ def testBugOnReprOf1SGTUMesh(self):
+ """ Non reg bug on repr of empty MEDCoupling1SGTUMesh instance """
+ m = MEDCoupling1SGTUMesh()
+ m.simpleRepr()
+ str(m)
+ m.advancedRepr()
+ repr(m)
+ m = MEDCoupling1DGTUMesh()
+ m.simpleRepr()
+ str(m)
+ m.advancedRepr()
+ repr(m)
+
+ def testCheckGeomConsistency0(self):
+ """Test of MEDCouplingUMesh.checkGeomConsistency"""
+ m = MEDCouplingUMesh("",2)
+ m.setCoords(DataArrayDouble([(0,0),(1,0),(2,0),(0,1),(1,1),(2,1)]))
+ m.allocateCells()
+ m.insertNextCell(NORM_TRI6,[0,1,2,3,4,5])
+ m.insertNextCell(NORM_TRI6,[0,1,3,3,4,5])
+ m.checkConsistency()
+ self.assertRaises(InterpKernelException,m.checkGeomConsistency) # cell1 is incorrect because node 3 is repeated twice
+ m.getNodalConnectivity()[10]=2 # replace 3 by 2 for cell#1 to fix the problem
+ m.checkConsistency()
+ m.checkGeomConsistency() # now m is OK
+
+ def testInt32Int64Arr0(self):
+ n = 30
+ arr = DataArrayInt32(n)
+ arr.iota() ; arr.rearrange(3)
+ comps = ["a","bb","ccc"]
+ name = "aaa"
+ arr.setName(name) ; arr.setInfoOnComponents(comps)
+ self.assertEqual(arr.accumulate(),[135, 145, 155])
+ arr2 = arr.convertToInt64Arr() # test is here
+ self.assertEqual(arr2.accumulate(),[135, 145, 155])
+ self.assertTrue(isinstance(arr2,DataArrayInt64))
+ self.assertEqual(arr2.getName(),name)
+ self.assertEqual(arr2.getInfoOnComponents(),comps)
+ arr3 = arr2.convertToInt32Arr() # test is here
+ self.assertEqual(arr3.accumulate(),[135, 145, 155])
+ self.assertTrue(isinstance(arr3,DataArrayInt32))
+ self.assertEqual(arr3.getName(),name)
+ self.assertEqual(arr3.getInfoOnComponents(),comps)
+ self.assertTrue(arr3.isEqual(arr))
+
+ def testComputeMeshCenterOfMass0(self):
+ #2D
+ arr = DataArrayDouble(5) ; arr.iota()
+ m = MEDCouplingCMesh() ; m.setCoords(arr,arr) ; m=m.buildUnstructured()
+ self.assertTrue( m.computeMeshCenterOfMass().isEqual(DataArrayDouble([2,2],1,2),1e-12) )
+ #3D
+ m = MEDCouplingCMesh() ; m.setCoords(arr,arr,arr) ; m=m.buildUnstructured()
+ self.assertTrue( m.computeMeshCenterOfMass().isEqual(DataArrayDouble([2,2,2],1,3),1e-12) )
+
+ def testBugPenta15_0(self):
+ """
+ Non regression test from Roberto Da Via pointing error in connectivity of 5th sub face of penta15 cell.
+ """
+ coo=DataArrayDouble([
+ (0,1,1),(0,0,1),(1,0,1),
+ (0,1,0),(0,0,0),(1,0,0),
+ (0,0.5,1),(0.5,0,1),(0.5,0.5,1),
+ (0,0.5,0),(0.5,0,0),(0.5,0.5,0),
+ (0,1,0.5),(0,0,0.5),(1,0,0.5)
+ ])
+
+ m = MEDCouplingUMesh("penta15",3)
+ m.setCoords(coo)
+ m.allocateCells()
+ m.insertNextCell(NORM_PENTA15,list(range(15)))
+ bm = m.buildBoundaryMesh(True)
+ bm.writeVTK("boundary.vtu")
+ conn_expected = [
+ [6, 0, 1, 2, 6, 7, 8],
+ [6, 3, 5, 4, 11, 10, 9],
+ [8, 0, 3, 4, 1, 12, 9, 13, 6],
+ [8, 1, 4, 5, 2, 13, 10, 14, 7],
+ [8, 2, 5, 3, 0, 14, 11, 12, 8] # old = [8, 2, 4, 5, 0, 14, 11, 12, 8]
+ ]
+ self.assertTrue( bm.getNodalConnectivity().isEqual(DataArrayInt(sum(conn_expected,[]))) )
+
+ def testBugWithPolyhedInterpWithMoreThan255Nodes(self):
+ """
+ [EDF25207] : Check interpolation containing polyhedron with more than 255 nodes is OK at bbox computation stage
+ """
+ n = 8
+ arr = DataArrayDouble(n) ; arr.iota()
+ m = MEDCouplingCMesh()
+ m.setCoords(arr,arr,arr)
+ m = m.buildUnstructured()
+ skin = m.computeSkin()
+ skin.zipCoords()
+ # check that skin contains more than 2**8-1 node to reveal bug
+ self.assertTrue(skin.getNumberOfNodes()>255)
+ # Build a single polyhedron cell from skin
+ skin1 = MEDCoupling1SGTUMesh(skin)
+ conn = skin1.getNodalConnectivity()
+ conn.rearrange( MEDCouplingUMesh.GetNumberOfNodesOfGeometricType(skin1.getCellModelEnum()) )
+ connPolyhed = conn.changeNbOfComponents(MEDCouplingUMesh.GetNumberOfNodesOfGeometricType(skin1.getCellModelEnum())+1,-1)
+ connPolyhed.rearrange(1)
+ connPolyhed.popBackSilent()
+ meshSinglePolyhed = MEDCouplingUMesh("",3)
+ meshSinglePolyhed.allocateCells()
+ meshSinglePolyhed.insertNextCell(NORM_POLYHED,connPolyhed.getValues())
+ meshSinglePolyhed.setCoords( skin1.getCoords() )
+
+ rem = MEDCouplingRemapper()
+ rem.prepare(meshSinglePolyhed,m,"P0P0")
+ res = rem.getCrudeMatrix()
+ self.assertTrue( all([len(elt)==1 for elt in res]) )
+ self.assertTrue( all([elt[0]>0.99 and elt[0]<1.01 for elt in res]) )
+
+ @unittest.skipUnless(MEDCouplingHasNumPyBindings(),"requires numpy")
+ def testShapeFuncAndDerivative0(self):
+ """
+ Test values returned by MEDCoupling on HEXA27 element of shape function and its derivatives.
+ See https://www.code-aster.org/V2/doc/v10/fr/man_r/r3/r3.01.01.pdf
+ """
+ import numpy as np
+
+ ref_coords_hexa27_med = [[-1.0, -1.0, -1.0], [-1.0, 1.0, -1.0], [1.0, 1.0, -1.0], [1.0, -1.0, -1.0], [-1.0, -1.0, 1.0], [-1.0, 1.0, 1.0], [1.0, 1.0, 1.0], [1.0, -1.0, 1.0], [-1.0, 0.0, -1.0], [0.0, 1.0, -1.0], [1.0, 0.0, -1.0], [0.0, -1.0, -1.0], [-1.0, 0.0, 1.0], [0.0, 1.0, 1.0], [1.0, 0.0, 1.0], [0.0, -1.0, 1.0], [-1.0, -1.0, 0.0], [-1.0, 1.0, 0.0], [1.0, 1.0, 0.0], [1.0, -1.0, 0.0], [0.0, 0.0, -1.0], [-1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [1.0, 0.0, 0.0], [0.0, -1.0, 0.0], [0.0, 0.0, 1.0], [0.0, 0.0, 0.0]]
+
+ def coor2index(coor):
+ zeMap = {-1.0 : 0, 0.0 : 2 , 1.0 : 1}
+ return zeMap[coor]
+
+ vcoor2index = np.vectorize( coor2index )
+ node2ijk_hexa27_med = vcoor2index( np.array(ref_coords_hexa27_med) )
+
+ def N_1d_quad(x):
+ return np.array([-0.5*x*(1-x), 0.5*x*(x+1), 1.-x*x])
+
+ def N_3d_hexa27(x, i, j, k):
+ return N_1d_quad(x[0])[i]*N_1d_quad(x[1])[j]*N_1d_quad(x[2])[k]
+
+ def N_hexa27(x):
+ return np.array([N_3d_hexa27(x, *node2ijk_hexa27_med[node,:]) for node in range(27)])
+
+ # Implementing shape function derivatives
+ def diff_N_1d_quad(x):
+ return np.array([x-0.5, x+0.5, -2.*x])
+
+ def diff_N_3d_hexa27(x, i, j, k):
+ return np.array([diff_N_1d_quad(x[0])[i]*N_1d_quad(x[1])[j] *N_1d_quad(x[2])[k],
+ N_1d_quad(x[0])[i] *diff_N_1d_quad(x[1])[j]*N_1d_quad(x[2])[k],
+ N_1d_quad(x[0])[i] *N_1d_quad(x[1])[j] *diff_N_1d_quad(x[2])[k]])
+
+ def diff_N_hexa27(x):
+ return np.array([diff_N_3d_hexa27(x, *node2ijk_hexa27_med[node,:]) for node in range(27)])
+ # computation of ref values
+ posInRefCoord = [-0.85685375,-0.90643355,-0.90796825]
+ ref = N_hexa27( np.array(posInRefCoord) )
+ ref2 = diff_N_hexa27( np.array(posInRefCoord) )
+ # computation using MEDCoupling
+ gl = MEDCouplingGaussLocalization(NORM_HEXA27,sum(ref_coords_hexa27_med,[]),posInRefCoord,[1])
+ mcShapeFunc = gl.getShapeFunctionValues()
+ mcShapeFunc.rearrange(1)
+ self.assertTrue( mcShapeFunc.isEqual(DataArrayDouble(ref),1e-12) )
+
+ mvDevOfShapeFunc = gl.getDerivativeOfShapeFunctionValues()
+ mvDevOfShapeFunc.rearrange(1)
+ ref2_mc = DataArrayDouble(ref2) ; ref2_mc.rearrange(1)
+ self.assertTrue( mvDevOfShapeFunc.isEqual( ref2_mc, 1e-12) )
+
+ def testShapeFuncAndDerivative1(self):
+ """
+ This test focus
+ """
+ def GetShapeFunc(ref_coord,vec):
+ gl3 = MEDCouplingGaussLocalization(gt,sum(ref_coord,[]), vec, [1])
+ funVal = gl3.getShapeFunctionValues()
+ funVal.rearrange(1)
+ return funVal
+
+ def GetDerivative(ref_coord,vec):
+ gl3 = MEDCouplingGaussLocalization(gt,sum(ref_coord,[]), vec, [1])
+ funVal = gl3.getDerivativeOfShapeFunctionValues()
+ return funVal
+ vec = [-0.85685375,-0.90643355,-0.90796825]
+ eps = 1e-6
+ # 3D cells
+ for gt in [NORM_TETRA4,NORM_TETRA10,NORM_HEXA8,NORM_PENTA6,NORM_PYRA5,NORM_PYRA13,NORM_PENTA15,NORM_PENTA6,NORM_PENTA18,NORM_HEXA20,NORM_HEXA27]: # type of cell for which derivatives are implemented
+ ref_coord = [list(elt) for elt in MEDCouplingGaussLocalization.GetDefaultReferenceCoordinatesOf(gt).getValuesAsTuple()]
+
+ der_computed = GetDerivative(ref_coord,vec)
+ der_computed.rearrange(3)
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps,vec[1],vec[2]])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_X = der_computed[:,0]-der_deduced
+ delta_X.abs()
+ self.assertTrue(delta_X.findIdsNotInRange(-1e-4,+1e-4).empty())
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1]+eps,vec[2]])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_Y = der_computed[:,1]-der_deduced
+ delta_Y.abs()
+ self.assertTrue(delta_Y.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1],vec[2]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_Z = der_computed[:,2]-der_deduced
+ delta_Z.abs()
+ self.assertTrue(delta_Z.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ for gt,ref_coord in [(NORM_TETRA4,[[0.0, 1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 1.0], [1.0, 0.0, 0.0]]),(NORM_TETRA10,[[0.0, 1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 1.0], [1.0, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5], [0.0, 0.5, 0.5], [0.5, 0.5, 0.0], [0.5, 0.0, 0.0], [0.5, 0.0, 0.5]]),(NORM_HEXA8,[[-1.0, -1.0, -1.0], [-1.0, 1.0, -1.0], [1.0, 1.0, -1.0], [1.0, -1.0, -1.0], [-1.0, -1.0, 1.0], [-1.0, 1.0, 1.0], [1.0, 1.0, 1.0], [1.0, -1.0, 1.0]]),(NORM_HEXA8,[[-1.0, 1.0, 0.0], [-1.0, -1.0, 0.0], [1.0, -1.0, 0.0], [1.0, 1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]),(NORM_HEXA8,[[-1.0, -1.0, 0.0], [-1.0, 1.0, 0.0], [1.0, 1.0, 0.0], [1.0, -1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]),(NORM_PENTA6,[[-1.0, 1.0, 0.0], [-1.0, 0.0, 0.0], [-1.0, -0.0, 1.0], [1.0, 1.0, 0.0], [1.0, 0.0, 0.0], [1.0, 0.0, 1.0]]),(NORM_PENTA6,[[-1.0, 1.0, 0.0], [-1.0, -1.0, 0.0], [1.0, -1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]),(NORM_PENTA6,[[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]),(NORM_PYRA5,[[1.0, 0.0, 0.0], [0.0, -1.0, 0.0], [-1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]),(NORM_PYRA13, [[1.0, 0.0, 0.0], [0.0, -1.0, 0.0], [-1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0], [0.5, -0.5, 0.0], [-0.5, -0.5, 0.0], [-0.5, 0.5, 0.0], [0.5, 0.5, 0.0], [0.5, 0.0, 0.5], [0.0, -0.5, 0.5], [-0.5, 0.0, 0.5], [0.0, 0.5, 0.5]]),(NORM_PENTA15,[[-1.0, 1.0, 0.0], [-1.0, 0.0, 0.0], [-1.0, -0.0, 1.0], [1.0, 1.0, 0.0], [1.0, 0.0, 0.0], [1.0, 0.0, 1.0], [-1.0, 0.5, 0.0], [-1.0, 0.0, 0.5], [-1.0, 0.5, 0.5], [1.0, 0.5, 0.0], [1.0, 0.0, 0.5], [1.0, 0.5, 0.5], [0.0, 1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 1.0]]),(NORM_PENTA18,[[-1.0, 1.0, 0.0], [-1.0, 0.0, 0.0], [-1.0, -0.0, 1.0], [1.0, 1.0, 0.0], [1.0, 0.0, 0.0], [1.0, 0.0, 1.0], [-1.0, 0.5, 0.0], [-1.0, 0.0, 0.5], [-1.0, 0.5, 0.5], [1.0, 0.5, 0.0], [1.0, 0.0, 0.5], [1.0, 0.5, 0.5], [0.0, 1.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 1.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5], [0.0, 0.5, 0.5]]),(NORM_HEXA20,[[-1.0, -1.0, -1.0], [-1.0, 1.0, -1.0], [1.0, 1.0, -1.0], [1.0, -1.0, -1.0], [-1.0, -1.0, 1.0], [-1.0, 1.0, 1.0], [1.0, 1.0, 1.0], [1.0, -1.0, 1.0], [-1.0, 0.0, -1.0], [0.0, 1.0, -1.0], [1.0, 0.0, -1.0], [0.0, -1.0, -1.0], [-1.0, 0.0, 1.0], [0.0, 1.0, 1.0], [1.0, 0.0, 1.0], [0.0, -1.0, 1.0], [-1.0, -1.0, 0.0], [-1.0, 1.0, 0.0], [1.0, 1.0, 0.0], [1.0, -1.0, 0.0]])]: # type of cell for which derivatives are implemented
+
+ der_computed = GetDerivative(ref_coord,vec)
+ der_computed.rearrange(3)
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps,vec[1],vec[2]])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_X = der_computed[:,0]-der_deduced
+ delta_X.abs()
+ self.assertTrue(delta_X.findIdsNotInRange(-1e-4,+1e-4).empty())
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1]+eps,vec[2]])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_Y = der_computed[:,1]-der_deduced
+ delta_Y.abs()
+ self.assertTrue(delta_Y.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1],vec[2]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_Z = der_computed[:,2]-der_deduced
+ delta_Z.abs()
+ self.assertTrue(delta_Z.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ # 2D cells
+ vec = [0.64,0.2]
+
+ for gt in [NORM_QUAD4,NORM_QUAD8,NORM_QUAD9,NORM_TRI3,NORM_TRI6,NORM_TRI7]:
+ ref_coord = [list(elt) for elt in MEDCouplingGaussLocalization.GetDefaultReferenceCoordinatesOf(gt).getValuesAsTuple()]
+
+ der_computed = GetDerivative(ref_coord,vec)
+ der_computed.rearrange(2)
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps,vec[1]])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_X = der_computed[:,0]-der_deduced
+ delta_X.abs()
+ self.assertTrue(delta_X.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_Y = der_computed[:,1]-der_deduced
+ delta_Y.abs()
+ self.assertTrue(delta_Y.findIdsNotInRange(-1e-4,+1e-4).empty())
+
+ # B version of TRI6, QUAD4 and QUAD8
+ for gt,ref_coord in [(NORM_TRI3,[[0.0, 0.0], [1.0, 0.0], [0.0, 1.0]]),(NORM_TRI6,[[0., 0.], [1., 0.], [0., 1.], [0.5, 0.], [0.5, 0.5], [0., 0.5]]),
+ (NORM_QUAD4,[[-1., -1.], [1., -1.], [1., 1.], [-1., 1.]]),(NORM_QUAD4,[[-1., -1.], [-1., 1.], [1., 1.], [1., -1.]]),(NORM_QUAD4,[[-1., 0.], [1., 0.], [0., 0.], [0., 0.]]),(NORM_QUAD8,[[-1., -1.], [1., -1.], [1., 1.], [-1., 1.], [0., -1.], [1., 0.], [0., 1.], [-1., 0.]])]:
+ der_computed = GetDerivative(ref_coord,vec)
+ der_computed.rearrange(2)
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps,vec[1]])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_X = der_computed[:,0]-der_deduced
+ delta_X.abs()
+ self.assertTrue(delta_X.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_Y = der_computed[:,1]-der_deduced
+ delta_Y.abs()
+ self.assertTrue(delta_Y.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ # 1D cells
+ vec = [0.64]
+
+ for gt in [NORM_SEG2,NORM_SEG3,NORM_SEG4]:
+ ref_coord = [list(elt) for elt in MEDCouplingGaussLocalization.GetDefaultReferenceCoordinatesOf(gt).getValuesAsTuple()]
+
+ der_computed = GetDerivative(ref_coord,vec)
+ der_computed.rearrange(1)
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_X = der_computed[:,0]-der_deduced
+ delta_X.abs()
+ self.assertTrue(delta_X.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+ #B version of SEG2
+ for gt,ref_coord in [(NORM_SEG2,[[0.], [1.]])]:
+ der_computed = GetDerivative(ref_coord,vec)
+ der_computed.rearrange(1)
+
+ der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+ delta_X = der_computed[:,0]-der_deduced
+ delta_X.abs()
+ self.assertTrue(delta_X.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+
+ def testComputeTriangleHeight0(self):
+ arr = DataArrayDouble([0,1])
+ m = MEDCouplingCMesh() ; m.setCoords(arr,arr)
+ m = m.buildUnstructured() ; m.simplexize(0) ; m = MEDCoupling1SGTUMesh(m)
+ res = m.computeTriangleHeight()
+ expected = DataArrayDouble([(1.0, 1.0, sqrt(2)/2.0), (sqrt(2)/2.0, 1.0, 1.0)])
+ self.assertTrue( res.isEqual(expected,1e-12) )
+ m.changeSpaceDimension(3,100)
+ res2 = m.computeTriangleHeight()
+ self.assertTrue( res2.isEqual(expected,1e-12) )
+ expected2 = DataArrayDouble([sqrt(2)/2.0, sqrt(2)/2.0])
+ self.assertTrue( res2.minPerTuple().isEqual(expected2,1e-12) )
+
+ def testDAILocateComponentId0(self):
+ arr = DataArrayInt( [(0, 1, 2), (3, 4, 5), (6, 2, 3), (7, 8, 9), (9, 0, 10), (11, 12, 13), (14, 5, 11), (15, 16, 17)] )
+ valToSearchIntoTuples = DataArrayInt( [1, 4, 6, 8, 10, 12, 14, 16, 17] )
+ tupleIdHint = DataArrayInt( [0, 1, 2, 3, 4, 5, 6, 7, 7] )
+ ret = arr.locateComponentId( valToSearchIntoTuples, tupleIdHint )
+ self.assertTrue( ret.isEqual(DataArrayInt([1, 1, 0, 1, 2, 1, 0, 1, 2]) ) )
+ pass
+
+ def testMeasureOnGaussPtMeshDimNotEqualSpaceDim0(self):
+ """
+ [EDF26877] : This test focuses on computation of measure field on field on Gauss Point in the special case where SpaceDim
+ are not eqaul to the meshDim.
+ """
+ seg2 = MEDCouplingUMesh("mesh",1)
+ seg2.setCoords(DataArrayDouble([(3,3),(4,4)]))
+ seg2.allocateCells()
+ seg2.insertNextCell(NORM_SEG2,[0,1])
+ fff=MEDCouplingFieldDouble.New(ON_GAUSS_PT) ; fff.setName("CH1RB") ; fff.setNature(IntensiveMaximum)
+ fff.setMesh(seg2)
+ fff.setGaussLocalizationOnCells([0], [0.,1.], [0.333333333333333], [1.0])
+ disc = fff.getDiscretization()
+ # spaceDim = 2 meshDim = 1
+ self.assertTrue( disc.getMeasureField(seg2,True).getArray().isEqual(DataArrayDouble([sqrt(2.0)]),1e-10) )
+ # spaceDim = 3 meshDim = 1
+ seg2.setCoords(DataArrayDouble([(3,3,3),(4,4,4)]))
+ self.assertTrue( disc.getMeasureField(seg2,True).getArray().isEqual(DataArrayDouble([sqrt(3.0)]),1e-10) )
+ # spaceDim = 3 meshDim = 2
+ tri = MEDCouplingUMesh("mesh",2)
+ tri.setCoords( DataArrayDouble([(0,0,0),(1,1,0),(2,2,2)]) )
+ tri.allocateCells()
+ tri.insertNextCell(NORM_TRI3,[0,1,2])
+ fff=MEDCouplingFieldDouble.New(ON_GAUSS_PT) ; fff.setName("CH1RB") ; fff.setNature(IntensiveMaximum)
+ fff.setMesh(tri)
+ fff.setGaussLocalizationOnCells(list(range(0, 1)), [0., 0., 1., 0., 0., 1.], [0.3333333333333333, 0.3333333333333333], [0.5])
+ disc = fff.getDiscretization()
+ self.assertTrue( disc.getMeasureField(tri,True).getArray().isEqual( tri.getMeasureField(True).getArray(), 1e-10) )
+ pass
+
+ def testUMeshExplodeMeshTo(self):
+ """
+ [EDF27988] : implementation of reduceToCells implies implementation of MEDCouplingUMesh.explodeMeshTo
+ """
+ arr = DataArrayDouble(5) ; arr.iota()
+ m = MEDCouplingCMesh() ; m.setCoords(arr,arr,arr)
+ m = m.buildUnstructured()
+ m1 = m[::2] ; m2 = m[1::2]
+ m1.simplexize(PLANAR_FACE_5)
+ m = MEDCouplingUMesh.MergeUMeshesOnSameCoords([m1,m2])
+ mE1 = m.explodeMeshTo(-1)
+ ref1 = m.buildDescendingConnectivity()
+ mE2 = m.explodeMeshTo(-2)
+ ref2 = m.explode3DMeshTo1D()
+ mE3 = m.explodeMeshTo(-3)
+ self.assertTrue( len(mE1) ==5 )
+ self.assertTrue( mE1[0].getNodalConnectivity().isEqual(ref1[0].getNodalConnectivity()) )
+ self.assertTrue( mE1[0].getNodalConnectivityIndex().isEqual(ref1[0].getNodalConnectivityIndex()) )
+ self.assertTrue( mE1[0].getCoords().getHiddenCppPointer() == m.getCoords().getHiddenCppPointer() )
+ for i in range(1,5):
+ self.assertTrue( mE1[i].isEqual(ref1[i]) )
+ #
+ self.assertTrue( len(mE2) ==5 )
+ self.assertTrue( mE2[0].getNodalConnectivity().isEqual(ref2[0].getNodalConnectivity()) )
+ self.assertTrue( mE2[0].getNodalConnectivityIndex().isEqual(ref2[0].getNodalConnectivityIndex()) )
+ self.assertTrue( mE2[0].getCoords().getHiddenCppPointer() == m.getCoords().getHiddenCppPointer() )
+ for i in range(1,5):
+ self.assertTrue( mE2[i].isEqual(ref2[i]) )
+ #
+ self.assertTrue( mE3[0].getMeshDimension() == 0 )
+ self.assertTrue( mE3[0].getNumberOfCells() == mE3[0].getNumberOfNodes() )
+ a,b = m.getReverseNodalConnectivity()
+ self.assertTrue( mE3[3].isEqual(a) and mE3[4].isEqual(b) )
+ ref3_2 = (m.getNodalConnectivityIndex().deltaShiftIndex()-1) ; ref3_2.computeOffsetsFull()
+ self.assertTrue( ref3_2.isEqual(mE3[2]) )
+ tmp = m.getNodalConnectivityIndex().deepCopy() ; tmp.popBackSilent() ; tmp = tmp.buildComplement( len(m.getNodalConnectivity()) ) ; ref3_1 = m.getNodalConnectivity()[tmp]
+ self.assertTrue( ref3_1.isEqual(mE3[1]) )
+ #
+ cellsInPolyh = [37,160]
+ polyh = m[cellsInPolyh]
+ polyh.convertAllToPoly()
+ m[cellsInPolyh] = polyh
+ pE3 = m.explodeMeshTo(-3)
+ self.assertTrue( pE3[0].getMeshDimension() == 0 )
+ self.assertTrue( pE3[0].getNumberOfCells() == pE3[0].getNumberOfNodes() )
+ a,b = m.getReverseNodalConnectivity()
+ self.assertTrue( pE3[3].isEqual(a) and pE3[4].isEqual(b) )
+ self.assertTrue( pE3[2].isEqual(mE3[2]) ) # indexed arrays are the same
+
+ ref_a,ref_b = DataArrayInt.ExtractFromIndexedArrays( DataArrayInt(cellsInPolyh).buildComplement(m.getNumberOfCells()), mE3[1], mE3[2] )
+ a,b = DataArrayInt.ExtractFromIndexedArrays( DataArrayInt(cellsInPolyh).buildComplement(m.getNumberOfCells()), pE3[1], pE3[2] )
+ self.assertTrue( ref_a.isEqual(a) )
+ self.assertTrue( ref_b.isEqual(b) )
+ for cell in cellsInPolyh:
+ ref_c,ref_d = DataArrayInt.ExtractFromIndexedArrays( cell, mE3[1], mE3[2] ) ; ref_c.sort()
+ c,d = DataArrayInt.ExtractFromIndexedArrays( cell, pE3[1], pE3[2] )
+ self.assertTrue( ref_c.isEqual(c) )
+ self.assertTrue( ref_d.isEqual(d) )
+
+ def testGetCellContainingPointRelativeEps(self):
+ """
+ See EDF27860 : This test checks that detection of point inside a cell works by normalizing cell around origin with factor equal to the max delta of bbox along axis X, Y or Z.
+ """
+ # in this test cell is vuluntary far from origin {15260.775604514516, 11197.646906189217, 14187.820484060947}
+ # and caracteritic size is ~ 1500
+ coo = DataArrayDouble( [(14724.199858870656, 11928.888084722483, 14442.32726944039), (14788.407409534622, 11992.60694822231, 14453.86181555231), (15572.175148726046, 10798.586790270576, 14471.54225356788), (15643.898717334796, 10853.094666047728, 14477.233802854305), (15005.31495255754, 11573.261110174888, 13933.313698681504), (15070.29423166349, 11636.377758513776, 13946.650959030132), (15797.351350158377, 10466.40572765595, 13965.524190108257), (15869.808770928525, 10519.99285973948, 13972.419352086607), (15273.866774426142, 11216.458197414971, 13433.169979717744), (15340.421031616577, 11277.882145177837, 13446.53598386297), (16013.382514001762, 10132.795887638129, 13465.184281842226), (16086.979064572806, 10184.802292369684, 13472.147425473782)] )
+ m = MEDCouplingUMesh("",3)
+ m.setCoords(coo)
+ m.allocateCells()
+ m.insertNextCell(NORM_TETRA4,[0,5,4,6])
+ m.insertNextCell(NORM_TETRA4,[4,5,9,7])
+
+ ##### See EDF2760 pt is outside cell 0 (6e-4) and 1 (8e-4)
+ pt = DataArrayDouble([(15263.41200205526, 11314.957094727113, 13950.0)])
+ a,b = m.getCellsContainingPoints(pt,1e-3)
+ self.assertTrue(a.isEqual(DataArrayInt([0,1])))
+ self.assertTrue(b.isEqual(DataArrayInt([0,2])))
+
+ # by shifting pt by 10 along Z pt in only inside cell # 0
+ pt += [0,0,10]
+ a1,b1 = m.getCellsContainingPoints(pt,1e-3)
+ self.assertTrue(a1.isEqual(DataArrayInt([0])))
+ self.assertTrue(b1.isEqual(DataArrayInt([0,1])))
+
pass
if __name__ == '__main__':
