doc/tutorial/atestMEDCouplingUMesh1.rst

   1
   2 .. _python_testMEDCouplingumesh1_solution:
   3
   4 Playing with unstructured mesh
   5 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
   6
   7 ::
   8
   9         import MEDCoupling as mc
  10
  11         # Build a 3D mesh from scratch mixing HEXA8 and POLYHED
  12         coords=[0.,0.,0., 1.,1.,0., 1.,1.25,0., 1.,0.,0., 1.,1.5,0., 2.,0.,0., 2.,1.,0., 1.,2.,0., 0.,2.,0., 3.,1.,0.,
  13                 3.,2.,0., 0.,1.,0., 1.,3.,0., 2.,2.,0., 2.,3.,0.,
  14                 0.,0.,1., 1.,1.,1., 1.,1.25,1., 1.,0.,1., 1.,1.5,1., 2.,0.,1., 2.,1.,1., 1.,2.,1., 0.,2.,1., 3.,1.,1.,
  15                 3.,2.,1., 0.,1.,1., 1.,3.,1., 2.,2.,1., 2.,3.,1.,
  16                 0.,0.,2., 1.,1.,2., 1.,1.25,2., 1.,0.,2., 1.,1.5,2., 2.,0.,2., 2.,1.,2., 1.,2.,2., 0.,2.,2., 3.,1.,2.,
  17                 3.,2.,2., 0.,1.,2., 1.,3.,2., 2.,2.,2., 2.,3.,2.,
  18                 0.,0.,3., 1.,1.,3., 1.,1.25,3., 1.,0.,3., 1.,1.5,3., 2.,0.,3., 2.,1.,3., 1.,2.,3., 0.,2.,3., 3.,1.,3.,
  19                 3.,2.,3., 0.,1.,3., 1.,3.,3., 2.,2.,3., 2.,3.,3.]
  20         conn=[0,11,1,3,15,26,16,18,   1,2,4,7,13,6,-1,1,16,21,6,-1,6,21,28,13,-1,13,7,22,28,-1,7,4,19,22,-1,4,2,17,19,-1,2,1,16,17,-1,16,21,28,22,19,17,
  21               1,6,5,3,16,21,20,18,   13,10,9,6,28,25,24,21, 11,8,7,4,2,1,-1,11,26,16,1,-1,1,16,17,2,-1,2,17,19,4,-1,4,19,22,7,-1,7,8,23,22,-1,8,11,26,23,-1,26,16,17,19,22,23,
  22               7,12,14,13,22,27,29,28,  15,26,16,18,30,41,31,33, 16,17,19,22,28,21,-1,16,31,36,21,-1,21,36,43,28,-1,28,22,37,43,-1,22,19,34,37,-1,19,17,32,34,-1,17,16,31,32,-1,31,36,43,37,34,32,
  23               16,21,20,18,31,36,35,33,   28,25,24,21,43,40,39,36, 26,23,22,19,17,16,-1,26,41,31,16,-1,16,31,32,17,-1,17,32,34,19,-1,19,34,37,22,-1,22,23,38,37,-1,23,26,41,38,-1,41,31,32,34,37,38,
  24               22,27,29,28,37,42,44,43, 30,41,31,33,45,56,46,48,  31,32,34,37,43,36,-1,31,46,51,36,-1,36,51,58,43,-1,43,37,52,58,-1,37,34,49,52,-1,34,32,47,49,-1,32,31,46,47,-1,46,51,58,52,49,47,
  25               31,36,35,33,46,51,50,48,  43,40,39,36,58,55,54,51, 41,38,37,34,32,31,-1,41,56,46,31,-1,31,46,47,32,-1,32,47,49,34,-1,34,49,52,37,-1,37,38,53,52,-1,38,41,56,53,-1,56,46,47,49,52,53,
  26               37,42,44,43,52,57,59,58]
  27         mesh3D = mc.MEDCouplingUMesh("mesh3D",3);
  28         mesh3D.allocateCells(18);
  29         mesh3D.insertNextCell(mc.NORM_HEXA8,conn[0:8]); mesh3D.insertNextCell(mc.NORM_POLYHED,conn[8:51]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[51:59]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[59:67]); mesh3D.insertNextCell(mc.NORM_POLYHED,conn[67:110]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[110:118]);
  30         mesh3D.insertNextCell(mc.NORM_HEXA8,conn[118:126]); mesh3D.insertNextCell(mc.NORM_POLYHED,conn[126:169]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[169:177]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[177:185]); mesh3D.insertNextCell(mc.NORM_POLYHED,conn[185:228]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[228:236]);
  31         mesh3D.insertNextCell(mc.NORM_HEXA8,conn[236:244]); mesh3D.insertNextCell(mc.NORM_POLYHED,conn[244:287]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[287:295]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[295:303]); mesh3D.insertNextCell(mc.NORM_POLYHED,conn[303:346]); mesh3D.insertNextCell(mc.NORM_HEXA8,conn[346:354]);
  32         myCoords = mc.DataArrayDouble(coords,60,3);
  33         myCoords.setInfoOnComponents(["X [m]","Y [m]","Z [m]"])
  34         mesh3D.setCoords(myCoords);
  35         mesh3D.orientCorrectlyPolyhedrons()
  36         mesh3D.sortCellsInMEDFileFrmt()
  37         mesh3D.checkCoherency()
  38         renum = mc.DataArrayInt(60) ; renum[:15]=range(15,30) ; renum[15:30]=range(15) ; renum[30:45]=range(45,60) ; renum[45:]=range(30,45)
  39         mesh3D.renumberNodes(renum,60)
  40         # Scale coordinates from meters to centimeters
  41         mesh3D.getCoords()[:] *= 100.
  42         mesh3D.getCoords().setInfoOnComponents(["X [cm]","Y [cm]","Z [cm]"])
  43         # Identify unique Z values
  44         zLev = mesh3D.getCoords()[:,2]
  45         zLev = zLev.getDifferentValues(1e-12)
  46         zLev.sort()
  47         # Extract cells from a given Z level - Solution 1
  48         tmp,cellIdsSol1 = mesh3D.buildSlice3D([0.,0.,(zLev[1]+zLev[2])/2],[0.,0.,1.],1e-12)
  49         # Idem - Solution 2
  50         bary = mesh3D.getBarycenterAndOwner()
  51         baryZ = bary[:,2]
  52         cellIdsSol2 = baryZ.getIdsInRange(zLev[1],zLev[2])
  53         # Idem - Solution 3
  54         nodeIds = mesh3D.findNodesOnPlane([0.,0.,zLev[0]],[0.,0.,1.],1e-10)
  55         mesh2D = mesh3D.buildFacePartOfMySelfNode(nodeIds,True)
  56         extMesh = mc.MEDCouplingExtrudedMesh(mesh3D,mesh2D,0)
  57         n_cells = mesh2D.getNumberOfCells()
  58         cellIdsSol3 = extMesh.getMesh3DIds()[n_cells:2*n_cells]
  59         # Compare the 3 methods
  60         print cellIdsSol1.getValues()
  61         print cellIdsSol2.getValues()
  62         print cellIdsSol3.getValues()
  63         # Extract part of the mesh
  64         mesh3DPart = mesh3D[cellIdsSol2] # equivalent to mesh3DPart = mesh3D.buildPartOfMySelf(cellIdsSol2,True)
  65         mesh3DPart.zipCoords()
  66         # Check geometric type ordering
  67         #print mesh3DPart.advancedRepr()
  68         print mesh3DPart.checkConsecutiveCellTypesAndOrder([mc.NORM_HEXA8,mc.NORM_POLYHED])
  69         print mesh3DPart.checkConsecutiveCellTypes()
  70         #print mesh3DPart.advancedRepr()
  71         # Extract cells along a line - Solution 1
  72         baryXY = bary[:,[0,1]]
  73         baryXY -= [250.,150.]
  74         magn = baryXY.magnitude()
  75         cellIds2Sol1 = magn.getIdsInRange(0.,1e-12)
  76         # Extract cells along a line - Solution 2
  77         bary2 = mesh2D.getBarycenterAndOwner()[:,[0,1]]
  78         bary2 -= [250.,150.]
  79         magn = bary2.magnitude()
  80         ids = magn.getIdsInRange(0.,1e-12)
  81         idStart = int(ids) # ids is assumed to contain only one value, if not an exception is thrown
  82         ze_range = range(idStart,mesh3D.getNumberOfCells(),mesh2D.getNumberOfCells())
  83         cellIds2Sol2 = extMesh.getMesh3DIds()[ze_range]
  84         # Construct the final sub-part
  85         mesh3DSlice2 = mesh3D[cellIds2Sol1]
  86         mesh3DSlice2.zipCoords()
  87         # Aggregate two meshes, one being the translated version of the original
  88         mesh3DSlice2bis = mesh3DSlice2.deepCpy()
  89         mesh3DSlice2bis.translate([0.,1000.,0.])
  90         mesh3DSlice2All = mc.MEDCouplingUMesh.MergeUMeshes([mesh3DSlice2,mesh3DSlice2bis])
  91         mesh3DSlice2All.writeVTK("mesh3DSlice2All.vtu")
  92         # Discover descending connectivity
  93         mesh3DSurf,desc,descIndx,revDesc,revDescIndx = mesh3D.buildDescendingConnectivity()
  94         numberOf3DCellSharing = revDescIndx.deltaShiftIndex()
  95         cellIds = numberOf3DCellSharing.getIdsNotEqual(1)
  96         mesh3DSurfInside = mesh3DSurf[cellIds]
  97         mesh3DSurfInside.writeVTK("mesh3DSurfInside.vtu")
  98