Salome HOME
Add examples to salome test. And add some checks to test the generated mesh
authorChristophe Bourcier <christophe.bourcier@cea.fr>
Fri, 21 Apr 2023 12:27:50 +0000 (14:27 +0200)
committerChristophe Bourcier <christophe.bourcier@cea.fr>
Fri, 21 Apr 2023 12:27:50 +0000 (14:27 +0200)
doc/salome/examples/examples.set
doc/salome/examples/ghs3dSetParametersDemo.py
doc/salome/examples/ghs3d_enfmesh.py
doc/salome/examples/ghs3d_enfvert.py
doc/salome/examples/ghs3d_optimization.py
doc/salome/examples/ghs3ddemo.py

index 11f137dab25db2428808b2683942fed6ef952fd0..ce8575174a645b0459cca6f81cc2c488cedf1b4c 100644 (file)
@@ -19,4 +19,8 @@
 
 SET(EXAMPLE_NAMES
   ghs3dSetParametersDemo
-)
\ No newline at end of file
+  ghs3ddemo
+  ghs3d_enfmesh
+  ghs3d_enfvert
+  ghs3d_optimization
+)
index b2c6f4cfa6dcd01cd6eaa49723eb559e0fc71c2f..82afb521be9faf750027853b914c8aceda96662c 100644 (file)
@@ -13,6 +13,8 @@ smesh =  smeshBuilder.New()
 box = geompy.MakeBoxDXDYDZ(200., 200., 200.)
 geompy.addToStudy(box, "box")
 
+expected_volume = 200**3
+
 # create a mesh on the box
 mgtetraMesh = smesh.Mesh(box,"box: MG-Tetra and NETGEN_1D_2D mesh")
 
@@ -30,6 +32,9 @@ MG_Tetra_Parameters_1.SetPthreadMode( 1 )       # 0 - none, 1 - aggressive, 2 -
 status = mgtetraMesh.Compute()
 assert( status )
 
+volume = smesh.GetVolume(mgtetraMesh)
+assert (volume-expected_volume)/expected_volume < 1e-12
+
 mgtetraHPCMesh       = smesh.Mesh(box,"box: MG-Tetra HPC and NETGEN_1D_2D mesh")
 status               = mgtetraHPCMesh.AddHypothesis(NETGEN_2D_Parameters_1)
 NETGEN_1D_2D_1       = mgtetraHPCMesh.Triangle(algo=smeshBuilder.NETGEN_1D2D)
@@ -44,6 +49,9 @@ MG_Tetra_Parameters_2.SetParallelMode( 1 )      # 0 - none, 1 - reproducible_giv
 status = mgtetraHPCMesh.Compute()
 assert( status )
 
+volume_HPC = smesh.GetVolume(mgtetraHPCMesh)
+assert (volume_HPC-expected_volume)/expected_volume < 1e-12
+
 
 # End of script
 
index 9752289431b863cc2ebb3500837a5d2d704e89a6..d05d718c46255f508a1c2f807c970004e6c8f90a 100644 (file)
@@ -62,13 +62,18 @@ Mesh_cylindre.AddHypothesis( MG_CADSurf_Parameters2 )
 face_cyl_faces = Mesh_cylindre.GroupOnGeom(face_cyl,'group_face_cyl', SMESH.FACE)
 face_cyl_edges = Mesh_cylindre.GroupOnGeom(face_cyl,'group_edge_cyl', SMESH.EDGE)
 face_cyl_nodes = Mesh_cylindre.GroupOnGeom(face_cyl,'group_node_cyl', SMESH.NODE)
-Mesh_cylindre.Compute()
+ok = Mesh_cylindre.Compute()
+if not ok:
+  raise Exception("Error when computing Mesh_cylindre with MG_CADSurf")
 
 # Create the mesh on the cylinder
 Mesh_box_tri = smesh.Mesh(box,"Mesh_box_tri")
 Mesh_box_tri.AddHypothesis( MG_CADSurf )
 Mesh_box_tri.AddHypothesis( MG_CADSurf_Parameters )
-Mesh_box_tri.Compute()
+ok = Mesh_box_tri.Compute()
+if not ok:
+  raise Exception("Error when computing Mesh_box_tri with MG_CADSurf")
+
 
 # Create 4 copies of the 2D mesh to test the 3 types of contraints (NODE, EDGE, FACE)
 # from the whole mesh and from groups of elements.
