Lecture, Êcriture d'un fichier MED grâce à l'API avancÊe de MEDLoader
---------------------------------------------------------------------

L'API avancÊe de MEDLoader est reprÊsentÊe par les classes ``MEDFile*`` de la bibliothèque MEDLoader.

* Au plus haut niveau, pour l'ensemble du fichier: ``MEDFileData``,
* Pour l'ensemble des maillages du fichier : ``MEDFileMeshes``,
* Pour chacun des maillages : ``MEDFileMeshMultiTS``, ``MEDFileMesh``, ``MEDFileUMesh``, ``MEDFileCMesh``,  
* Pour l'ensemble des champs du fichier : ``MEDFileFields``, ``MEDFileFieldGlobs``, 
* Et enfin pour chacun des champs : ``MEDFileField1TS``, ``MEDFileFieldMultiTS``


Objectif
~~~~~~~~

Ecrire un maillage et un champ Ă  partir de rien, les relire et comparer les rĂŠsultats.

Points abordĂŠs : en utilisant l'API avancĂŠe de MEDLoader,

* Ecrire un fichier 
* Lire un fichier

DĂŠbut d'implĂŠmentation
~~~~~~~~~~~~~~~~~~~~~~

Cet exercice repose comme tous les autres sur le language de script Python. On charge 
le module Python ``MEDLoader``.

Pour information, le module ``MEDCoupling`` complet est inclus dans ``MEDLoader``. Pas besoin de l'importer
si ``MEDLoader`` a ĂŠtĂŠ chargĂŠ. ::

	import MEDLoader as ml

Lecture, ĂŠcriture d'un maillage
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Nous crĂŠons tout d'abord le mĂŞme maillage ``targetMesh`` que pour l'API simple. ::

	targetCoords = [-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 ]
	targetConn = [0,3,4,1, 1,4,2, 4,5,2, 6,7,4,3, 7,8,5,4]
	targetMesh = ml.MEDCouplingUMesh("MyMesh",2)
	targetMesh.allocateCells(5)
	targetMesh.insertNextCell(ml.NORM_TRI3,3,targetConn[4:7])
	targetMesh.insertNextCell(ml.NORM_TRI3,3,targetConn[7:10])
	targetMesh.insertNextCell(ml.NORM_QUAD4,4,targetConn[0:4])
	targetMesh.insertNextCell(ml.NORM_QUAD4,4,targetConn[10:14])
	targetMesh.insertNextCell(ml.NORM_QUAD4,4,targetConn[14:18])
	myCoords = ml.DataArrayDouble(targetCoords,9,2)
	myCoords.setInfoOnComponents(["X [km]","YY [mm]"])
	targetMesh.setCoords(myCoords)        

.. note:: Le maillage ``targetMesh`` est ordonnĂŠ par type gĂŠomĂŠtrique.

Nous construisons ensuite ``targetMesh1`` reprĂŠsentant les sous-constituants (*faces*) du maillage
``targetMesh``, et nous en extrayons seulement les cellules (donc ici des surfaces) [3,4,7,8]. 
Pour plus de dĂŠtails sur la connectivitĂŠ descendante, 
consulter la section :ref:`exo-umesh-desc-connec` du deuxième exercise.
Cet ensemble peut par exemple reprĂŠsenter un ensemble d'intĂŠrĂŞt pour un calcul : ::

	targetMeshConsti, _, _, _, _ = targetMesh.buildDescendingConnectivity()
	targetMesh1 = targetMeshConsti[[3,4,7,8]]
	targetMesh1.setName(targetMesh.getName())

.. note:: En Python, le underscore ``_`` signifie que l'on attend une valeur de retour, mais qu'on n'en aura pas l'usage 
	(on ne la *bind* pas).
.. note:: ``targetMesh1`` sera sauvĂŠ comme ĂŠtant une partie du mĂŞme maillage global dans le fichier MED. 
	Il doit donc avoir le mĂŞme nom. C'est lĂ  qu'on voit qu'un maillage au sens MED fichier peut mĂŠlanger les dimensions. 

On peut alors ĂŠcrire les deux maillages dans le fichier "TargetMesh2.med". ::

	meshMEDFile = ml.MEDFileUMesh()
	meshMEDFile.setMeshAtLevel(0,targetMesh)
	meshMEDFile.setMeshAtLevel(-1,targetMesh1)
	meshMEDFile.write("TargetMesh2.med",2)         # 2 stands for 'write from scratch'

Lecture, ĂŠcriture de groupes de mailles
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

CrĂŠons deux groupes de cellules sur le maillage 2D, c'est Ă  dire au niveau relatif 0 (ici, le niveau relatif 0 correspond
Ă  la 2D, le niveau -1 
correspond Ă  la 1D,  etc ...). Le premier groupe ``grp0_Lev0`` contient les cellules [0,1,3] 
le second ``grp1_Lev0`` les cellules [1,2,3,4] : ::

	grp0_0 = ml.DataArrayInt([0,1,3]) 
	grp0_0.setName("grp0_Lev0")
	grp1_0 = ml.DataArrayInt([1,2,3,4])
	grp1_0.setName("grp1_Lev0")
	meshMEDFile.setGroupsAtLevel(0, [grp0_0,grp1_0])

.. note:: On voit ĂŠvidemment ici l'importance de nommer les tableaux : c'est le nom qui sera utilisĂŠ pour le groupe. 

