# -*- coding: iso-8859-1 -*-
-# Copyright (C) 2007-2013 CEA/DEN, EDF R&D
+# Copyright (C) 2007-2014 CEA/DEN, EDF R&D
#
# 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.
+# 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
self.assertAlmostEqual(0.3521445110626687 ,m[6][170],12)
pass
- @unittest.skipUnless(MEDCouplingHasNumPyBindings(),"requires numpy")
+ @unittest.skipUnless(MEDCouplingHasNumPyBindings() and MEDCouplingHasSciPyBindings(),"requires numpy AND scipy")
def testGetCrudeCSRMatrix1(self):
""" testing CSR matrix output using numpy/scipy.
"""
- from scipy.sparse import diags
+ from scipy.sparse import spdiags #diags
import scipy
from numpy import array
arr=DataArrayDouble(3) ; arr.iota()
pass
m2=rem.getCrudeCSRMatrix()
diff=m-m2
- assert(isinstance(m,scipy.sparse.csr.csr_matrix))
- assert(m.getnnz()==7)
+ self.assertTrue(isinstance(m,scipy.sparse.csr.csr_matrix))
+ self.assertEqual(m.getnnz(),7)
self.assertAlmostEqual(m[0,0],0.25,12)
self.assertAlmostEqual(m[1,0],0.1,12)
self.assertAlmostEqual(m[1,1],0.15,12)
self.assertEqual(diff.getnnz(),0)
# IntegralGlobConstraint (division by sum of cols)
colSum=m.sum(axis=0)
- m_0=m*diags(array(1/colSum),[0])
+ # version 0.12.0 # m_0=m*diags(array(1/colSum),[0])
+ m_0=m*spdiags(array(1/colSum),[0],colSum.shape[1],colSum.shape[1])
del colSum
self.assertAlmostEqual(m_0[0,0],0.625,12)
self.assertAlmostEqual(m_0[1,0],0.25,12)
self.assertAlmostEqual(m_0[2,1],0.3333333333333333,12)
self.assertAlmostEqual(m_0[2,2],1.,12)
self.assertAlmostEqual(m_0[2,3],1.,12)
- assert(m_0.getnnz()==7)
+ self.assertEqual(m_0.getnnz(),7)
# ConservativeVolumic (division by sum of rows)
rowSum=m.sum(axis=1)
- m_1=diags(array(1/rowSum.transpose()),[0])*m
+ # version 0.12.0 # m_1=diags(array(1/rowSum.transpose()),[0])*m
+ m_1=spdiags(array(1/rowSum.transpose()),[0],rowSum.shape[0],rowSum.shape[0])*m
del rowSum
self.assertAlmostEqual(m_1[0,0],1.,12)
self.assertAlmostEqual(m_1[1,0],0.4,12)
self.assertAlmostEqual(m_1[2,1],0.3,12)
self.assertAlmostEqual(m_1[2,2],0.2,12)
self.assertAlmostEqual(m_1[2,3],0.3,12)
- assert(m_1.getnnz()==7)
+ self.assertEqual(m_1.getnnz(),7)
+ pass
+
+ @unittest.skipUnless(MEDCouplingHasNumPyBindings() and MEDCouplingHasSciPyBindings(),"requires numpy AND scipy")
+ def testP0P1Bary_1(self):
+ a=MEDCouplingUMesh("a",2)
+ a.allocateCells()
+ conna=[0,1,3,2,1,4,5,3,4,6,7,5,6,8,9,7,8,10,11,9,10,12,13,11,12,14,15,13,14,16,17,15,16,18,19,17,18,20,21,19,20,22,23,21,22,24,25,23,24,26,27,25]
+ a.setCoords(DataArrayDouble([1.54,0,-0.01,1.54,0.02,-0.01,1.54,0,0.01,1.54,0.02,0.01,1.54,0.04,-0.01,1.54,0.04,0.01,1.54,0.06,-0.01,1.54,0.06,0.01,1.54,0.08,-0.01,1.54,0.08,0.01,1.54,0.1,-0.01,1.54,0.1,0.01,1.54,0.12,-0.01,1.54,0.12,0.01,1.54,0.14,-0.01,1.54,0.14,0.01,1.54,0.16,-0.01,1.54,0.16,0.01,1.54,0.18,-0.01,1.54,0.18,0.01,1.54,0.2,-0.01,1.54,0.2,0.01,1.54,0.22,-0.01,1.54,0.22,0.01,1.54,0.24,-0.01,1.54,0.24,0.01,1.54,0.26,-0.01,1.54,0.26,0.01],28,3))
+ for i in xrange(13):
+ a.insertNextCell(NORM_QUAD4,conna[4*i:4*(i+1)])
+ pass
+ a.finishInsertingCells() ; a.simplexize(0)
+ #
+ connb=[0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,0,2,39,3,5,40,6,8,41,9,11,42,12,14,43,15,17,44,18,20,45,21,23,46,24,26,47,27,29,48,30,32,49,33,35,50,36,38,51,52,2,39,53,5,40,54,8,41,55,11,42,56,14,43,57,17,44,58,20,45,59,23,46,60,26,47,61,29,48,62,32,49,63,35,50,64,38,51,52,2,65,53,5,66,54,8,67,55,11,68,56,14,69,57,17,70,58,20,71,59,23,72,60,26,73,61,29,74,62,32,75,63,35,76,64,38,77,53,2,65,54,5,66,55,8,67,56,11,68,57,14,69,58,17,70,59,20,71,60,23,72,61,26,73,62,29,74,63,32,75,64,35,76,78,38,77,53,2,40,54,5,41,55,8,42,56,11,43,57,14,44,58,17,45,59,20,46,60,23,47,61,26,48,62,29,49,63,32,50,64,35,51,78,38,79,3,2,40,6,5,41,9,8,42,12,11,43,15,14,44,18,17,45,21,20,46,24,23,47,27,26,48,30,29,49,33,32,50,36,35,51,80,38,79,3,2,1,6,5,4,9,8,7,12,11,10,15,14,13,18,17,16,21,20,19,24,23,22,27,26,25,30,29,28,33,32,31,36,35,34,80,38,37]
+ b=MEDCouplingUMesh("b",2)
+ b.allocateCells()
+ for i in xrange(104):
+ b.insertNextCell(NORM_TRI3,connb[3*i:3*(i+1)])
+ pass
+ b.setCoords(DataArrayDouble([1.54,0,-0.01,1.54,0.01,-0.01,1.54,0.01,0,1.54,0.02,-0.01,1.54,0.03,-0.01,1.54,0.03,0,1.54,0.04,-0.01,1.54,0.05,-0.01,1.54,0.05,0,1.54,0.06,-0.01,1.54,0.07,-0.01,1.54,0.07,0,1.54,0.08,-0.01,1.54,0.09,-0.01,1.54,0.09,0,1.54,0.1,-0.01,1.54,0.11,-0.01,1.54,0.11,0,1.54,0.12,-0.01,1.54,0.13,-0.01,1.54,0.13,0,1.54,0.14,-0.01,1.54,0.15,-0.01,1.54,0.15,0,1.54,0.16,-0.01,1.54,0.17,-0.01,1.54,0.17,0,1.54,0.18,-0.01,1.54,0.19,-0.01,1.54,0.19,0,1.54,0.2,-0.01,1.54,0.21,-0.01,1.54,0.21,0,1.54,0.22,-0.01,1.54,0.23,-0.01,1.54,0.23,0,1.54,0.24,-0.01,1.54,0.25,-0.01,1.54,0.25,0,1.54,0,0,1.54,0.02,0,1.54,0.04,0,1.54,0.06,0,1.54,0.08,0,1.54,0.1,0,1.54,0.12,0,1.54,0.14,0,1.54,0.16,0,1.54,0.18,0,1.54,0.2,0,1.54,0.22,0,1.54,0.24,0,1.54,0,0.01,1.54,0.02,0.01,1.54,0.04,0.01,1.54,0.06,0.01,1.54,0.08,0.01,1.54,0.1,0.01,1.54,0.12,0.01,1.54,0.14,0.01,1.54,0.16,0.01,1.54,0.18,0.01,1.54,0.2,0.01,1.54,0.22,0.01,1.54,0.24,0.01,1.54,0.01,0.01,1.54,0.03,0.01,1.54,0.05,0.01,1.54,0.07,0.01,1.54,0.09,0.01,1.54,0.11,0.01,1.54,0.13,0.01,1.54,0.15,0.01,1.54,0.17,0.01,1.54,0.19,0.01,1.54,0.21,0.01,1.54,0.23,0.01,1.54,0.25,0.01,1.54,0.26,0.01,1.54,0.26,0,1.54,0.26,-0.01],81,3))
+ #
+ rem=MEDCouplingRemapper() ; rem.setIntersectionType(Barycentric)
+ rem.prepare(a,b,"P1P0")
+ m0=rem.getCrudeCSRMatrix()
+ self.assertEqual(m0.nnz,312)
+ #
+ ids=4*[None] ; vs=4*[None]
+ ids[0]=DataArrayInt([0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48,51,54,57,60,63,66,69,72,75,158,161,164,167,170,173,176,179,182,185,188,191,194,197,200,203,206,209,212,215,218,221,224,227,230,233]) ; vs[0]=10./3.
+ ids[1]=DataArrayInt([1,2,4,5,7,8,10,11,13,14,16,17,19,20,22,23,25,26,28,29,31,32,34,35,37,38,40,41,43,44,46,47,49,50,52,53,55,56,58,59,61,62,64,65,67,68,70,71,73,74,76,77,80,83,86,89,92,95,98,101,104,107,110,113,116,117,120,123,126,129,132,135,138,141,144,147,150,153,156,157,159,160,162,163,165,166,168,169,171,172,174,175,177,178,180,181,183,184,186,187,189,190,192,193,195,196,198,199,201,202,204,205,207,208,210,211,213,214,216,217,219,220,222,223,225,226,228,229,231,232,234,237,240,243,246,249,252,255,258,261,264,267,270,275,278,281,284,287,290,293,296,299,302,305,308,311]) ; vs[1]=5./6.
+ ids[2]=DataArrayInt([78,81,84,87,90,93,96,99,102,105,108,111,114,119,122,125,128,131,134,137,140,143,146,149,152,155,236,239,242,245,248,251,254,257,260,263,266,269,272,273,276,279,282,285,288,291,294,297,300,303,306,309]) ; vs[2]=5./3.
+ ids[3]=DataArrayInt([79,82,85,88,91,94,97,100,103,106,109,112,115,118,121,124,127,130,133,136,139,142,145,148,151,154,235,238,241,244,247,250,253,256,259,262,265,268,271,274,277,280,283,286,289,292,295,298,301,304,307,310]) ; vs[3]=2.5
+ vals=DataArrayDouble(312,1)
+ for idd,v in zip(ids,vs):
+ vals[idd]=v
+ pass
+ vals*=1e-5
+ eps0=DataArrayDouble(m0.data)-vals ; eps0.abs()
+ self.assertTrue(eps0.getIdsInRange(1e-17,1e300).empty())
+ self.assertTrue(DataArrayInt(m0.indices).isEqual(DataArrayInt([0,1,3,1,4,5,4,6,7,6,8,9,8,10,11,10,12,13,12,14,15,14,16,17,16,18,19,18,20,21,20,22,23,22,24,25,24,26,27,0,2,3,1,3,5,4,5,7,6,7,9,8,9,11,10,11,13,12,13,15,14,15,17,16,17,19,18,19,21,20,21,23,22,23,25,24,25,27,0,2,3,1,3,5,4,5,7,6,7,9,8,9,11,10,11,13,12,13,15,14,15,17,16,17,19,18,19,21,20,21,23,22,23,25,24,25,27,0,2,3,1,3,5,4,5,7,6,7,9,8,9,11,10,11,13,12,13,15,14,15,17,16,17,19,18,19,21,20,21,23,22,23,25,24,25,27,0,2,3,1,3,5,4,5,7,6,7,9,8,9,11,10,11,13,12,13,15,14,15,17,16,17,19,18,19,21,20,21,23,22,23,25,24,25,27,0,1,3,1,4,5,4,6,7,6,8,9,8,10,11,10,12,13,12,14,15,14,16,17,16,18,19,18,20,21,20,22,23,22,24,25,24,26,27,0,1,3,1,4,5,4,6,7,6,8,9,8,10,11,10,12,13,12,14,15,14,16,17,16,18,19,18,20,21,20,22,23,22,24,25,24,26,27,0,1,3,1,4,5,4,6,7,6,8,9,8,10,11,10,12,13,12,14,15,14,16,17,16,18,19,18,20,21,20,22,23,22,24,25,24,26,27])))
+ self.assertTrue(DataArrayInt(m0.indptr).isEqual(DataArrayInt([0,3,6,9,12,15,18,21,24,27,30,33,36,39,42,45,48,51,54,57,60,63,66,69,72,75,78,81,84,87,90,93,96,99,102,105,108,111,114,117,120,123,126,129,132,135,138,141,144,147,150,153,156,159,162,165,168,171,174,177,180,183,186,189,192,195,198,201,204,207,210,213,216,219,222,225,228,231,234,237,240,243,246,249,252,255,258,261,264,267,270,273,276,279,282,285,288,291,294,297,300,303,306,309,312])))
+ #
+ rem2=MEDCouplingRemapper() ; rem2.setIntersectionType(Barycentric)
+ rem2.prepare(b,a,"P0P1")
+ m1=rem2.getCrudeCSRMatrix()
+ self.assertEqual(m1.nnz,312)
+ #
+ m1=rem2.getCrudeCSRMatrix()
+ m1t=m1.transpose()
+ delta=m0-m1t
+ self.assertTrue(DataArrayDouble(delta.data).isUniform(0.,1e-17))
+ pass
+
+ @unittest.skipUnless(MEDCouplingHasNumPyBindings() and MEDCouplingHasSciPyBindings(),"requires numpy AND scipy")
+ def testNonConformWithRemapper_1(self):
+ coo=DataArrayDouble([-0.396700000780411,-0.134843245350081,-0.0361311386958691,-0.407550009429364,-0.13484324535008,-0.0361311386958923,-0.396700000780411,-0.132191446077668,-0.0448729493559049,-0.407550009429364,-0.132191446077666,-0.0448729493559254,-0.396700000780411,-0.128973582738749,-0.0534226071577727,-0.407550009429364,-0.128973582738747,-0.0534226071577904,-0.396700000780411,-0.128348829636458,-0.0346583696473619,-0.407550009429364,-0.128348829636457,-0.0346583696473822,-0.396700000780411,-0.125874740261886,-0.0430683597970123,-0.407550009429364,-0.125874740261885,-0.0430683597970302,-0.396700000780411,-0.122905344829122,-0.051310216195766,-0.407550009429364,-0.12290534482912,-0.0513102161957814],12,3)
+ conn=DataArrayInt([2,9,3,11,2,3,5,11,2,8,9,11,2,10,8,11,2,5,4,11,2,4,10,11,3,0,1,6,3,1,7,6,3,2,0,6,3,8,2,6,3,7,9,6,3,9,8,6])
+ m=MEDCoupling1SGTUMesh("mesh",NORM_TETRA4)
+ m.setNodalConnectivity(conn)
+ m.setCoords(coo)
+ # m is ready
+ m1,d,di,rd,rdi=m.buildUnstructured().buildDescendingConnectivity()
+ rdi2=rdi.deltaShiftIndex()
+ cellIds=rdi2.getIdsEqual(1)
+ skinAndNonConformCells=m1[cellIds]
+ skinAndNonConformCells.zipCoords() # at this point skinAndNonConformCells contains non conform cells and skin cells. Now trying to split them in two parts.
+ #
+ rem=MEDCouplingRemapper()
+ rem.setMaxDistance3DSurfIntersect(1e-12)
+ rem.setMinDotBtwPlane3DSurfIntersect(0.99)# this line is important it is to tell to remapper to select only cells with very close orientation
+ rem.prepare(skinAndNonConformCells,skinAndNonConformCells,"P0P0")
+ mat=rem.getCrudeCSRMatrix()
+ indptr=DataArrayInt(mat.indptr)
+ indptr2=indptr.deltaShiftIndex()
+ cellIdsOfNonConformCells=indptr2.getIdsNotEqual(1)
+ cellIdsOfSkin=indptr2.getIdsEqual(1)
+ self.assertTrue(cellIdsOfSkin.isEqual(DataArrayInt([1,2,3,5,6,7,8,9,10,11,12,13,14,15,16,17,19,20,21,23])))
+ self.assertTrue(cellIdsOfNonConformCells.isEqual(DataArrayInt([0,4,18,22])))
+ pass
+
+ def test3D1DOnP1P0_1(self):
+ """ This test focused on P1P0 interpolation with a source with meshDim=1 spaceDim=3 and a target with meshDim=3.
+ This test has revealed a bug in remapper. A reverse matrix is computed so a reverse method should be given in input.
+ """
+ target=MEDCouplingCMesh()
+ arrX=DataArrayDouble([0,1]) ; arrY=DataArrayDouble([0,1]) ; arrZ=DataArrayDouble(11) ; arrZ.iota()
+ target.setCoords(arrX,arrY,arrZ)
+ target=target.buildUnstructured() ; target.setName("TargetSecondaire")
+ #
+ sourceCoo=DataArrayDouble([(0.5,0.5,0.1),(0.5,0.5,1.2),(0.5,0.5,1.6),(0.5,0.5,1.8),(0.5,0.5,2.43),(0.5,0.5,2.55),(0.5,0.5,4.1),(0.5,0.5,4.4),(0.5,0.5,4.9),(0.5,0.5,5.1),(0.5,0.5,7.6),(0.5,0.5,7.7),(0.5,0.5,8.2),(0.5,0.5,8.4),(0.5,0.5,8.6),(0.5,0.5,8.8),(0.5,0.5,9.2),(0.5,0.5,9.6),(0.5,0.5,11.5)])
+ source=MEDCoupling1SGTUMesh("SourcePrimaire",NORM_SEG2)
+ source.setCoords(sourceCoo)
+ source.allocateCells()
+ for i in xrange(len(sourceCoo)-1):
+ source.insertNextCell([i,i+1])
+ pass
+ source=source.buildUnstructured()
+ fsource=MEDCouplingFieldDouble(ON_NODES) ; fsource.setName("field")
+ fsource.setMesh(source)
+ arr=DataArrayDouble(len(sourceCoo)) ; arr.iota(0.7) ; arr*=arr
+ fsource.setArray(arr)
+ fsource.setNature(ConservativeVolumic)
+ #
+ rem=MEDCouplingRemapper()
+ rem.setIntersectionType(PointLocator)
+ rem.prepare(source,target,"P1P0")
+ f2Test=rem.transferField(fsource,-27)
+ self.assertEqual(f2Test.getName(),fsource.getName())
+ self.assertEqual(f2Test.getMesh().getHiddenCppPointer(),target.getHiddenCppPointer())
+ expArr=DataArrayDouble([0.49,7.956666666666667,27.29,-27,59.95666666666667,94.09,-27,125.69,202.89,296.09])
+ self.assertTrue(f2Test.getArray().isEqual(expArr,1e-12))
+ f2Test=rem.reverseTransferField(f2Test,-36)
+ self.assertEqual(f2Test.getName(),fsource.getName())
+ self.assertEqual(f2Test.getMesh().getHiddenCppPointer(),source.getHiddenCppPointer())
+ expArr2=DataArrayDouble([0.49,7.956666666666667,7.956666666666667,7.956666666666667,27.29,27.29,59.95666666666667,59.95666666666667,59.95666666666667,94.09,125.69,125.69,202.89,202.89,202.89,202.89,296.09,296.09,-36.])
+ self.assertTrue(f2Test.getArray().isEqual(expArr2,1e-12))
+ pass
+
+ def testRemapperAMR1(self):
+ """ This test is the origin of the ref values for MEDCouplingBasicsTest.testAMR2"""
+ coarse=DataArrayDouble(35) ; coarse.iota(0) #X=5,Y=7
+ fine=DataArrayDouble(3*2*4*4) ; fine.iota(0) #X=3,Y=2 refined by 4
+ MEDCouplingIMesh.CondenseFineToCoarse([5,7],fine,[(1,4),(2,4)],[4,4],coarse)
+ #
+ m=MEDCouplingCartesianAMRMesh("mesh",2,[6,8],[0.,0.],[1.,1.])
+ trgMesh=m.buildUnstructured()
+ m.addPatch([(1,4),(2,4)],[4,4])
+ srcMesh=m[0].getMesh().buildUnstructured()
+ srcField=MEDCouplingFieldDouble(ON_CELLS)
+ fine2=DataArrayDouble(3*2*4*4) ; fine2.iota(0) ; srcField.setArray(fine2)
+ srcField.setMesh(srcMesh) ; srcField.setNature(Integral)
+ #
+ trgField=MEDCouplingFieldDouble(ON_CELLS)
+ coarse2=DataArrayDouble(35) ; coarse2.iota(0) ; trgField.setArray(coarse2)
+ trgField.setMesh(trgMesh) ; trgField.setNature(Integral)
+ #
+ rem=MEDCouplingRemapper()
+ rem.prepare(srcMesh,trgMesh,"P0P0")
+ rem.partialTransfer(srcField,trgField)
+ #
+ self.assertTrue(coarse.isEqual(trgField.getArray(),1e-12))
pass
def build2DSourceMesh_1(self):
