From 9353560dfc959c78b47fbabf6cc5c4a493f974e7 Mon Sep 17 00:00:00 2001 From: geay Date: Thu, 6 Mar 2014 17:16:16 +0100 Subject: [PATCH] bug correction concerning configuration b of PYRA13 --- .../GaussPoints/InterpKernelGaussCoords.cxx | 77 ++++++++----------- src/MEDCoupling_Swig/MEDCouplingBasicsTest.py | 25 ++++++ 2 files changed, 58 insertions(+), 44 deletions(-) diff --git a/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx b/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx index 13f90bc69..14f8af006 100644 --- a/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx +++ b/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx @@ -1372,9 +1372,9 @@ void GaussInfo::pyra13bInit() coords[1] = 0.0; coords[2] = 0.0; break; - case 3: + case 1: coords[0] = 0.0; - coords[1] = 1.0; + coords[1] = -1.0; coords[2] = 0.0; break; case 2: @@ -1382,9 +1382,9 @@ void GaussInfo::pyra13bInit() coords[1] = 0.0; coords[2] = 0.0; break; - case 1: + case 3: coords[0] = 0.0; - coords[1] = -1.0; + coords[1] = 1.0; coords[2] = 0.0; break; case 4: @@ -1392,24 +1392,24 @@ void GaussInfo::pyra13bInit() coords[1] = 0.0; coords[2] = 1.0; break; - case 8: + case 5: coords[0] = 0.5; - coords[1] = 0.5; + coords[1] = -0.5; coords[2] = 0.0; break; - case 7: + case 6: coords[0] = -0.5; - coords[1] = 0.5; + coords[1] = -0.5; coords[2] = 0.0; break; - case 6: + case 7: coords[0] = -0.5; - coords[1] = -0.5; + coords[1] = 0.5; coords[2] = 0.0; break; - case 5: + case 8: coords[0] = 0.5; - coords[1] = -0.5; + coords[1] = 0.5; coords[2] = 0.0; break; case 9: @@ -1417,9 +1417,9 @@ void GaussInfo::pyra13bInit() coords[1] = 0.0; coords[2] = 0.5; break; - case 12: + case 10: coords[0] = 0.0; - coords[1] = 0.5; + coords[1] = -0.5; coords[2] = 0.5; break; case 11: @@ -1427,42 +1427,31 @@ void GaussInfo::pyra13bInit() coords[1] = 0.0; coords[2] = 0.5; break; - case 10: + case 12: coords[0] = 0.0; - coords[1] = -0.5; + coords[1] = 0.5; coords[2] = 0.5; break; LOCAL_COORD_MACRO_END; SHAPE_FUN_MACRO_BEGIN; - funValue[0] = 0.5*(-gc[0] + gc[1] + gc[2] - 1.0)*(-gc[0] - gc[1] + gc[2] - 1.0)* - (gc[0] - 0.5)/(1.0 - gc[2]); - funValue[3] = 0.5*(-gc[0] - gc[1] + gc[2] - 1.0)*(+gc[0] - gc[1] + gc[2] - 1.0)* - (gc[1] - 0.5)/(1.0 - gc[2]); - funValue[2] = 0.5*(+gc[0] - gc[1] + gc[2] - 1.0)*(+gc[0] + gc[1] + gc[2] - 1.0)* - (-gc[0] - 0.5)/(1.0 - gc[2]); - funValue[1] = 0.5*(+gc[0] + gc[1] + gc[2] - 1.0)*(-gc[0] + gc[1] + gc[2] - 1.0)* - (-gc[1] - 0.5)/(1.0 - gc[2]); - - funValue[4] = 2.0*gc[2]*(gc[2] - 0.5); - - funValue[8] = 0.5*(-gc[0] + gc[1] + gc[2] - 1.0)*(-gc[0] - gc[1] + gc[2] - 1.0)* - (gc[0] - gc[1] + gc[2] - 1.0)/(1.0 - gc[2]); - funValue[7] = 0.5*(-gc[0] - gc[1] + gc[2] - 1.0)*(gc[0] - gc[1] + gc[2] - 1.0)* - (gc[0] + gc[1] + gc[2] - 1.0)/(1.0 - gc[2]); - funValue[6] = 0.5*(gc[0] - gc[1] + gc[2] - 1.0)*(gc[0] + gc[1] + gc[2] - 1.0)* - (-gc[0] + gc[1] + gc[2] - 1.0)/(1.0 - gc[2]); - funValue[5] = 0.5*(gc[0] + gc[1] + gc[2] - 1.0)*(-gc[0] + gc[1] + gc[2] - 1.0)* - (-gc[0] - gc[1] + gc[2] - 1.0)/(1.0 - gc[2]); - - funValue[9] = 0.5*gc[2]*(-gc[0] + gc[1] + gc[2] - 1.0)*(-gc[0] - gc[1] + gc[2] - 1.0)/ - (1.0 - gc[2]); - funValue[12] = 0.5*gc[2]*(-gc[0] - gc[1] + gc[2] - 1.0)*(gc[0] - gc[1] + gc[2] - 1.0)/ - (1.0 - gc[2]); - funValue[11] = 0.5*gc[2]*(gc[0] - gc[1] + gc[2] - 1.0)*(gc[0] + gc[1] + gc[2] - 1.0)/ - (1.0 - gc[2]); - funValue[10] = 0.5*gc[2]*(gc[0] + gc[1] + gc[2] - 1.0)*(-gc[0] + gc[1] + gc[2] - 1.0)/ - (1.0 - gc[2]); + funValue[0] =0.5*(-gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]-gc[1]+gc[2]-1.0)*(gc[0]-0.5)/(1.0-gc[2]); + funValue[1] =0.5*(+gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]+gc[1]+gc[2]-1.0)*(-gc[1]-0.5)/(1.0-gc[2]); + funValue[2] =0.5*(+gc[0]-gc[1]+gc[2]-1.0)*(+gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]-0.5)/(1.0-gc[2]); + funValue[3] =0.5*(-gc[0]-gc[1]+gc[2]-1.0)*(+gc[0]-gc[1]+gc[2]-1.0)*(gc[1]-0.5)/(1.0-gc[2]); + + funValue[4] =2.*gc[2]*(gc[2]-0.5); + + funValue[5] =-0.5*(gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]-gc[1]+gc[2]-1.0)/(1.0-gc[2]); + funValue[6] =-0.5*(gc[0]-gc[1]+gc[2]-1.0)*(gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]+gc[1]+gc[2]-1.0)/(1.0-gc[2]); + funValue[7] =-0.5*(-gc[0]-gc[1]+gc[2]-1.0)*(gc[0]-gc[1]+gc[2]-1.0)*(gc[0]+gc[1]+gc[2]-1.0)/(1.0-gc[2]); + funValue[8] =-0.5*(-gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]-gc[1]+gc[2]-1.0)*(gc[0]-gc[1]+gc[2]-1.0)/(1.0-gc[2]); + + funValue[9] =gc[2]*(-gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]-gc[1]+gc[2]-1.0)/(1.0-gc[2]); + funValue[10]=gc[2]*(gc[0]+gc[1]+gc[2]-1.0)*(-gc[0]+gc[1]+gc[2]-1.0)/(1.0-gc[2]); + funValue[11]=gc[2]*(gc[0]-gc[1]+gc[2]-1.0)*(gc[0]+gc[1]+gc[2]-1.0)/(1.0-gc[2]); + funValue[12]=gc[2]*(-gc[0]-gc[1]+gc[2]-1.0)*(gc[0]-gc[1]+gc[2]-1.0)/(1.0-gc[2]); + SHAPE_FUN_MACRO_END; } diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py index 6bc2ddd4f..de8113d3f 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py @@ -14466,6 +14466,31 @@ class MEDCouplingBasicsTest(unittest.TestCase): self.assertTrue(arrOfDisc2.isEqual(coo,1e-12)) pass + def testSwig2Pyra13GP1(self): + coo=DataArrayDouble([[0.,2.,0.],[2.,2.,0.],[2.,0.,0.],[0.,0.,0.],[1.,1.,2.],[1.,2.,0.],[2.,1.,0.],[1.,0.,0.],[0.,1.,0.],[0.5,1.5,1.],[1.5,1.5,1.],[1.5,0.5,1.],[0.5,0.5,1.]]) + m=MEDCouplingUMesh("mesh",3) ; m.setCoords(coo) + m.allocateCells() + # the cell description is exactly those described in the description of PYRA13 in MED file 3.0.7 documentation + m.insertNextCell(NORM_PYRA13,[0,1,2,3,4,5,6,7,8,9,10,11,12]) + refCoords=[1.,0.,0.,0.,-1.,0.,-1.,0.,0.,0.,1.,0.,0.,0.,1.,0.5,-0.5,0.,-0.5,-0.5,0.,-0.5,0.5,0.,0.5,0.5,0.,0.5,0.,0.5,0.,-0.5,0.5,-0.5,0.,0.5,0.,0.5,0.5] + gaussCoords=[0.,0.,0.5,0.21210450275,0.21210450275,0.5,-0.21210450275,0.21210450275,0.5,-0.21210450275,-0.21210450275,0.5,0.21210450275,-0.21210450275,0.5,0.,0.,0.07579099449999999,0.,0.,0.9242090055000001,0.5394929090572634,0.,0.17359176399999998,0.,0.5394929090572634,0.17359176399999998,-0.5394929090572634,0.,0.17359176399999998,0.,-0.5394929090572634,0.17359176399999998,0.1133235629427366,0.,0.826408236,0.,0.1133235629427366,0.826408236,-0.1133235629427366,0.,0.826408236,0.,-0.1133235629427366,0.826408236,0.5826406005183961,0.5826406005183961,-0.053206449499999975,-0.5826406005183961,0.5826406005183961,-0.053206449499999975,-0.5826406005183961,-0.5826406005183961,-0.053206449499999975,0.5826406005183961,-0.5826406005183961,-0.053206449499999975,0.5532064495,0.,0.5,0.,0.5532064495,0.5,-0.5532064495,0.,0.5,0.,-0.5532064495,0.5,-0.029434151018396033,-0.029434151018396033,1.0532064495,0.029434151018396033,-0.029434151018396033,1.0532064495,0.029434151018396033,0.029434151018396033,1.0532064495,-0.029434151018396033,0.029434151018396033,1.0532064495] + weights=[0.0492545926875,0.031210562625,0.031210562625,0.031210562625,0.031210562625,0.10663554205740113,0.0007171281994273535,0.0816994048010844,0.0816994048010844,0.0816994048010844,0.0816994048010844,0.0036048554264914074,0.0036048554264914074,0.0036048554264914074,0.0036048554264914074,0.008958181586640837,0.008958181586640837,0.008958181586640837,0.008958181586640837,0.002018983875,0.002018983875,0.002018983875,0.002018983875,2.286237794882217e-05,2.286237794882217e-05,2.286237794882217e-05,2.286237794882217e-05] + fGauss=MEDCouplingFieldDouble(ON_GAUSS_PT) ; fGauss.setName("fGauss") + fGauss.setMesh(m) + fGauss.setGaussLocalizationOnType(NORM_PYRA13,refCoords,gaussCoords,weights) + arr=DataArrayDouble(fGauss.getNumberOfTuplesExpected()) ; arr.iota() + fGauss.setArray(arr) + arrOfDisc=fGauss.getLocalizationOfDiscr() + # the test is here + self.assertTrue(arrOfDisc.isEqual(DataArrayDouble([1.,1.,1.,0.5757909945,1.,1.,1.,0.5757909945,1.,1.4242090055,1.,1.,1.,1.4242090055,1.,1.,1.,0.151581989,1.,1.,1.848418011,0.4605070909427367,1.5394929090572635,0.347183528,0.4605070909427367,0.4605070909427367,0.347183528,1.5394929090572638,0.4605070909427366,0.347183528,1.5394929090572635,1.5394929090572638,0.347183528,0.8866764370572636,1.1133235629427367,1.652816472,0.8866764370572636,0.8866764370572636,1.652816472,1.1133235629427367,0.8866764370572636,1.652816472,1.1133235629427365,1.1133235629427367,1.652816472,-0.16528120103679209,1.,-0.106412899,1.,-0.1652812010367921,-0.106412899,2.1652812010367914,1.,-0.106412899,1.,2.165281201036791,-0.106412899,0.4467935505,1.5532064495,1.,0.4467935505,0.4467935505,1.,1.5532064495,0.4467935505,1.,1.5532064495,1.5532064495,1.,1.0588683020367922,1.,2.106412899,1.,1.0588683020367922,2.106412899,0.9411316979632077,1.,2.106412899,1.,0.9411316979632078,2.106412899],27,3),1e-12)) + # + weights=13*[1] + gaussCoords=refCoords[:] ; gaussCoords[14]=0.9999999999999 # change z of point #4 0.999... instead of 1. because with shape function it leads to division by 0. ! + fGauss.setGaussLocalizationOnType(NORM_PYRA13,refCoords,gaussCoords,weights) + arrOfDisc2=fGauss.getLocalizationOfDiscr() + self.assertTrue(arrOfDisc2.isEqual(coo,1e-10)) # be less exigent 1e-10 instead of 1e-12 due to shape function sensitivity arount 0.,0.,1. ! + pass + def setUp(self): pass pass -- 2.39.2