From 27653b3419580d3711578991ae0d4bc8bdda25b3 Mon Sep 17 00:00:00 2001 From: vsr Date: Wed, 5 Nov 2008 10:13:21 +0000 Subject: [PATCH] Add new example python script (from CEA) --- src/SMESH_SWIG/ex29_refine.py | 205 ++++++++++++++++++++++++++++++++++ 1 file changed, 205 insertions(+) create mode 100644 src/SMESH_SWIG/ex29_refine.py diff --git a/src/SMESH_SWIG/ex29_refine.py b/src/SMESH_SWIG/ex29_refine.py new file mode 100644 index 000000000..448cbf833 --- /dev/null +++ b/src/SMESH_SWIG/ex29_refine.py @@ -0,0 +1,205 @@ +# CEA/LGLS-LETR 2008, Francis KLOSS (OCC) +# ======================================= + +# Procedure that take a triangulation and split all triangles in 4 others triangles + +import geompy +import smesh + +# Values +# ------ + +# Path for ".med" files +path = "/tmp/ex29_" + +# Name of the shape and the mesh +name = "Carre" + +# Add a node and needed edges +# --------------------------- + +def node(m, f, n1, n2, lnv): + x1, y1, z1 = m.GetNodeXYZ(n1) + x2, y2, z2 = m.GetNodeXYZ(n2) + + x = (x1 + x2) / 2.0 + y = (y1 + y2) / 2.0 + z = (z1 + z2) / 2.0 + + i = m.AddNode(x, y, z) + + in1 = m.GetShapeID(n1) + in2 = m.GetShapeID(n2) + + if (in1==f) or (in2==f): + m.SetNodeOnFace(i, f, 0, 0) + + else: + e1 = m.AddEdge([ n1, i ]) + e2 = m.AddEdge([ i , n2 ]) + + if n1 in lnv: + e = in2 + else: + e = in1 + + m.SetMeshElementOnShape(e1, e) + m.SetMeshElementOnShape(e2, e) + m.SetNodeOnEdge(i, e, 0) + + return i + +# Add a triangle and associate to the CAD face +# -------------------------------------------- + +def triangle(m, f, n1, n2, n3): + i = m.AddFace([ n1, n2, n3 ]) + m.SetMeshElementOnShape(i, f) + +# Split all triangles in 4 triangles +# ---------------------------------- + +def SplitTrianglesIn4(m): + # Get all triangles + triangles = m.GetElementsByType(smesh.FACE) + + # Remove all edges + m.RemoveElements(m.GetElementsByType(smesh.EDGE)) + + # Get the list of nodes (ids) associated with the CAD vertices + shape = m.GetShape() + lnv = [] + for v in geompy.SubShapeAll(shape, geompy.ShapeType["VERTEX"]): + lnv = lnv + m.GetSubMeshNodesId(v, True) + + # Split every triangle + for t in triangles: + noeud_1, noeud_2, noeud_3 = m.GetElemNodes(t) + + face = m.GetShapeIDForElem(t) + + noeud_12 = node(m, face, noeud_1, noeud_2, lnv) + noeud_23 = node(m, face, noeud_2, noeud_3, lnv) + noeud_13 = node(m, face, noeud_1, noeud_3, lnv) + + triangle(m, face, noeud_1 , noeud_12, noeud_13) + triangle(m, face, noeud_2 , noeud_23, noeud_12) + triangle(m, face, noeud_3 , noeud_13, noeud_23) + triangle(m, face, noeud_12, noeud_23, noeud_13) + + # Remove all initial triangles + m.RemoveElements(triangles) + + # Merge all identical nodes + m.MergeNodes(m.FindCoincidentNodes(0.0001)) + +# Build a CAD square +# ------------------ + +x0 = 0.0 ; y0 = 0.0 ; z0 = 0.0 +x1 = 1.0 ; y1 = 0.0 ; z1 = 0.0 +x2 = 1.0 ; y2 = 1.0 ; z2 = 0.0 +x3 = 0.0 ; y3 = 1.0 ; z3 = 0.0 + +P0 = geompy.MakeVertex(x0, y0, z0) +P1 = geompy.MakeVertex(x1, y1, z1) +P2 = geompy.MakeVertex(x2, y2, z2) +P3 = geompy.MakeVertex(x3, y3, z3) + +square = geompy.MakeQuad4Vertices(P0, P1, P2, P3) +geompy.addToStudy(square, name) + +# Refine edges and create group of mesh +# ------------------------------------- + +def refine(m, p1, p2, n, k, name): + s = m.GetShape() + + g = geompy.CreateGroup(s, geompy.ShapeType["EDGE"]) + e = geompy.GetEdge(s, p1, p2) + i = geompy.GetSubShapeID(s, e) + geompy.AddObject(g, i) + m.Group(g, name) + + a = m.Segment(e) + a.NumberOfSegments(n, k) + +# Mesh the square +# --------------- + +MyMesh = smesh.Mesh(square) + +refine(MyMesh, P1, P2, 8, 7, "Droite") +refine(MyMesh, P3, P0, 9, 10, "Gauche") +refine(MyMesh, P0, P1, 7, 9, "Bas" ) +refine(MyMesh, P2, P3, 12, 14, "Haut" ) + +algo2D = MyMesh.Triangle() +algo2D.MaxElementArea(0.07) + +MyMesh.Compute() + +MyMesh.ExportMED(path+"110_triangles.med", 0) + +# Disturb the mesh +# ---------------- + +MyMesh.MoveNode( 37, 0.05 , 0.368967 , 0 ) +MyMesh.MoveNode( 38, 0.34 , 0.0762294, 0 ) +MyMesh.MoveNode( 40, 0.8 , 0.42 , 0 ) +MyMesh.MoveNode( 42, 0.702662, 0.74 , 0 ) +MyMesh.MoveNode( 46, 0.4 , 0.374656 , 0 ) +MyMesh.MoveNode( 47, 0.13 , 0.63 , 0 ) +MyMesh.MoveNode( 49, 0.222187, 0.3 , 0 ) +MyMesh.MoveNode( 54, 0.557791, 0.05 , 0 ) +MyMesh.MoveNode( 55, 0.7 , 0.2 , 0 ) +MyMesh.MoveNode( 56, 0.73 , 0.52 , 0 ) +MyMesh.MoveNode( 58, 0.313071, 0.31 , 0 ) +MyMesh.MoveNode( 59, 0.8 , 0.56 , 0 ) +MyMesh.MoveNode( 62, 0.592703, 0.95 , 0 ) +MyMesh.MoveNode( 63, 0.28 , 0.5 , 0 ) +MyMesh.MoveNode( 65, 0.49 , 0.93 , 0 ) +MyMesh.MoveNode( 68, 0.501038, 0.65 , 0 ) +MyMesh.MoveNode( 69, 0.37 , 0.63 , 0 ) +MyMesh.MoveNode( 70, 0.597025, 0.52 , 0 ) +MyMesh.MoveNode( 72, 0.899 , 0.878589 , 0 ) +MyMesh.MoveNode( 73, 0.92 , 0.85 , 0 ) +MyMesh.MoveNode( 74, 0.820851, 0.75 , 0 ) + +NbCells1 = 110 +MyMesh.ExportMED(path+"110_triangles_2.med", 0) + +# First mesh refining +# ------------------- + +SplitTrianglesIn4(MyMesh) + +NbCells2 = NbCells1*4 +print("Mesh with "+str(NbCells2)+" cells computed.") + +MyMesh.ExportMED(path+str(NbCells2)+"_triangles.med", 0) + +# Second mesh refining +# -------------------- + +SplitTrianglesIn4(MyMesh) + +NbCells3 = NbCells2*4 +print("Mesh with "+str(NbCells3)+" cells computed.") + +MyMesh.ExportMED(path+str(NbCells3)+"_triangles.med",0) + +# Third mesh refining +# ------------------- + +SplitTrianglesIn4(MyMesh) + +NbCells4 = NbCells3*4 +print("Mesh with "+str(NbCells4)+" cells computed.") + +MyMesh.ExportMED(path+str(NbCells4)+"_triangles.med", 0) + +# Update the object browser +# ------------------------- + +salome.sg.updateObjBrowser(1) -- 2.30.2