X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=test%2Ftest_import_1D2D_with_tetras_and_pyras.py;fp=test%2Ftest_import_1D2D_with_tetras_and_pyras.py;h=e1b4d07ddd0f0b0c4fb78082b8e0d6d7e20a0445;hp=0000000000000000000000000000000000000000;hb=e0c1dd6f368b19ac7e2f2841a096c06ef2f82eed;hpb=2f549d1c9a2a2791550d3f617354ade97d588a78 diff --git a/test/test_import_1D2D_with_tetras_and_pyras.py b/test/test_import_1D2D_with_tetras_and_pyras.py new file mode 100644 index 000000000..e1b4d07dd --- /dev/null +++ b/test/test_import_1D2D_with_tetras_and_pyras.py @@ -0,0 +1,179 @@ +#!/usr/bin/env python + +# test fix bos #33557 Bad pyramids generated + +import sys +import salome + + +#nb_segs_right = 30 +#nb_segs_right = 100 +nb_segs_right = 150 +#nb_segs_right = 200 + +algo = "MG-Tetra" +#algo = "Netgen" +#algo = "GMSH" + +algo_gmsh = "Delaunay" +#algo_gmsh = "HXT" + +salome.salome_init() + + + +### +### GEOM component +### + +import GEOM +from salome.geom import geomBuilder +import math +import SALOMEDS + + +geompy = geomBuilder.New() + +O = geompy.MakeVertex(0, 0, 0) +OX = geompy.MakeVectorDXDYDZ(1, 0, 0) +OY = geompy.MakeVectorDXDYDZ(0, 1, 0) +OZ = geompy.MakeVectorDXDYDZ(0, 0, 1) +Vertex_1 = geompy.MakeVertex(-0.62375, 0.0575, 0.02) +Vertex_2 = geompy.MakeVertex(-0.62375, -0.0575, 0.02) +Vertex_3 = geompy.MakeVertex(-0.365, -0.0575, 0.02) +Vertex_4 = geompy.MakeVertex(-0.365, 0.0575, 0.02) +geomObj_1 = geompy.MakeMarker(0, 0, 0, 1, 0, 0, 0, 1, 0) +Curve_1 = geompy.MakePolyline([Vertex_1, Vertex_4, Vertex_3, Vertex_2], True) +Face_1 = geompy.MakeFaceWires([Curve_1], 1) +[right_side] = geompy.SubShapes(Face_1, [6]) +[right_side] = geompy.GetExistingSubObjects(Face_1, False) +path = geompy.MakePrismDXDYDZ(Vertex_3, 0, 0, 0.1) +Translation_1 = geompy.MakeTranslation(Vertex_3, 0, -0.02, 0) +Translation_2 = geompy.MakeTranslation(Vertex_4, 0, 0.02, 0) +Translation_3 = geompy.MakeTranslation(Vertex_4, 0, 0.02, 0) +Line_1 = geompy.MakeLineTwoPnt(Translation_3, Translation_1) +Extrusion_1 = geompy.MakePrismDXDYDZ(Line_1, 0.2, 0, 0) +Partition_1 = geompy.MakePartition([Extrusion_1], [Face_1], [], [], geompy.ShapeType["FACE"], 0, [], 0) +Extrusion_2 = geompy.MakePrismDXDYDZ(Partition_1, 0, 0, 0.1) +[Face_to_enforce] = geompy.SubShapes(Extrusion_2, [30]) +geompy.addToStudy( O, 'O' ) +geompy.addToStudy( OX, 'OX' ) +geompy.addToStudy( OY, 'OY' ) +geompy.addToStudy( OZ, 'OZ' ) +geompy.addToStudy( Vertex_1, 'Vertex_1' ) +geompy.addToStudy( Vertex_2, 'Vertex_2' ) +geompy.addToStudy( Vertex_3, 'Vertex_3' ) +geompy.addToStudy( Vertex_4, 'Vertex_4' ) +geompy.addToStudy( Curve_1, 'Curve_1' ) +geompy.addToStudy( Face_1, 'Face_1' ) +geompy.addToStudyInFather( Face_1, right_side, 'right_side' ) +geompy.addToStudy( path, 'path' ) +geompy.addToStudy( Translation_2, 'Translation_2' ) +geompy.addToStudy( Translation_1, 'Translation_1' ) +geompy.addToStudy( Translation_3, 'Translation_3' ) +geompy.addToStudy( Line_1, 'Line_1' ) +geompy.addToStudy( Extrusion_1, 'Extrusion_1' ) +geompy.addToStudy( Partition_1, 'Partition_1' ) +geompy.addToStudy( Extrusion_2, 'Extrusion_2' ) +geompy.addToStudyInFather( Extrusion_2, Face_to_enforce, 'Face_to_enforce' ) + + +# Create a group with other faces than Face_to_enforce +all_faces = geompy.SubShapeAll(Extrusion_2, geompy.ShapeType["FACE"]) +gr_other_faces = geompy.CreateGroup(Extrusion_2, geompy.ShapeType["FACE"]) +geompy.UnionList(gr_other_faces, all_faces) +geompy.DifferenceList(gr_other_faces, [Face_to_enforce]) +geompy.addToStudyInFather( Extrusion_2, gr_other_faces, 'other_faces' ) + +geom_volume = geompy.BasicProperties(Extrusion_2)[2] + +### +### SMESH component +### + +import SMESH, SALOMEDS +from salome.smesh import smeshBuilder + +smesh = smeshBuilder.New() + +Mesh_1 = smesh.Mesh(Face_1,'Mesh_1') + +Regular_1D = Mesh_1.Segment() +Number_of_Segments_1 = Regular_1D.NumberOfSegments(10) + +Quadrangle_2D = Mesh_1.Quadrangle(algo=smeshBuilder.QUADRANGLE) + +right_side_1 = Mesh_1.GroupOnGeom(right_side,'right_side',SMESH.EDGE) +Regular_1D_1 = Mesh_1.Segment(geom=right_side) + + +Number_of_Segments_2 = Regular_1D_1.NumberOfSegments(nb_segs_right) +Propagation_of_1D_Hyp = Regular_1D_1.Propagation() + +isDone = Mesh_1.Compute() + +Mesh_path = smesh.Mesh(path,'Mesh_path') + +Regular_1D_2 = Mesh_path.Segment() +Local_Length_1 = Regular_1D_2.LocalLength(0.01,None,1e-07) + +isDone = Mesh_path.Compute() + +([ right_side_extruded, right_side_top ], error) = Mesh_1.ExtrusionAlongPathObjects( [], [ Mesh_1 ], [ Mesh_1 ], Mesh_path, None, 1, 0, [ ], 0, 0, [ 0, 0, 0 ], 1, [ ], 0 ) + +Mesh_2 = smesh.Mesh(Extrusion_2,'Mesh_2') + +# Enforce elements from Mesh_1 +Import_1D2D = Mesh_2.UseExisting2DElements(geom=Face_to_enforce) +Source_Faces_1 = Import_1D2D.SourceFaces([ right_side_extruded ],0,0) + +if algo != "GMSH": + # CADSurf is global mesh + MG_CADSurf = Mesh_2.Triangle(algo=smeshBuilder.MG_CADSurf) +else: + # CADSurf is a submesh on other faces + MG_CADSurf = Mesh_2.Triangle(algo=smeshBuilder.MG_CADSurf, geom=gr_other_faces) +MG_CADSurf_Parameters_1 = MG_CADSurf.Parameters() +MG_CADSurf_Parameters_1.SetGeometricMesh( 0 ) +MG_CADSurf_Parameters_1.SetPhySize( 0.005 ) +MG_CADSurf_Parameters_1.SetMinSize( 0.005 ) +MG_CADSurf_Parameters_1.SetMaxSize( 0.005 ) +MG_CADSurf_Parameters_1.SetCorrectSurfaceIntersection( False ) +MG_CADSurf_Parameters_1.SetElementType( 1 ) # Quadrangle preference +if algo == "MG-Tetra": + Mesh_2.Tetrahedron(algo=smeshBuilder.MG_Tetra) +elif algo == "Netgen": + Mesh_2.Tetrahedron() +elif algo == "GMSH": + GMSH = Mesh_2.Tetrahedron(algo=smeshBuilder.GMSH) + Gmsh_Parameters = GMSH.Parameters() + Gmsh_Parameters.Set2DAlgo( 0 ) + Gmsh_Parameters.SetMinSize( 0 ) + Gmsh_Parameters.SetMaxSize( 0.005 ) + if algo_gmsh == "HXT": + Gmsh_Parameters.Set3DAlgo( 4 ) + else: + Gmsh_Parameters.Set3DAlgo( 0 ) + Gmsh_Parameters.SetIs2d( 0 ) + +Face_to_enforce_1 = Mesh_2.GroupOnGeom(Face_to_enforce,'Face_to_enforce',SMESH.FACE) + + +isDone = Mesh_2.Compute() + +Mesh_2.MakeGroup("pyramids", SMESH.VOLUME, CritType=SMESH.FT_ElemGeomType, Threshold=SMESH.Geom_PYRAMID) + +if salome.sg.hasDesktop(): + salome.sg.updateObjBrowser() + +if not isDone: + raise Exception("Error on mesh compute") + +mesh_volume = Mesh_2.GetVolume() +assert abs(mesh_volume-geom_volume)/geom_volume < 1e-7, "Wrong mesh volume" + +min_aspect_ratio, max_aspect_ratio = Mesh_2.GetMinMax(SMESH.FT_AspectRatio3D) +assert max_aspect_ratio < 200, "Bad aspect ratio 3D: %.1f"%max_aspect_ratio + +min_volume, max_volume = Mesh_2.GetMinMax(SMESH.FT_Volume3D) +assert min_volume > 1e-10, "Bad min volume: %s"%min_volume