]> SALOME platform Git repositories - tools/medcoupling.git/blob - doc/tutorial/medloader_SplitAndMerge1_fr.rst
Salome HOME
Merge branch 'master' of https://codev-tuleap.cea.fr/plugins/git/salome/medcoupling
[tools/medcoupling.git] / doc / tutorial / medloader_SplitAndMerge1_fr.rst
1
2 Spliter et fusionner un fichier MED grâce à l'API avancée de MEDLoader
3 ----------------------------------------------------------------------
4
5 Objectif
6 ~~~~~~~~
7
8 Cet exercise présente un cas complet et avancé d'utilisation de l'API avancée de MEDLoader.
9 Le but est de créer un maillage multi-type à partir de rien avec 2 champs :
10
11 * un champ aux cellules "CellField"
12 * un champ aux noeuds "NodeField"
13  
14 Nous allons ensuite couper ces champs en deux parties (dans le but d'un traitement en parallèle par un code par exemple)
15 et aussi montrer comment re-fusionner deux champs à partir de morceaux disjoints. 
16
17 Début de l'implémentation
18 ~~~~~~~~~~~~~~~~~~~~~~~~~
19
20 Créer un unstructured mesh ``m0`` issu d'un maillage structuré (meshDim=2, spaceDim=2) de 30*30.
21 Chacune des cellules paires du maillage sera *simplexisée* (i.e. coupée en triangle - méthode ``MEDCouplingUMesh.simplexize(0)``) ::
22
23         import MEDLoader as ml
24         
25         m0 = ml.MEDCouplingCMesh()
26         arr = ml.DataArrayDouble(31,1) ; arr.iota(0.)
27         m0.setCoords(arr,arr)
28         m0 = m0.buildUnstructured()
29         m00 = m0[::2]                # Extract even cells
30         m00.simplexize(0) 
31         m01 = m0[1::2]
32         m0 = ml.MEDCouplingUMesh.MergeUMeshes([m00,m01])
33         m0.getCoords()[:] *= 1/15.   # Illustrate how to quickly rescale a mesh
34         m0.setName("mesh")
35
36 .. note:: Le ``setName()`` sur "m0" est obligatoire. Ne pas oublier que dans le contexte MED fichier 
37         le nommage correct des maillages est fondamental.
38
39 Créer les champs ``cellField`` et ``nodeField`` au pas de temps identifié à (5,6) et au pas de temps 5.6 s. ::
40
41         # Cell field
42         cellField = ml.MEDCouplingFieldDouble(ml.ON_CELLS, ml.ONE_TIME) 
43         cellField.setTime(5.6,5,6)
44         cellField.setMesh(m0)
45         cellField.setName("CellField")
46         cellField.fillFromAnalytic(1,"exp(-((x-1)*(x-1)+(y-1)*(y-1)))")
47         cellField.getArray().setInfoOnComponent(0,"powercell [W]")
48         # Node field
49         nodeField = ml.MEDCouplingFieldDouble(ml.ON_NODES,ml.ONE_TIME) 
50         nodeField.setTime(5.6,5,6)
51         nodeField.setMesh(m0)
52         nodeField.setName("NodeField")
53         nodeField.fillFromAnalytic(1,"exp(-((x-1)*(x-1)+(y-1)*(y-1)))")
54         nodeField.getArray().setInfoOnComponent(0,"powernode [W]")
55
56 On obtient par exemple pour "CellField" ceci :
57
58 .. image:: images/SplitAndMergeCell1.jpg        
59
60
61 Partitionnement de maillage
62 ~~~~~~~~~~~~~~~~~~~~~~~~~~~
63
64 Couper ``m0`` en deux parties distinctes. Les deux parties seront nommées ``proc0`` et ``proc1``. 
65 ``proc0`` sera la partie dans la boîte englobante (``MEDCouplingUMesh.getCellsInBoundingBox()``) ``[(0.,0.4),(0.,0.4)]``
66 à 1e-10 près. ``proc1`` sera le complémentaire (``DataArrayInt.buildComplement()``). ::
67
68         proc0 = m0.getCellsInBoundingBox([(0.,0.4),(0.,0.4)],1e-10)
69         proc1 = proc0.buildComplement(m0.getNumberOfCells())
70
71 .. image:: images/SplitAndMerge2.jpg
72
73 Ecriture dans 2 fichiers MED séparés
74 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
75
76 En partant du partitionnement ``proc0`` et ``proc1`` créer 2 fichiers MED appelés "proc0.med" et "proc1.med" : ::
77
78         nodeField0 = nodeField[proc0] ; cellField0 = cellField[proc0] ; cellField0.setMesh(nodeField0.getMesh())
79         nodeField1 = nodeField[proc1] ; cellField1 = cellField[proc1] ; cellField1.setMesh(nodeField1.getMesh())
80         
81         proc0_fname = "proc0.med"
82         ml.WriteField(proc0_fname, nodeField0, True)
83         ml.WriteFieldUsingAlreadyWrittenMesh(proc0_fname, cellField0)
84         
85         proc1_fname = "proc1.med"
86         ml.WriteField(proc1_fname,nodeField1,True)
87         ml.WriteFieldUsingAlreadyWrittenMesh(proc1_fname,cellField1)
88
89 Lecture et fusion des 2 fichiers MED séparés (non optimal)
90 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
91
92 Partant de "proc0.med" et de "proc1.med" lire leur "CellField" respectif avec l'API basique, 
93 agréger les deux et mettre le résultat dans ``cellField_read`` : ::
94
95         cellField0_read = ml.ReadFieldCell("proc0.med","mesh",0,"CellField",5,6)
96         cellField1_read = ml.ReadFieldCell("proc1.med","mesh",0,"CellField",5,6)
97         cellField_read = ml.MEDCouplingFieldDouble.MergeFields([cellField0_read,cellField1_read])
98
99 .. note:: On peut avoir l'impression que l'information Cell (méthode ``ReadFieldCell``) est répétée de manière abusive
100         (effectivement le champ "CellField" a été créé aux cellules), 
101         mais ne pas oublier que dans la norme MED fichier rien n'interdit qu'un champ repose sur des cellules mais 
102         aussi simultanément sur des noeuds, ou des points de Gauss ...
103
104 Comparer ``cellField_read`` et ``cellField0``. Problème, à cause de la contrainte sur la numérotation MED fichier, 
105 on a perdu la numérotation initiale. Ou plus exactement il n'y a pas
106 de moyen standard de retrouver la numérotation originale. Donc un ``MEDCouplingFieldDouble.isEqual()`` 
107 n'est pas suffisant. Utilisons un ``MEDCouplingFieldDouble.substractInPlaceDM()``
108 qui opère pour nous une renumérotation suivant une politique particulière (*policy*, voir doc html). 
109 Pour ce faire une copie profonde (*deep copy*) de ``cellField`` vers ``cellFieldCpy`` et opérer sur cette copie
110 un ``substractInPlaceDM`` (DM pour "Different Meshes", contrairement à ``substract`` qui ne marche que 
111 s'ils partagent le même maillage): ::
112
113         cellFieldCpy = cellField.deepCopy()
114         cellFieldCpy.substractInPlaceDM(cellField_read,10,1e-12)
115         cellFieldCpy.getArray().abs()
116         print cellFieldCpy.getArray().isUniform(0.,1e-12)
117
118 Opérons le même travail sur "NodeField" que celui réalisé plus haut sur "CellField".
119 La différence ici c'est qu'il va y avoir duplication de l'information à la frontière, car les noeuds limites sont partagés
120 des deux côtés : ::
121
122         nodeField0_read = ml.ReadFieldNode("proc0.med","mesh",0,"NodeField",5,6)
123         nodeField1_read = ml.ReadFieldNode("proc1.med","mesh",0,"NodeField",5,6)
124         nodeField_read = ml.MEDCouplingFieldDouble.MergeFields([nodeField0_read, nodeField1_read])
125
126 .. note:: Dans cette partie, on a donc relu le maillage une deuxième fois ce qui peut être pénalisant ...
127
128 Invoquer ``MEDCouplingUMesh.mergeNodes()`` sur ``nodeField_read`` pour lui retirer les noeuds dupliqués. 
129 Faire une deep copy appelée ``nodeFieldCpy`` de ``nodeField``
130 et supprimer encore les doublons : ::
131
132         nodeField_read.mergeNodes(1e-10)
133         nodeFieldCpy = nodeField.deepCopy()
134         nodeFieldCpy.mergeNodes(1e-10)
135
136 .. note:: A noter que ``mergeNodes()`` possède deux paramètres de précisions (*epsilons*), le premier, 
137         classique, sur la distance absolue entre les noeuds, et l'autre sur la tolérance acceptée sur les valeurs du champ. 
138         Si la valeur du champ de deux noeuds à fusionner dépasse ce deuxième epsilon, une exception est levée.
139
140 Comparer ``nodeFieldCpy`` et ``nodeField_read`` toujours en utilisant ``MEDCouplingFieldDouble.substractInPlaceDM()`` : ::
141
142         nodeFieldCpy.substractInPlaceDM(nodeField_read,10,1e-12)
143         print nodeFieldCpy.getArray().isUniform(0.,1e-12)
144
145
146 Lecture et merge des 2 fichiers MED séparés (moins facile, mais plus optimal)
147 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
148
149 Il s'agit ici de faire une méthode plus systématique et potentiellement plus générale de fusion de fichiers.
150 Pour de gros fichiers cette approche est à préférer.
151 Outre la performance, cette approche a l'avantage de pouvoir rajouter des infos.
152
153 Avec l'API avancée lire les maillages des deux fichiers "proc0.med" et "proc1.med" et agréger le résultat 
154 dans une instance ``mergeMLMesh`` de ``MEDFileUMesh``.
155 Traiter tous les niveaux de dimension (même si ici il n'y en a qu'un seul) en utilisant la méthode ``MEDFileUMesh.getNonEmptyLevels()`` 
156 sur l'instance venant de "proc0.med".
157
158 La solution donnée ci-dessous est la plus générique possible, car elle traite aussi les différents pas de temps et les
159 différents types géométriques : ::
160
161         fileNames = ["proc0.med","proc1.med"]
162         msML = [ml.MEDFileMesh.New(fname) for fname in fileNames]
163         fsML = [ml.MEDFileFields.New(fname) for fname in fileNames]
164         mergeMLMesh = ml.MEDFileUMesh()
165         mergeMLFields = ml.MEDFileFields()
166         for lev in msML[0].getNonEmptyLevels():
167                 o2nML = len(msML[0].getNonEmptyLevels())*[None]
168                 cs = [mML.getCoords() for mML in msML]
169                 mergeMLMesh.setCoords(ml.DataArrayDouble.Aggregate(cs))
170                 ms = [mML.getMeshAtLevel(lev) for mML in msML]
171                 m = ml.MEDCouplingUMesh.MergeUMeshes(ms) ; m.setCoords(mergeMLMesh.getCoords())
172                 o2nML[lev] = m.sortCellsInMEDFileFrmt()
173                 mergeMLMesh.setMeshAtLevel(lev,m)
174                 pass
175         
176         for fieldName in fsML[0].getFieldsNames():
177                 fmts = [fML[fieldName] for fML in fsML]
178                 mergeField = ml.MEDFileFieldMultiTS()
179                 for dt,it,tim in fmts[0].getTimeSteps():
180                         fts = [fmt[dt,it] for fmt in fmts]
181                         arrs = len(fts)*[None]
182                         for typp in fts[0].getTypesOfFieldAvailable():
183                                 arr1s = []
184                                 if typp == ml.ON_CELLS:
185                                         for ft in fts:
186                                                 for geoTyp,smth in ft.getFieldSplitedByType():
187                                                         if geoTyp != ml.NORM_ERROR:
188                                                                 smth1 = filter(lambda x:x[0] == ml.ON_CELLS,smth)
189                                                                 arr2s = [ft.getUndergroundDataArray()[elt[1][0]:elt[1][1]] for elt in smth1]
190                                                                 arr1s.append(ml.DataArrayDouble.Aggregate(arr2s))
191                                                                 pass
192                                                         pass
193                                                 pass
194                                         pass
195                                 else:
196                                         for ft in fts:
197                                                 smth = filter(lambda x:x[0] == ml.NORM_ERROR,ft.getFieldSplitedByType())
198                                                 arr2 = ml.DataArrayDouble.Aggregate([ft.getUndergroundDataArray()[elt[1][0][1][0]:elt[1][0][1][1]] for elt in smth])
199                                                 arr1s.append(arr2)
200                                                 pass
201                                         pass
202                                 arr = ml.DataArrayDouble.Aggregate(arr1s)
203                                 if typp == ml.ON_CELLS:
204                                      arr.renumberInPlace(o2nML[lev])
205                                 mcf = ml.MEDCouplingFieldDouble(typp,ml.ONE_TIME) ; mcf.setName(fieldName) ; mcf.setTime(tim,dt,it) ; mcf.setArray(arr)
206                                 mcf.setMesh(mergeMLMesh.getMeshAtLevel(lev)) ; mcf.checkConsistencyLight()
207                                 mergeField.appendFieldNoProfileSBT(mcf)
208                                 pass
209                         pass
210                 mergeMLFields.pushField(mergeField)
211                 pass
212         mergeMLMesh.write("merge.med",2)
213         mergeMLFields.write("merge.med",0)
214
215
216 Solution
217 ~~~~~~~~
218
219 :ref:`python_testMEDLoaderSplitAndMerge1_solution`