CrĂŠons trois groupes de niveau -1, c'est Ă  dire des groupes de faces. Le premier appelĂŠ 
``grp0_LevM1`` aux cellules [0,1], le second appelÊ ``grp1_LevM1`` aux cellules [0,1,2], et le 3ème ``grp2_LevM1``
aux cellules [1,2,3] : ::

	grp0_M1 = ml.DataArrayInt([0,1])
	grp0_M1.setName("grp0_LevM1")
	grp1_M1 = ml.DataArrayInt([0,1,2])
	grp1_M1.setName("grp1_LevM1")
	grp2_M1 = ml.DataArrayInt([1,2,3])
	grp2_M1.setName("grp2_LevM1")
	meshMEDFile.setGroupsAtLevel(-1,[grp0_M1,grp1_M1,grp2_M1])
	
Ecrivons le tout : ::
	
	meshMEDFile.write("TargetMesh2.med",2)         # 2 stands for 'write from scratch'
	
Nous pouvons ensuite re-lire le fichier MED : ::

	meshMEDFileRead = ml.MEDFileMesh.New("TargetMesh2.med") # a new is needed because it returns a MEDFileUMesh (MEDFileMesh is abstract)
	meshRead0 = meshMEDFileRead.getMeshAtLevel(0)
	meshRead1 = meshMEDFileRead.getMeshAtLevel(-1)
	print "Is level 0 in the file equal to 'targetMesh'?", meshRead0.isEqual(targetMesh,1e-12)
	print "Is level 0 in the file equal to 'targetMesh1'?", meshRead1.isEqual(targetMesh1,1e-12)

Affichons les niveaux disponibles pour le groupe ``grp0_Lev0`` : ::

	print meshMEDFileRead.getGrpNonEmptyLevels("grp0_Lev0")

Et rĂŠcupĂŠrons enfin les identifiants de cellules contenus dans le groupe ``grp0_Lev0`` : ::

	grp0_0_read = meshMEDFileRead.getGroupArr(0,"grp0_Lev0")
	print "Is group 'grp0_Lev0' equal to what is read in the file?" , grp0_0_read.isEqual(grp0_0)

Lire/ĂŠcrire des champs avec l'API avancĂŠe
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

CrĂŠons un champ de vecteurs simple, aux cellules (P0), avec un seul pas de temps, appelĂŠ ``f``. ::

	f = ml.MEDCouplingFieldDouble(ml.ON_CELLS, ml.ONE_TIME)
	f.setTime(5.6,7,8)
	f.setArray(targetMesh.computeCellCenterOfMass())
	f.setMesh(targetMesh)
	f.setName("AFieldName")

Stocker ``f`` dans un object ``MEDFileField1TS`` (un champ avec un seul pas de temps -- *one time-step, 1TS*) 
pour prĂŠparer l'ĂŠcriture MED ::

	fMEDFile = ml.MEDFileField1TS()
	fMEDFile.setFieldNoProfileSBT(f)     # No profile desired on the field, Sort By Type

Ajouter le champ au fichier "TargetMesh2.med" ::

	fMEDFile.write("TargetMesh2.med",0) # 0 is paramount to indicate that we *append* (and no overwrite) to the MED file

.. note:: Noter l'utilisation du 0 pour indiquer que nous dĂŠsirons ajouter au fichier existant.

Lire le champ : ::

	fMEDFileRead = ml.MEDFileField1TS("TargetMesh2.med",f.getName(),7,8)
	fRead1 = fMEDFileRead.getFieldOnMeshAtLevel(ml.ON_CELLS,0,meshMEDFileRead) # Quickest way, not re-reading mesh in the file.
	fRead2 = fMEDFileRead.getFieldAtLevel(ml.ON_CELLS,0)                       # Like above, but this time the mesh is read!
	print "Does the field remain OK with the quick method?", fRead1.isEqual(f,1e-12,1e-12)
	print "Does the field remain OK with the slow method?", fRead2.isEqual(f,1e-12,1e-12)
	
Lire/ĂŠcrire un champ sur un "profil"
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Nous allons maintenant voir un concept avancĂŠ des fichiers MED, Ă  savoir la possibilitĂŠ d'ĂŠcrire un champ sur seulement
une *partie* du maillage. La technique habituellement utilisÊe est plutôt de mettre des valeurs particulières (e.g. +infinity
soit 1e+300) sur les zones oĂš le champ n'a pas de sens, permettant ainsi de repĂŠrer en plus des bugs ĂŠventuels lors du calcul.

Le mode de fonctionnement avec les profils reste donc peu courant.

Construisons une rĂŠduction aux cellules [1,2,3] de ``f`` et appelons la ``fPart`` : ::

	pfl = ml.DataArrayInt([1,2,3]) 
	pfl.setName("My1stPfl")
	fPart = f.buildSubPart(pfl)
	fPart.setName("fPart")

La stocker dans la structure ``MEDFileField1TS`` et invoquer ``setFieldProfile()``. ::

	fMEDFile2 = ml.MEDFileField1TS()
	fMEDFile2.setFieldProfile(fPart,meshMEDFileRead,0,pfl) # 0 is the relative level (here 0 means 2D)
	fMEDFile2.write("TargetMesh2.med",0) # 0 is paramount to indicate that we *append* (and no overwrite) to the MED file

Lire le champ ``fPart`` du fichier "TargetMesh2.med" et les identifiants de cellules correspondant. ::

	fMEDFileRead2 = ml.MEDFileField1TS("TargetMesh2.med",fPart.getName(),7,8)
	fPartRead, pflRead = fMEDFileRead2.getFieldWithProfile(ml.ON_CELLS,0,meshMEDFileRead)
	print "Is the partial field correclty read?", fPartRead.isEqualWithoutConsideringStr(fPart.getArray(),1e-12)
	print "Is the list of cell identifiers matching?", pflRead.isEqualWithoutConsideringStr(pfl)

Solution
~~~~~~~~

:ref:`python_testMEDLoaderAdvancedAPI1_solution`