Salome HOME
Make joints easy to use with python
[tools/medcoupling.git] / src / MEDLoader / Swig / SauvLoaderTest.py
1 #  -*- coding: iso-8859-1 -*-
2 # Copyright (C) 2007-2016  CEA/DEN, EDF R&D
3 #
4 # This library is free software; you can redistribute it and/or
5 # modify it under the terms of the GNU Lesser General Public
6 # License as published by the Free Software Foundation; either
7 # version 2.1 of the License, or (at your option) any later version.
8 #
9 # This library is distributed in the hope that it will be useful,
10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
12 # Lesser General Public License for more details.
13 #
14 # You should have received a copy of the GNU Lesser General Public
15 # License along with this library; if not, write to the Free Software
16 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
17 #
18 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 #
20 # Author : Edward AGAPOV (eap)
21
22 from MEDLoader import *
23 import unittest, os
24 from MEDLoaderDataForTest import MEDLoaderDataForTest
25
26 class SauvLoaderTest(unittest.TestCase):
27
28     def __getResourcesDirectory(self):
29         med_root_dir = os.getenv("MEDCOUPLING_ROOT_DIR")
30         if med_root_dir:
31             pth = os.path.join( os.getenv("MEDCOUPLING_ROOT_DIR"), "share","resources","med")
32             if os.path.exists(pth):
33               return pth
34         current_dir = os.path.dirname(os.path.realpath(__file__))
35         pth = os.path.join(current_dir, "..", "..", "..", "resources")
36         if not os.path.exists(pth):
37           raise Exception("SauvLoaderTest: Unable to get resource directory")
38         return pth
39         pass
40
41     def testSauv2Med(self):
42         # get a file containing all types of readable piles
43         sauvFile = os.path.join( self.__getResourcesDirectory(),"allPillesTest.sauv")
44         self.assertTrue( os.access( sauvFile, os.F_OK))
45
46         # read SAUV and write MED
47         medFile = "SauvLoaderTest.med"
48         sr=SauvReader(sauvFile);
49         d2=sr.loadInMEDFileDS();
50         d2.write(medFile,0);
51
52         # check
53         self.assertEqual(1,d2.getNumberOfMeshes())
54         self.assertEqual(8+97,d2.getNumberOfFields())
55         mm = d2.getMeshes()
56         m = mm.getMeshAtPos(0)
57         self.assertEqual(17,len(m.getGroupsNames()))
58
59         os.remove( medFile )
60         pass
61
62     def testMed2Sauv(self):
63         # read pointe.med
64         medFile = os.path.join(self.__getResourcesDirectory(),"pointe.med")
65         self.assertTrue( os.access( medFile, os.F_OK))
66         pointeMed = MEDFileData.New( medFile )
67
68         # add 3 faces to pointeMed
69         pointeMedMesh = pointeMed.getMeshes().getMeshAtPos(0)
70         pointeM1D = MEDCouplingUMesh.New()
71         pointeM1D.setCoords( pointeMedMesh.getCoords() )
72         pointeM1D.setMeshDimension( 2 )
73         pointeM1D.allocateCells( 3 )
74         pointeM1D.insertNextCell( NORM_TRI3, 3, [0,1,2])
75         pointeM1D.insertNextCell( NORM_TRI3, 3, [0,1,3])
76         pointeM1D.insertNextCell( NORM_QUAD4, 4, [10,11,12,13])
77         pointeM1D.finishInsertingCells()
78         pointeMedMesh.setMeshAtLevel( -1, pointeM1D )
79         pointeMed.getMeshes().setMeshAtPos( 0, pointeMedMesh )
80
81         # add a field on 2 faces to pointeMed
82         ff1=MEDFileFieldMultiTS.New()
83         f1=MEDCouplingFieldDouble.New(ON_GAUSS_NE,ONE_TIME)
84         #f1.setMesh( pointeM1D )
85         f1.setName("Field on 2 faces")
86         d=DataArrayDouble.New()
87         d.alloc(3+4,2)
88         d.setInfoOnComponent(0,"sigX [MPa]")
89         d.setInfoOnComponent(1,"sigY [GPa]")
90         d.setValues([311,312,321,322,331,332,411,412,421,422,431,432,441,442],3+4,2)
91         f1.setArray(d)
92         da=DataArrayInt.New()
93         da.setValues([0,2],2,1)
94         da.setName("sup2")
95         ff1.appendFieldProfile(f1,pointeMedMesh,-1,da)
96         pointeMed.getFields().pushField( ff1 )
97
98         #remove fieldnodeint
99         pointeFields = pointeMed.getFields()
100         for i in range( pointeFields.getNumberOfFields() ):
101             if pointeFields.getFieldAtPos(i).getName() == "fieldnodeint":
102                 pointeFields.destroyFieldAtPos( i )
103                 break
104
105         # write pointeMed to SAUV
106         sauvFile = "SauvLoaderTest.sauv"
107         sw=SauvWriter();
108         sw.setMEDFileDS(pointeMed);
109         sw.write(sauvFile);
110
111         # read SAUV and check
112         sr=SauvReader.New(sauvFile);
113         d2=sr.loadInMEDFileDS();
114         self.assertEqual(1,d2.getNumberOfMeshes())
115         self.assertEqual(4,d2.getNumberOfFields())
116         m = d2.getMeshes().getMeshAtPos(0)
117         self.assertEqual("maa1",m.getName())
118         self.assertEqual(6,len(m.getGroupsNames()))
119         self.assertEqual(3,m.getMeshDimension())
120         groups = m.getGroupsNames()
121         self.assertTrue( "groupe1" in groups )
122         self.assertTrue( "groupe2" in groups )
123         self.assertTrue( "groupe3" in groups )
124         self.assertTrue( "groupe4" in groups )
125         self.assertTrue( "groupe5" in groups )
126         self.assertTrue( "maa1" in groups )
127         self.assertEqual(16,m.getSizeAtLevel(0))
128         um0 = m.getMeshAtLevel(0)
129         self.assertEqual(12, um0.getNumberOfCellsWithType( NORM_TETRA4 ))
130         self.assertEqual(2, um0.getNumberOfCellsWithType( NORM_PYRA5 ))
131         self.assertEqual(2, um0.getNumberOfCellsWithType( NORM_HEXA8 ))
132         um1 = m.getMeshAtLevel(-1)
133         #self.assertEqual(2, um1.getNumberOfCellsWithType( NORM_TRI3 ))
134         self.assertEqual(1, um1.getNumberOfCellsWithType( NORM_TRI3 ))
135         pointeUM0 = pointeMedMesh.getMeshAtLevel(0)
136         self.assertTrue(m.getCoords().isEqualWithoutConsideringStr(pointeMedMesh.getCoords(),1e-12))
137         self.assertEqual( um0.getMeasureField(False).accumulate(0),
138                           pointeUM0.getMeasureField(False).accumulate(0),1e-12)
139         # check fields
140         # fieldnodedouble
141         fieldnodedoubleTS1 = pointeMed.getFields().getFieldWithName("fieldnodedouble")
142         fieldnodedoubleTS2 = d2.getFields().getFieldWithName("fieldnodedouble")
143         self.assertEqual( fieldnodedoubleTS1.getInfo(), fieldnodedoubleTS2.getInfo())
144         self.assertEqual( fieldnodedoubleTS1.getNumberOfTS(), fieldnodedoubleTS2.getNumberOfTS())
145         io1 = fieldnodedoubleTS1.getIterations()
146         io2 = fieldnodedoubleTS2.getIterations()
147         for i in range(fieldnodedoubleTS1.getNumberOfTS() ):
148             fnd1 = fieldnodedoubleTS1.getFieldOnMeshAtLevel(ON_NODES, io1[i][0],io1[i][1],pointeUM0)
149             fnd2 = fieldnodedoubleTS2.getFieldOnMeshAtLevel(ON_NODES, io2[i][0],io2[i][1],um0)
150             self.assertTrue( fnd1.getArray().isEqual( fnd2.getArray(), 1e-12 ))
151         # fieldcelldoublevector
152         fieldnodedoubleTS1 = pointeMed.getFields().getFieldWithName("fieldcelldoublevector")
153         fieldnodedoubleTS2 = d2.getFields().getFieldWithName("fieldcelldoublevector")
154         self.assertEqual( fieldnodedoubleTS1.getInfo(), fieldnodedoubleTS2.getInfo())
155         self.assertEqual( fieldnodedoubleTS1.getNumberOfTS(), fieldnodedoubleTS2.getNumberOfTS())
156         io1 = fieldnodedoubleTS1.getIterations()
157         io2 = fieldnodedoubleTS2.getIterations()
158         for i in range(fieldnodedoubleTS1.getNumberOfTS() ):
159             fnd1 = fieldnodedoubleTS1.getFieldOnMeshAtLevel(ON_CELLS, io1[i][0],io1[i][1],pointeUM0)
160             fnd2 = fieldnodedoubleTS2.getFieldOnMeshAtLevel(ON_CELLS, io2[i][0],io2[i][1],um0)
161             self.assertAlmostEqual( fnd1.accumulate(0), fnd2.accumulate(0) )
162             self.assertAlmostEqual( fnd1.accumulate(1), fnd2.accumulate(1) )
163             self.assertAlmostEqual( fnd1.accumulate(2), fnd2.accumulate(2) )
164             # Field on 2 faces
165             fieldOnFaces = d2.getFields().getFieldWithName(f1.getName())
166             io1 = fieldOnFaces.getIterations()
167             fof = fieldOnFaces.getFieldOnMeshAtLevel(f1.getTypeOfField(),io1[i][0],io1[i][1],um1)
168             self.assertTrue( d.isEqual( fof.getArray(), 1e-12 ))
169             pass
170         del sr
171         os.remove( sauvFile )
172         pass
173
174     def testSauv2MedWONodeFamilyNum(self):
175         """test for issue 0021673: [CEA 566] Bug in SauvWriter when writing meshes
176         having no family ids on nodes."""
177
178         myCoords=DataArrayDouble.New([-0.3,-0.3, 0.2,-0.3, 0.7,-0.3, -0.3,0.2, 0.2,0.2, 0.7,0.2, -0.3,0.7, 0.2,0.7, 0.7,0.7 ],9,2)
179         targetConn=[0,3,4,1, 1,4,2, 4,5,2, 6,7,4,3, 7,8,5,4];
180         targetMesh=MEDCouplingUMesh.New("BugInSauvWriter",2);
181         targetMesh.allocateCells(5);
182         targetMesh.insertNextCell(NORM_TRI3,3,targetConn[4:7]);
183         targetMesh.insertNextCell(NORM_TRI3,3,targetConn[7:10]);
184         targetMesh.insertNextCell(NORM_QUAD4,4,targetConn[0:4]);
185         targetMesh.insertNextCell(NORM_QUAD4,4,targetConn[10:14]);
186         targetMesh.insertNextCell(NORM_QUAD4,4,targetConn[14:18]);
187         targetMesh.finishInsertingCells();
188         targetMesh.setCoords(myCoords);
189         #
190         m=MEDFileUMesh.New()
191         m.setMeshAtLevel(0,targetMesh)
192         # start of bug
193         fam=DataArrayInt.New(targetMesh.getNumberOfNodes())
194         fam[:]=0
195         #m.setFamilyFieldArr(1,fam)
196         #end of bug
197
198         ms=MEDFileMeshes.New()
199         ms.setMeshAtPos(0,m)
200         meddata=MEDFileData.New()
201         meddata.setMeshes(ms)
202
203         medFile = "BugInSauvWriter.sauv"
204         sw=SauvWriter.New();
205         sw.setMEDFileDS(meddata);
206         sw.write(medFile);
207
208         os.remove( medFile )
209         pass
210
211     def testSauv2MedOnPipe1D(self):
212         """test for issue 0021745: [CEA 600] Some missing groups in mesh after reading a SAUV file with SauvReader."""
213         sauvFile="Test_sauve_1D.sauv"
214         # Make a sauve file with a qudratic 1D mesh
215         m=MEDCouplingUMesh.New("pipe1D",1)
216         m.allocateCells(2);
217         targetConn=[0,2,1, 2,4,3]
218         m.insertNextCell(NORM_SEG3,3,targetConn[0:3])
219         m.insertNextCell(NORM_SEG3,3,targetConn[3:6])
220         m.finishInsertingCells();
221         # coords
222         coords=[ 0.,1.,2.,4.,5. ];
223         c=DataArrayDouble.New()
224         c.setValues(coords,5,1)
225         m.setCoords(c)
226         # MEDFileUMesh
227         mm=MEDFileUMesh.New()
228         mm.setName(m.getName())
229         mm.setDescription("1D mesh")
230         mm.setCoords(c)
231         mm.setMeshAtLevel(0,m);
232         # MEDFileData
233         mfd1 = MEDFileData.New()
234         ms=MEDFileMeshes.New(); ms.setMeshAtPos(0,mm)
235         mfd1.setMeshes(ms)
236         # write
237         sw=SauvWriter.New()
238         sw.setMEDFileDS(mfd1)
239         sw.write(sauvFile)
240         # Check connectivity read from the sauv file
241         sr = SauvReader.New(sauvFile)
242         mfd2 = sr.loadInMEDFileDS()
243         mfMesh = mfd2.getMeshes()[0]
244         mesh = mfMesh.getMeshAtLevel(0)
245         self.assertTrue(mesh.getNodalConnectivity().isEqual(m.getNodalConnectivity()))
246         #
247         del sr
248         os.remove(sauvFile)
249         pass
250
251     @unittest.skipUnless(HasXDR(),"requires XDR")
252     def testMissingGroups(self):
253         """test for issue 0021749: [CEA 601] Some missing groups in mesh after reading a SAUV file with SauvReader."""
254         sauvFile = os.path.join(self.__getResourcesDirectory(),"BDC-714.sauv")
255         self.assertTrue( os.access( sauvFile, os.F_OK))
256         name_of_group_on_cells='Slice10:ABSORBER'
257         name_of_group_on_cells2='Slice10:00LR'
258         sr=SauvReader.New(sauvFile)
259         mfd2=sr.loadInMEDFileDS()
260         mfMesh=mfd2.getMeshes()[0]
261         #
262         self.assertTrue(name_of_group_on_cells in mfMesh.getGroupsNames())
263         self.assertTrue(name_of_group_on_cells2 in mfMesh.getGroupsNames())
264         self.assertEqual(270,len(mfMesh.getGroupsNames()))
265         #
266         ids1=mfMesh.getGroupArr(0,name_of_group_on_cells)
267         ids2=mfMesh.getGroupArr(0,name_of_group_on_cells2)
268         ids1.setName("")
269         ids2.setName("")
270         self.assertTrue(ids1.isEqual(ids2))
271         pass
272
273     def testGaussPt(self):
274         """issue 22321: [CEA 933] Bug when reading a sauve file containing field on Gauss Pt.
275         The problem was that a field ON_GAUSS_PT was created but no Gauss Localization
276         was defined"""
277
278         # create a MEDFileData with a field ON_GAUSS_PT: 9 Gauss points, on 4 QUAD8 elements
279         f=MEDCouplingFieldDouble(ON_GAUSS_PT)
280         m=MEDCouplingUMesh("mesh",2) ; m.allocateCells()
281         m.insertNextCell(NORM_QUAD8,[0,2,4,6,1,3,5,7])
282         m.insertNextCell(NORM_QUAD8,[2,9,11,4,8,10,12,3])
283         m.insertNextCell(NORM_QUAD8,[6,4,14,16,5,13,15,17])
284         m.insertNextCell(NORM_QUAD8,[4,11,19,14,12,18,20,13])
285         m.setCoords(DataArrayDouble([(0,0),(0,0.25),(0,0.5),(0.25,0.5),(0.5,0.5),(0.5,0.25),(0.5,0),(0.25,0),(0,0.75),(0,1),(0.25,1),(0.5,1),(0.5,0.75),(0.75,0.5),(1,0.5),(1,0.25),(1,0),(0.75,0),(0.75,1),(1,1),(1,0.75)],21,2))
286         f.setMesh(m)
287         arr=DataArrayDouble(4*9*2) ; arr.iota() ; arr.rearrange(2) ; arr.setInfoOnComponents(["YOUN []","NU []"])
288         f.setArray(arr)
289         refCoo=[-1,-1,1,-1,1,1,-1,1,0,-1,1,0,0,1,-1,0]
290         gpCoo=[-0.7,-0.7,0.7,-0.7,0.7,0.7,-0.7,0.7,0,-0.7,0.7,0,0,0.7,-0.7,0,0,0]
291         wgt=[0.3,0.3,0.3,0.3,0.4,0.4,0.4,0.4,0.7]
292         f.setGaussLocalizationOnType(NORM_QUAD8,refCoo,gpCoo,wgt)
293         f.setName("SIGT")
294         f.checkConsistencyLight()
295         #
296         mm=MEDFileUMesh()
297         mm.setMeshAtLevel(0,m)
298         mfm = MEDFileMeshes()
299         mfm.pushMesh( mm )
300         ff=MEDFileField1TS()
301         ff.setFieldNoProfileSBT(f)
302         mfmts = MEDFileFieldMultiTS()
303         mfmts.pushBackTimeStep(ff)
304         mff = MEDFileFields()
305         mff.pushField( mfmts )
306         mfd = MEDFileData.New()
307         mfd.setFields( mff )
308         mfd.setMeshes( mfm )
309
310         # convert the MED file to a SAUV file
311         sauvFile = "SauvLoaderTest_testGaussPt.sauv"
312         sw=SauvWriter.New();
313         sw.setMEDFileDS(mfd);
314         sw.write(sauvFile);
315
316         # convert the SAUV file back to MED
317         sr=SauvReader.New(sauvFile);
318         d2=sr.loadInMEDFileDS();
319
320         self.assertEqual( 1, d2.getNumberOfFields() )
321         self.assertEqual( 1, d2.getNumberOfMeshes() )
322         mfm2 = d2.getMeshes()[0]
323         mff2 = d2.getFields()[0]
324         m2 = mfm2.getMeshAtLevel(0)
325         f2 = mff2.getTimeStepAtPos(0).getFieldOnMeshAtLevel(f.getTypeOfField(),0,mfm2)
326         f2.setGaussLocalizationOnType(NORM_QUAD8,refCoo,gpCoo,wgt) # not stored in SAUV
327         #f2.setOrder( f.getTime()[2] ) # not stored in SAUV
328         self.assertTrue( m2.isEqual( m, 1e-12 ))
329         self.assertTrue( f2.isEqual( f, 1e-12, 1e-12 ))
330         del sr
331         os.remove( sauvFile )
332         pass
333
334     def testSauvWriterGroupWithOneFamily(self):
335         """
336         This test checks an option for sauv writing. It is requested here to copy a group from a family if a group is lying on a single family.
337         """
338         import re
339         mfd=MEDLoaderDataForTest.buildAMEDFileDataWithGroupOnOneFamilyForSauv()
340         sauvFile = "mesh.sauv"
341         sw=SauvWriter.New()
342         sw.setMEDFileDS(mfd)
343         self.assertTrue(not sw.getCpyGrpIfOnASingleFamilyStatus())
344         sw.setCpyGrpIfOnASingleFamilyStatus(True)
345         self.assertTrue(sw.getCpyGrpIfOnASingleFamilyStatus())
346         sw.write(sauvFile)
347
348         f = open(sauvFile)
349         # String pattern for the header of the sub meshes record ("PILE" number, number of named objects, number of objects)
350         pattern_pile= re.compile(r'\sPILE\sNUMERO\s+(?P<number>[0-9]+)NBRE\sOBJETS\sNOMMES\s+(?P<nbnamed>[0-9]+)NBRE\sOBJETS\s+(?P<nbobjects>[0-9]+)')
351         # String pattern for a sub mesh header (cell type, number of components and three numbers)
352         pattern_header=re.compile(r'\s+(?P<type>[0-9]+)\s+(?P<nbsubs>[0-9]+)\s+[0-9]+\s+[0-9]+\s+[0-9]+')
353
354         nbobjects=0
355         line = f.readline()
356         while(line):
357             match_pile = pattern_pile.match(line)
358             if match_pile:
359                 number=int(match_pile.group("number"))
360                 if number == 1:
361                     nbnamed=int(match_pile.group("nbnamed"))
362                     nbobjects=int(match_pile.group("nbobjects"))
363                     break
364                 pass
365             line=f.readline()
366             pass
367
368         # Skipping the objects names
369         f.readline()
370         # Skipping the objects ids
371         f.readline()
372
373         # Looking for each sub-mesh header
374         line = f.readline()
375         cur_object=0
376         while(line and cur_object < nbobjects):
377             match_header=pattern_header.match(line)
378             if match_header:
379                 cell_type=int(match_header.group("type"))
380                 nb_subs=int(match_header.group("nbsubs"))
381                 # Looking for a compound object
382                 if cell_type == 0:
383                     # Testing if there is only one component
384                     self.assertTrue(nb_subs > 1)
385                 else:
386                     f.readline()
387                     f.readline()
388                     cur_object = cur_object + 1
389                     pass
390                 pass
391             line=f.readline()
392             pass
393         f.close()
394         os.remove(sauvFile)
395         pass
396
397     pass
398
399 unittest.main()