src/SMESH_SWIG/smesh_tools.py

   1 #!/usr/bin/env python3
   2
   3 import sys
   4 import salome
   5 import medcoupling as mc
   6 from math import pi
   7
   8 #salome.salome_init()
   9
  10 import GEOM
  11 from salome.geom import geomBuilder
  12
  13 geompy = geomBuilder.New()
  14
  15 import  SMESH, SALOMEDS
  16 from salome.smesh import smeshBuilder
  17
  18 smesh = smeshBuilder.New()
  19
  20 def smesh_create_dual_mesh(mesh_ior, output_file, adapt_to_shape=True, mesh_name="MESH"):
  21     """ Create a dual of the mesh in input_file into output_file
  22
  23     Args:
  24         mesh_ior (string): corba Id of the Tetrahedron mesh
  25         output_file (string): dual mesh file
  26     """
  27     # Import mesh from file
  28     mesh = salome.orb.string_to_object(mesh_ior)
  29     if not mesh:
  30         raise Exception("Could not find mesh using id: ", mesh_ior)
  31
  32     shape = mesh.GetShapeToMesh()
  33
  34     # We got a meshProxy so we need to convert pointer to MEDCoupling
  35     int_ptr = mesh.ExportMEDCoupling(True, True)
  36     dab = mc.FromPyIntPtrToDataArrayByte(int_ptr)
  37     tetras = mc.MEDFileMesh.New(dab)[0]
  38     # End of SMESH -> MEDCoupling part for dualmesh
  39
  40     tetras = mc.MEDCoupling1SGTUMesh(tetras)
  41     polyh = tetras.computeDualMesh()
  42     dual_volume_raw = polyh.getMeasureField(True).accumulate()[0]
  43
  44     # Getting list of new points added on the skin
  45     skin = tetras.buildUnstructured().computeSkin()
  46     skin_polyh = polyh.buildUnstructured().computeSkin()
  47     allNodesOnSkinPolyh = skin_polyh.computeFetchedNodeIds()
  48     allNodesOnSkin = skin.computeFetchedNodeIds()
  49     ptsAdded = allNodesOnSkinPolyh.buildSubstraction(allNodesOnSkin)
  50     ptsAddedMesh = mc.MEDCouplingUMesh.Build0DMeshFromCoords( skin_polyh.getCoords()[ptsAdded] )
  51
  52     if adapt_to_shape:
  53         ptsAddedCoo = ptsAddedMesh.getCoords()
  54         ptsAddedCooModified = ptsAddedCoo[:]
  55
  56         # Matching faces with their ids
  57         faces = geompy.ExtractShapes(shape, geompy.ShapeType["FACE"], True)
  58         id2face = {}
  59         for face in faces:
  60             id2face[face.GetSubShapeIndices()[0]] = face
  61         print(id2face)
  62
  63         ## Projecting each points added by the dual mesh on the surface it is
  64         # associated with
  65         for i, tup in enumerate(ptsAddedCooModified):
  66             vertex = geompy.MakeVertex(*tuple(tup))
  67             shapes = geompy.GetShapesNearPoint(shape, vertex,
  68                                                geompy.ShapeType["FACE"])
  69             prj = geompy.MakeProjection(vertex,
  70                                         id2face[shapes.GetSubShapeIndices()[0]])
  71             new_coor = geompy.PointCoordinates(prj)
  72             ptsAddedCooModified[i] = new_coor
  73
  74         polyh.getCoords()[ptsAdded] = ptsAddedCooModified
  75
  76     polyh.setName(mesh_name)
  77     polyh.write(output_file)