import salome
import medcoupling as mc
from math import pi
+import numpy as np
#salome.salome_init()
smesh = smeshBuilder.New()
+from salome.kernel.logger import Logger
+logger = Logger("salome.smesh.smesh_tools")
+logger.setLevel("DEBUG")
+
def smesh_create_dual_mesh(mesh_ior, output_file, adapt_to_shape=True, mesh_name="MESH"):
""" Create a dual of the mesh in input_file into output_file
shape = mesh.GetShapeToMesh()
+ # Creating output file
+ logger.debug("Creating file with mesh: "+mesh_name)
+ myfile = mc.MEDFileUMesh()
+ myfile.setName(mesh_name)
+
+
# We got a meshProxy so we need to convert pointer to MEDCoupling
int_ptr = mesh.ExportMEDCoupling(True, True)
dab = mc.FromPyIntPtrToDataArrayByte(int_ptr)
- tetras = mc.MEDFileMesh.New(dab)[0]
+ mc_mesh_file = mc.MEDFileMesh.New(dab)
+ tetras = mc_mesh_file[0]
# End of SMESH -> MEDCoupling part for dualmesh
tetras = mc.MEDCoupling1SGTUMesh(tetras)
polyh = tetras.computeDualMesh()
- dual_volume_raw = polyh.getMeasureField(True).accumulate()[0]
+
+ ## Adding skin + transfering groups on faces from tetras mesh
+ mesh2d = polyh.buildUnstructured().computeSkin()
+ mesh2d.setName(mesh_name)
+ myfile.setMeshAtLevel(-1, mesh2d)
+
+
+ for grp_name in mc_mesh_file.getGroupsOnSpecifiedLev(-1):
+ # This group is created by the export
+ if grp_name == "Group_Of_All_Faces":
+ logger.debug("Skipping group: "+ grp_name)
+ continue
+ logger.debug("Transferring group: "+ grp_name)
+
+ grp_tria = mc_mesh_file.getGroup(-1, grp_name)
+ # Retrieve the nodes in group
+ grp_nodes = grp_tria.computeFetchedNodeIds()
+ # Find all the cells lying on one of the nodes
+ id_grp_poly = mesh2d.getCellIdsLyingOnNodes(grp_nodes, False)
+
+ grp_poly = mesh2d[id_grp_poly]
+
+ # We use the interpolation to remove the element that are not really in
+ # the group (the ones that are next to a nodes nut not in the group
+ # will have the sum of their column in the enterpolation matrix equal
+ # to zero)
+ rem = mc.MEDCouplingRemapper()
+
+ rem.prepare(grp_poly, grp_tria, "P0P0")
+ m = rem.getCrudeCSRMatrix()
+ _, id_to_keep = np.where(m.sum(dtype=np.int64, axis=0) >= 1e-07)
+
+ id_grp_poly = id_grp_poly[id_to_keep.tolist()]
+ id_grp_poly.setName(grp_name)
+
+ myfile.addGroup(-1, id_grp_poly)
# Getting list of new points added on the skin
skin = tetras.buildUnstructured().computeSkin()
ptsAddedMesh = mc.MEDCouplingUMesh.Build0DMeshFromCoords( skin_polyh.getCoords()[ptsAdded] )
if adapt_to_shape:
+ logger.debug("Adapting to shape")
ptsAddedCoo = ptsAddedMesh.getCoords()
ptsAddedCooModified = ptsAddedCoo[:]
id2face = {}
for face in faces:
id2face[face.GetSubShapeIndices()[0]] = face
- print(id2face)
## Projecting each points added by the dual mesh on the surface it is
# associated with
polyh.getCoords()[ptsAdded] = ptsAddedCooModified
polyh.setName(mesh_name)
- polyh.write(output_file)
+ myfile.setMeshAtLevel(0, polyh)
+
+ logger.debug("Writting dual mesh in :"+output_file)
+ myfile.write(output_file, 2)