doc/tutorial/medcoupling_3Dcube.rst

   1
   2 Meshing a 3D cube
   3 -----------------
   4
   5 Objective
   6 ~~~~~~~~~
   7
   8 The meshing class of the SALOME MED module allows user to create a mesh from scratch.
   9 In this example we propose to build a mesh on a 3D cube by three methods (classical method, method by extrusion and grid method). Each cell of the mesh must be a hexaedron.
  10 We see also how creating a group.
  11 Then we create a field on all the 3D cube.
  12 Each result will be saved in a med file, and then visualized with the SMESH module of Salome.
  13
  14 In spite of a mesh in MEDCoupling has only one dimension, it's possible to save some meshes with different dimension in one med file. We will see this method.
  15
  16 .. image:: images/Mesh_cube3D.jpg
  17
  18
  19 Beginning of implementation
  20 ~~~~~~~~~~~~~~~~~~~~~~~~~~~
  21
  22 To implement this exercice we use the python language script and import the MEDCoupling and MEDLoader parts of the MED module. We need also mathematical functions, so we import the python math module::
  23
  24         from MEDCoupling import *
  25         from MEDLoader import *
  26         from math import *
  27
  28 You must define 3 variables for space dimension, number of nodes on each dimension and total number of nodes::
  29
  30         MeshDim3D = 3
  31         MeshDim2D = 2
  32         N = ...
  33         NbCell2D = (N-1)*(N-1)
  34         NbCell3D = NbCell2D*(N-1)
  35         NbNode2D = N*N
  36         NbNode3D = NbNode2D*N
  37
  38
  39 Classical method
  40 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  41
  42 First instantiate a meshing object. Therefore, we need to define :
  43
  44  * its name
  45  * its dimension
  46  * its number of cells
  47
  48 .. note:: All this initialisation are necessary. If one lacks, you will have a segmentation fault!.
  49
  50 ::
  51
  52         mesh=MEDCouplingUMesh.New()
  53         mesh.setMeshDimension(3)
  54         mesh.allocateCells(...)
  55         mesh.setName("3Dcube")
  56
  57 Definition of nodes coordinates
  58 ```````````````````````````````
  59
  60 Define the coordinates of the nodes of the 3D cube mesh, and then use the setCoords function to set it::
  61
  62         # Initialisation of coordinates
  63         coordinates = []
  64         for k in range(N):
  65                 for j in range(N):
  66                         for i in range(N):
  67                                 coordinates.append(...)
  68
  69         myCoords = DataArrayDouble.New()
  70         myCoords.setValues(coordinates,nbOfNodes,3)
  71         mesh.setCoords(myCoords)
  72
  73
  74 Definition of hexahedrons connectivity
  75 ``````````````````````````````````````
  76 For each hexahedron of the mesh, you have to give its connectivity: the list of the nodes which belong to the hexahedron. The order of the nodes in the connectivity array must respect the MEDCoupling format (see the following figure).
  77
  78 .. image:: images/cube.jpg
  79
  80 .. warning:: Connectivity elements begin to 0 from (n-1) elements
  81
  82 ::
  83
  84         connectivity = []
  85         for k in range(N-1):
  86                 for j in range(N-1):
  87                         for i in range(N-1):
  88                                 inode = ...
  89                                 connectivity.append(inode)
  90                                 connectivity.append(...)
  91
  92         # Adding cells in meshing
  93         for i in range(nbOfCells):
  94                 mesh.insertNextCell(NORM_HEXA8,8,connectivity[8*i:8*(i+1)])
  95                 pass
  96
  97         # Check mesh consistency:
  98         mesh.checkConsistencyLight()
  99
 100 Method by extrusion
 101 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 102
 103 In order to create a extruded mesh, we need one 2D mesh and one 1D mesh which define the vector of extrusion and the number of steps.
 104
 105 Definition of 2D mesh
 106 ``````````````````````````````````````
 107 ::
 108
 109         coordinates = []
 110         for j in range(N):
 111                 for i in range(N):
 112                         coordinates.append(...)
 113                         ...
 114         Connectivities = [...]
 115         myCoords = DataArrayDouble.New()
 116         myCoords.setValues(coordinates,NbNode2D,MeshDim2D)
 117
 118         m1 = MEDCouplingUMesh.New()
 119         m1.setMeshDimension(MeshDim2D)
 120         m1.allocateCells(NbCell2D)
 121         m1.setCoords(myCoords)
 122         m1.setName("2D_Support")
 123
 124         for i in range(NbCell2D):
 125                 m1.insertNextCell(NORM_QUAD4,4,Connectivities[4*i:4*(i+1)])
 126         m1.changeSpaceDimension(3)
 127
 128 Definition of 1D mesh
 129 ``````````````````````````````````````
 130 ::
 131
 132         coords = [ ... ]
 133         conn   = [ ... ]
 134         m2 = MEDCouplingUMesh.New()
 135         m2.setMeshDimension(1)
 136         m2.allocateCells(3)
 137         m2.insertNextCell(NORM_SEG2,2,conn[0:2])
 138         m2.insertNextCell(NORM_SEG2,2,conn[2:4])
 139         m2.insertNextCell(NORM_SEG2,2,conn[4:6])
 140         myCoords1D=DataArrayDouble.New()
 141         myCoords1D.setValues(coords,4,1)
 142         m2.setCoords(myCoords1D)
 143         m2.changeSpaceDimension(3)
 144
 145 Definition of extruded mesh
 146 ``````````````````````````````````````
 147
 148 Since 1D meshing isn't well oriented (along 0x vector), we need to imply a transformation on it.
 149 Then, we can extrude 2D meshing.
 150
 151 ::
 152
 153         center = [...]
 154         vector = [...]
 155         m2.rotate(...)
 156         m3 = m1.buildExtrudedMesh(m2,0)
 157         m3.setName("Extrusion")
 158
 159 Grid method
 160 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 161
 162 it's the easiest way to create a grid since you have no connectivity to set. They will be automatically setting. Incrementation of ids will be made first along Ox axis, then along Oy axis and finally along Oz axis.
 163
 164 ::
 165
 166         mesh=MEDCouplingCMesh.New()
 167         coordsX=DataArrayDouble.New()
 168         arrX=[ ... ]
 169         coordsX.setValues(arrX,4,1)
 170         coordsY=DataArrayDouble.New()
 171         arrY=[ ... ]
 172         coordsY.setValues(arrY,4,1)
 173         coordsZ=DataArrayDouble.New()
 174         arrZ=[ ... ]
 175         coordsZ.setValues(arrZ,4,1)
 176         mesh.setCoords(coordsX,coordsY,coordsZ)
 177
 178 Really in order to save this mesh, you need to transform this structured mesh to an unstructerd mesh.
 179 ::
 180
 181         meshU=mesh.buildUnstructured()
 182         meshU.setName("Grid")
 183
 184 Creation of groups
 185 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 186
 187 A group is a set of cells defining by their id. This set must the input for creating a group.
 188 Generally ids cells using in group are known. So you just need put these ids in a DataArray.
 189 ::
 190
 191         tabIdCells = DataArrayInt.New()
 192         IdCells = [ ... ]
 193         tabIdCells.setValues(IdCells,...)
 194
 195
 196 .. note:: It's also possible to retrieve ids cells from a submesh of the principal mesh.
 197
 198 ::
 199
 200         ret,tabIdCells = mesh.areCellsIncludedIn(subMesh,0)
 201
 202
 203 Once the DataArray is created, some initializations have to be done.
 204 ::
 205
 206         # Definition of the name group
 207         tabIdCells.setName("meshGroup")
 208
 209
 210 In order to add a group on a mesh, you need to transform your unstructured mesh in a file unstructured mesh.
 211 Moreover, we need to define:
 212
 213  * its name
 214  * its description
 215  * its coordinates
 216  * its dimension
 217  * its number of cells
 218
 219
 220 ::
 221
 222         # Passing MEDCoupling to MEDFile
 223         fmeshU = MEDFileUMesh.New()
 224         fmeshU.setName("Grid")
 225         fmeshU.setDescription("IHopeToConvinceLastMEDMEMUsers")
 226         myCoords = meshU.getCoords()
 227         fmeshU.setCoords(myCoords)
 228         fmeshU.setMeshAtLevel(0,meshU)
 229
 230 Then, you can
 231 Therefore, you need to define the level (ie. the dimension) of the group.
 232 This information is given by a number : 0,-1 or -2.
 233
 234  * 0 means the same level at mesh
 235
 236 ::
 237
 238         fmeshU.setGroupsAtLevel(0,[tabIdCells],False)
 239
 240 Create field on 3D cube
 241 ~~~~~~~~~~~~~~~~~~~~~~~
 242
 243 First you need to create a CouplingField and initialize some data:
 244
 245  * its name
 246  * its support (ie mesh)
 247  * its nature
 248  * its values
 249
 250
 251 The field will be a sin function dependent of distance of the barycenter of each cell from origin. So we need to create a barycenter field on the 3D mesh::
 252
 253         # Creation of field : with following definition
 254         # => Definition of the mesh support
 255         # => Definition of field name
 256         # => Definition of field nature
 257         field = MEDCouplingFieldDouble.New(ON_CELLS)
 258         field.setMesh(mesh)
 259         field.setName("field")
 260         field.setNature(ExtensiveMaximum)
 261
 262         # Computing and setting field values
 263         myCoords=DataArrayDouble.New()
 264         sampleTab=[]
 265         bar = mesh.computeCellCenterOfMass()
 266         print bar.getNbOfElems()
 267         for i in range(nbOfCells):
 268                 x = bar.getIJ(...)
 269                 y = bar.getIJ(...)
 270                 z = bar.getIJ(...)
 271                 d = sqrt(x*x+y*y+z*z)
 272                 sinus = sin(d)
 273         .       sampleTab.append(sinus)
 274
 275         myCoords.setValues(sampleTab,nbOfCells,1)
 276         field.setArray(myCoords)
 277
 278
 279 Saving the mesh in a med file
 280 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 281 General Case
 282 ````````````
 283 ::
 284
 285         medFileName = "MEDCoupling_Extrudedcube3D.med"
 286         MEDLoader.WriteUMesh(medFileName,meshU,True)
 287
 288 .. note:: True / False in Write* functions : True for overwriting existing file and False for adding in existing file
 289
 290 Multi mesh Case
 291 ````````````````
 292
 293 In spite of a MEDCoupling mesh has only one dimension, it's possible to genrate a file with multi dimension.
 294 Therefore, you need to create as meshes as necessary dimensions.
 295
 296 You have to give the connectivity of the faces on the bottom face of the 3D cube: the list of the nodes which belong to the face.
 297 The connectivity must respect following figure:
 298
 299 .. image:: images/face.jpg
 300
 301 ::
 302
 303         # Extraction of surfacic meshing
 304         pt=[0.,0.,0.]
 305         vec=[0.,0.,1.]
 306         nodes = mesh.findNodesOnPlane(pt,vec,1e-12)
 307         mesh2D = mesh.buildFacePartOfMySelfNode(nodes,True)
 308         #print mesh2D
 309         mesh2D.setName("3Dcube")
 310         mesh2D.checkConsistencyLight()
 311
 312         medFileName = "MEDCoupling_cube3D.med"
 313         meshes=[mesh2D,mesh]
 314         MEDLoader.WriteUMeshes(medFileName,meshes,True)
 315
 316 Group Case
 317 ````````````
 318 ::
 319
 320         medFileName = "MEDCoupling_Gridcube3D.med"
 321         fmeshU.write(medFileName,2)
 322
 323 Saving the fields in the med file
 324 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 325
 326 ::
 327
 328         MEDLoader.WriteField(medFileName,field,False)
 329
 330 Visualize the mesh with the SMESH module of Salome
 331 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 332
 333 Launch Salome platform, then select SMESH module and import your MED file. First You can see the number of elements of your mesh. For that, select your mesh in the object browser, set select Mesh menu and "Advanced Mesh Info" option. Verify the number of faces and the number of hexahedrons. To visualize your mesh: click right bottom on your mesh and select "Show" option. You can also visualize your groups. Select one group, click right bottom on your group and select "Show only" option.
 334
 335 Visualize the fields with the VISU module of Salome
 336 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 337
 338 Launch Salome platform, then select VISU module and import your MED file. You can see in the object browser the 2 fields you have created. Then you have to create a scalar map on each field to visualize them.
 339
 340 .. image:: images/Field_Cube3D.jpg
 341
 342 Solution
 343 ~~~~~~~~
 344
 345 :ref:`python_testMEDCouplingcube_solution`