@@ -97,9 +102,49 @@ MG_Tetra_Parameters_edge.SetEnforcedMeshWithGroup(face_cyl_edges,SMESH.EDGE,"edg
 MG_Tetra_Parameters_face.SetEnforcedMeshWithGroup(face_cyl_faces,SMESH.FACE,"faces from face_cyl_faces")
 
 #Compute the meshes
-mesh_node.Compute()
-mesh_edge.Compute()
-mesh_face.Compute()
-mesh_mesh.Compute()
+ok = mesh_node.Compute()
+if not ok:
+  raise Exception("Error when computing mesh_node")
+ok = mesh_edge.Compute()
+if not ok:
+  raise Exception("Error when computing mesh_edge")
+ok = mesh_face.Compute()
+if not ok:
+  raise Exception("Error when computing mesh_face")
+ok = mesh_mesh.Compute()
+if not ok:
+  raise Exception("Error when computing mesh_mesh")
+
+## Get the number of nodes of a mesh or a group of mesh
+def getNbNodes(mesh):
+  try:
+    #mesh
+    nb_nodes = mesh.NbNodes()
+  except:
+    # group
+    nb_nodes = mesh.GetNumberOfNodes()
+  return nb_nodes
+
+## Check enforced mesh is correct
+def checkEnforcedMesh(mesh, enforced_group, result_group, tol=1e-7):
+  name = "Check " + mesh.GetName()
+  # test nb nodes are the same between enforced group and the result group
+  if getNbNodes(enforced_group) != getNbNodes(enforced_group):
+    raise Exception("Enforce mesh failed for %s, wrong nummber of nodes"%name)
+  # test that all nodes are imposed at the same place by finding coincident nodes in a compound mesh
+  mesh_check = smesh.Concatenate([enforced_group, result_group], 0, name=name)
+  ll_coincident_nodes = mesh_check.FindCoincidentNodes(tol)
+  coincident_nodes = [item for sublist in ll_coincident_nodes for item in sublist]
+  mesh_check.MakeGroupByIds("coincident_nodes", SMESH.NODE, coincident_nodes)
+  mesh_nodes = mesh_check.GetNodesId()
+  if len(ll_coincident_nodes) != getNbNodes(enforced_group):
+    non_coincident_nodes = list(set(mesh_nodes) - set(coincident_nodes))
+    mesh_check.MakeGroupByIds("non_coincident_nodes", SMESH.NODE, non_coincident_nodes)
+    raise Exception("Enforce mesh failed for %s"%name)
+
+checkEnforcedMesh(mesh_node, face_cyl_nodes, mesh_node.GetGroupByName("nodes from face_cyl_nodes")[0])
+checkEnforcedMesh(mesh_edge, face_cyl_edges, mesh_edge.GetGroupByName("edges from face_cyl_edges")[0])
+checkEnforcedMesh(mesh_face, face_cyl_faces, mesh_face.GetGroupByName("faces from face_cyl_faces")[0])
+checkEnforcedMesh(mesh_mesh, Mesh_cylindre.GetMesh(), mesh_mesh.GetGroupByName("faces from cylinder")[0])
 
 # End of script
index 14c25cd5074bc93bdc99fbd32c70610693ff5fb8..c4cf16588eb7f345de0865ed3dd020ac286c7355 100644 (file)
@@ -10,6 +10,8 @@
 # Ex1: Add one enforced vertex with coordinates (50,50,100) 
 #      and physical size 2.
 
+import math
+
 import salome
 salome.salome_init()
 
@@ -20,6 +22,13 @@ import SMESH
 from salome.smesh import smeshBuilder
 smesh =  smeshBuilder.New()
 
+## check if a node is at the coords given the tolerance
+def isNodeAtCoord(node_id, x, y, z, tol=1e-12):
+  xn, yn, zn = mgtetraMesh.GetNodeXYZ(node_id)
+  dist = math.sqrt((x-xn)**2+(y-yn)**2+(z-zn)**2)
+  print("dist: ", dist)
+  return dist < tol
+
 # create a box
 box = geompy.MakeBoxDXDYDZ(200., 200., 200.)
 geompy.addToStudy(box, "box")
@@ -28,7 +37,9 @@ mgtetraMesh = smesh.Mesh(box,"box: MG-Tetra and MG-CADSurf mesh")
 # create a MG-CADSurf algorithm for faces
 mgtetraMesh.Triangle(algo=smeshBuilder.MG_CADSurf)
 # compute the mesh
-mgtetraMesh.Compute()
+ok  = mgtetraMesh.Compute()
+if not ok:
+  raise Exception("Error when computing MG-CADSurf mesh")
 
 # Make a copy of the 2D mesh
 mgtetraMesh_wo_geometry = smesh.CopyMesh( mgtetraMesh, 'MG-Tetra w/o geometry', 0, 0)
@@ -37,10 +48,19 @@ mgtetraMesh_wo_geometry = smesh.CopyMesh( mgtetraMesh, 'MG-Tetra w/o geometry',
 MG_Tetra = mgtetraMesh.Tetrahedron( smeshBuilder.MG_Tetra )
 MG_Tetra_Parameters = MG_Tetra.Parameters()
 # Create the enforced vertex
-MG_Tetra_Parameters.SetEnforcedVertex( 50, 50, 100, 2) # no group
+x1 = 50
+y1 = 50
+z1 = 100
+MG_Tetra_Parameters.SetEnforcedVertex( x1, y1, z1, 2) # no group
 # Compute the mesh
-mgtetraMesh.Compute()
+ok = mgtetraMesh.Compute()
+if not ok:
+  raise Exception("Error when computing MG_Tetra mesh")
+
+# Check that the enforced node is at the enforced coords
+node_closest = mgtetraMesh.FindNodeClosestTo(x1, y1, z1)
 
+assert isNodeAtCoord(node_closest, x1, y1, z1)
 
 # Ex2: Add one vertex enforced by a GEOM vertex at (50,50,100) 
 #      with physical size 5 and add it to a group called "My special nodes"
@@ -51,13 +71,26 @@ mgtetraMesh_wo_geometry.AddHypothesis( MG_Tetra )
 mgtetraMesh_wo_geometry.AddHypothesis( MG_Tetra_Parameters_wo_geometry )
 
 # Create the enforced vertex
-p1 = geompy.MakeVertex(150, 150, 100)
-geompy.addToStudy(p1, "p1")
-MG_Tetra_Parameters_wo_geometry.SetEnforcedVertexGeomWithGroup( p1, 5 , "My special nodes")
+x2 = 50
+y2 = 50
+z2 = 100
+p2 = geompy.MakeVertex(x2, y2, z2)
+geompy.addToStudy(p2, "p2")
+gr_enforced_name = "My special nodes"
+MG_Tetra_Parameters_wo_geometry.SetEnforcedVertexGeomWithGroup( p2, 5 , gr_enforced_name)
 #MG_Tetra_Parameters.SetEnforcedVertexGeom( p1, 5 ) # no group
 
 # compute the mesh
-mgtetraMesh_wo_geometry.Compute()
+ok = mgtetraMesh_wo_geometry.Compute()
+if not ok:
+  raise Exception("Error when computing MG_Tetra mesh without geometry")
+
+# Check that the enforced node is at the enforced coords
+gr_enforced_nodes = mgtetraMesh_wo_geometry.GetGroupByName(gr_enforced_name)[0]
+assert (gr_enforced_nodes.Size() == 1)
+node_enforced = gr_enforced_nodes.GetIDs()[0]
+
+assert isNodeAtCoord(node_enforced, x2, y2, z2)
 
 # Erase all enforced vertices
 MG_Tetra_Parameters.ClearEnforcedVertices()
index 34d9e6d58295af16c22e47946b93cbae654310d1..baa6d7feac6d186ca09c4a4a8b55250b2ebe0e3e 100644 (file)
@@ -2,10 +2,11 @@ import salome
 salome.salome_init()
 
 from salome.geom import geomBuilder
-geompy = geomBuilder.New(salome.myStudy)
+geompy = geomBuilder.New()
 
+import SMESH
 from salome.smesh import smeshBuilder
-smesh =  smeshBuilder.New(salome.myStudy)
+smesh =  smeshBuilder.New()
 
 # create a disk
 disk = geompy.MakeDiskR(100., 1, theName="disk")
@@ -14,7 +15,9 @@ disk = geompy.MakeDiskR(100., 1, theName="disk")
 mesh = smesh.Mesh( disk )
 cadsurf = mesh.Triangle( smeshBuilder.MG_CADSurf )
 cadsurf.SetQuadAllowed( True )
-mesh.Compute()
+ok = mesh.Compute()
+if not ok:
+  raise Exception("Error when computing mesh")
 
 # extrude the 2D mesh into a prismatic mesh
 mesh.ExtrusionSweepObject( mesh, [0,0,10], 7 )
@@ -32,10 +35,19 @@ mg_opt.SetSmoothOffSlivers( True )
 mg_opt.SetOptimizationLevel( smeshBuilder.Strong_Optimization )
 
 # run optimization
-optMesh.Compute()
+ok = optMesh.Compute()
+if not ok:
+  raise Exception("Error when computing optimization mesh")
 
 print("Nb tetra before optimization", mesh.NbTetras())
 print("Nb tetra after  optimization", optMesh.NbTetras())
 
+# Check that aspect ratio 3D of optimized mesh is better than original mesh
+min_aspectratio_orig, max_aspectratio_orig = mesh.GetMinMax(SMESH.FT_AspectRatio3D)
+min_aspectratio_optim, max_aspectratio_optim = optMesh.GetMinMax(SMESH.FT_AspectRatio3D)
+
+assert (min_aspectratio_orig - min_aspectratio_optim)/min_aspectratio_orig < 0.5
+assert (max_aspectratio_orig - max_aspectratio_optim)/max_aspectratio_orig < 0.5
+
 # End of script
 
index 46d098a84857754e80cdbc9b557bacbe41e7eb78..c25869235d2d0de6b5dea74a40d6dfd121d7dad0 100644 (file)
@@ -19,7 +19,15 @@ MG_CADSurf = mgtetraMesh.Triangle(algo=smeshBuilder.MG_CADSurf)
 MG_Tetra = mgtetraMesh.Tetrahedron(algo=smeshBuilder.MG_Tetra)
 
 # compute the mesh
-mgtetraMesh.Compute()
+ok = mgtetraMesh.Compute()
+
+if not ok:
+  raise Exception("Error when computing mgtetraMesh")
+
+volume = smesh.GetVolume(mgtetraMesh)
+expected_volume = 200**3
+
+assert (volume-expected_volume)/expected_volume < 1e-12
 
 # End of script