1 # -*- coding: utf-8 -*-
2 # Copyright (C) 2018-2019 CEA/DEN, EDF R&D
4 # This library is free software; you can redistribute it and/or
5 # modify it under the terms of the GNU Lesser General Public
6 # License as published by the Free Software Foundation; either
7 # version 2.1 of the License, or (at your option) any later version.
9 # This library is distributed in the hope that it will be useful,
10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 # Lesser General Public License for more details.
14 # You should have received a copy of the GNU Lesser General Public
15 # License along with this library; if not, write to the Free Software
16 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
31 from salome.geom import geomBuilder
36 geompy = geomBuilder.New()
48 O = geompy.MakeVertex(0, 0, 0)
49 OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
50 OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
51 OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
52 geompy.addToStudy( O, 'O' )
53 geompy.addToStudy( OX, 'OX' )
54 geompy.addToStudy( OY, 'OY' )
55 geompy.addToStudy( OZ, 'OZ' )
57 Cylinder_1 = geompy.MakeCylinderRH(r1, h1)
58 Cylinder_2 = geompy.MakeCylinderRH(r2, h2)
59 Rotation_1 = geompy.MakeRotation(Cylinder_2, OY, -90*math.pi/180.0)
60 Translation_1 = geompy.MakeTranslation(Rotation_1, 0, 0, length_piquage)
62 piquage = geompy.MakeFuseList([Cylinder_1, Translation_1], True, True)
63 geompy.addToStudy( piquage, 'piquage' )
65 Inlet_z = geompy.GetFaceNearPoint(piquage, O)
66 geompy.addToStudyInFather( piquage, Inlet_z, 'Inlet_z' )
68 p_inlet_x = geompy.MakeVertex(-h2, 0, length_piquage)
69 Inlet_x = geompy.GetFaceNearPoint(piquage, p_inlet_x)
70 geompy.addToStudyInFather( piquage, Inlet_x, 'Inlet_x' )
72 p_outlet = geompy.MakeVertex(0, 0, h1)
73 Outlet = geompy.GetFaceNearPoint(piquage, p_outlet)
74 geompy.addToStudyInFather( piquage, Outlet, 'Outlet' )
76 Wall = geompy.CreateGroup(piquage, geompy.ShapeType["FACE"])
77 faces = geompy.SubShapeAll(piquage, geompy.ShapeType["FACE"])
78 geompy.UnionList(Wall, faces)
79 geompy.DifferenceList(Wall, [Inlet_x, Inlet_z, Outlet])
80 geompy.addToStudyInFather( piquage, Wall, 'Wall' )
82 p_corner = geompy.MakeVertex(-r2, 0, length_piquage+r2)
83 corner = geompy.GetVertexNearPoint(piquage, p_corner)
84 geompy.addToStudyInFather( piquage, corner, 'corner' )
86 geom_groups = [Inlet_x, Inlet_z, Outlet, Wall]
93 import SMESH, SALOMEDS
94 from salome.smesh import smeshBuilder
96 from salome.StdMeshers import StdMeshersBuilder
97 from salome.HYBRIDPlugin import HYBRIDPluginBuilder
99 smesh = smeshBuilder.New()
101 # Viscous layers with imprinting
102 # ==============================
104 Mesh_1 = smesh.Mesh(piquage, "Mesh_with_imprinting")
106 NETGEN_2D = Mesh_1.Triangle(algo=smeshBuilder.NETGEN_1D2D)
107 NETGEN_2D_Parameters = NETGEN_2D.Parameters()
109 NETGEN_2D_Parameters.SetMinSize( 0.01 )
110 NETGEN_2D_Parameters.SetLocalSizeOnShape(corner, 0.01)
111 NETGEN_2D_Parameters.SetFineness( 5 )
112 NETGEN_2D_Parameters.SetGrowthRate( 0.1 )
113 NETGEN_2D_Parameters.SetNbSegPerEdge( 2 )
114 NETGEN_2D_Parameters.SetNbSegPerRadius( 3 )
116 wall_faces = geompy.SubShapeAll(Wall, geompy.ShapeType["FACE"])
117 wall_ids = geompy.GetSubShapesIDs(piquage, wall_faces)
118 Inlet_x_id, Inlet_z_id, Outlet_id = geompy.GetSubShapesIDs(piquage, [Inlet_x, Inlet_z, Outlet])
120 MG_Hybrid_2 = Mesh_1.Tetrahedron(algo=smeshBuilder.HYBRID)
121 MG_Hybrid_Parameters_2 = MG_Hybrid_2.Parameters()
122 MG_Hybrid_Parameters_2.SetLayersOnAllWrap( 0 )
123 MG_Hybrid_Parameters_2.SetElementGeneration( HYBRIDPluginBuilder.Generation_Tetra_Dominant )
124 MG_Hybrid_Parameters_2.SetHeightFirstLayer( 0.01 )
125 MG_Hybrid_Parameters_2.SetBoundaryLayersProgression( 1.1 )
126 MG_Hybrid_Parameters_2.SetNbOfBoundaryLayers( 3 )
127 MG_Hybrid_Parameters_2.SetFacesWithLayers( wall_ids )
128 MG_Hybrid_Parameters_2.SetFacesWithImprinting( [ Inlet_x_id, Inlet_z_id, Outlet_id ] )
130 isDone = Mesh_1.Compute()
133 raise Exception("Error when computing Mesh_with_imprinting")
135 assert Mesh_1.NbQuadrangles() > 0
139 for geom_group in geom_groups:
140 name = geom_group.GetName()
141 gr = Mesh_1.Group(geom_group)
142 d_groups_1[name] = gr
143 d_geom_groups[name] = geom_group
145 faces_imprinted = ["Inlet_x", "Inlet_z", "Outlet"]
147 # Check viscous layers with imprinting
148 for name in faces_imprinted:
149 geom_area = geompy.BasicProperties(d_geom_groups[name])[1]
150 gr_1 = d_groups_1[name]
151 #assert gr_1.Size() > 0
152 assert gr_1.Size() > 0
153 # Nb of quadrangles is in 7th index of mesh info
154 assert gr_1.GetMeshInfo()[7] > 0
156 # Compare mesh group and geom group
157 area_error_1 = abs(geom_area-smesh.GetArea(gr_1))/geom_area
158 assert area_error_1 < 0.025
160 # Compare whole mesh volume
161 shape_volume = geompy.BasicProperties(piquage)[2]
162 mesh_1_volume = Mesh_1.GetVolume()
163 volume_error_1 = abs(shape_volume-mesh_1_volume)/shape_volume
165 assert volume_error_1 < 0.02
170 Mesh_snap = smesh.Mesh(piquage, "Mesh_with_snapping")
172 HYBRID_3D = Mesh_snap.Tetrahedron(algo=smeshBuilder.HYBRID)
173 HYBRID_Parameters_2 = HYBRID_3D.Parameters()
174 HYBRID_Parameters_2.SetElementGeneration( 0 )
175 HYBRID_Parameters_2.SetHeightFirstLayer( 0.01 )
176 HYBRID_Parameters_2.SetBoundaryLayersProgression( 1.1 )
177 HYBRID_Parameters_2.SetNbOfBoundaryLayers( 3 )
178 HYBRID_Parameters_2.SetFacesWithLayers( wall_ids )
179 HYBRID_Parameters_2.SetFacesWithSnapping([ Inlet_x_id, Inlet_z_id, Outlet_id ])
180 HYBRID_Parameters_2.SetLayersOnAllWrap( 0 )
182 Import_1D2D = Mesh_snap.UseExisting2DElements(geom=Inlet_x)
183 Source_Faces_1 = Import_1D2D.SourceFaces([ d_groups_1["Inlet_x"] ],0,0)
185 Import_1D2D_1 = Mesh_snap.UseExisting2DElements(geom=Inlet_z)
186 Source_Faces_2 = Import_1D2D_1.SourceFaces([ d_groups_1["Inlet_z"] ],0,0)
188 Import_1D2D_2 = Mesh_snap.UseExisting2DElements(geom=Outlet)
189 Source_Faces_3 = Import_1D2D_2.SourceFaces([ d_groups_1["Outlet"] ],0,0)
191 Import_1D2D_3 = Mesh_snap.UseExisting2DElements(geom=Wall)
192 Source_Faces_4 = Import_1D2D_3.SourceFaces([ d_groups_1["Wall"] ],0,0)
194 isDone = Mesh_snap.Compute()
198 for geom_group in geom_groups:
199 name = geom_group.GetName()
200 gr = Mesh_snap.Group(geom_group)
201 d_groups_2[name] = gr
203 # Compare whole mesh volume
204 mesh_2_volume = Mesh_snap.GetVolume()
205 volume_error_2 = abs(mesh_1_volume-mesh_2_volume)/mesh_1_volume
207 assert volume_error_2 < 1e-4
209 # Check viscous layers with snapping
210 for name in faces_imprinted:
211 gr_1 = d_groups_1[name]
212 gr_2 = d_groups_2[name]
213 assert gr_2.Size() > 0
214 # Nb of quadrangles is in 7th index of mesh info
215 assert gr_1.GetMeshInfo()[7] == gr_1.GetMeshInfo()[7]
217 # Compare mesh group with imprinting set by ids and mesh group with imprinting set by geom group
218 area_error_2 = abs(smesh.GetArea(gr_1)-smesh.GetArea(gr_2))/smesh.GetArea(gr_1)
219 assert area_error_2 < 1e-10
224 Mesh_snap_on_group = smesh.Mesh(piquage, "Mesh_with_snapping_on_groups")
226 HYBRID_3D_on_group = Mesh_snap_on_group.Tetrahedron(algo=smeshBuilder.HYBRID)
227 HYBRID_Parameters_3 = HYBRID_3D_on_group.Parameters()
228 HYBRID_Parameters_3.SetElementGeneration( 0 )
229 HYBRID_Parameters_3.SetHeightFirstLayer( 0.01 )
230 HYBRID_Parameters_3.SetBoundaryLayersProgression( 1.1 )
231 HYBRID_Parameters_3.SetNbOfBoundaryLayers( 3 )
232 HYBRID_Parameters_3.SetFacesWithLayers( wall_ids )
233 HYBRID_3D_on_group.SetFacesWithSnapping([ Inlet_x, Inlet_z, Outlet ])
234 HYBRID_Parameters_3.SetLayersOnAllWrap( 0 )
236 Import_1D2D = Mesh_snap_on_group.UseExisting2DElements(geom=Inlet_x)
237 Source_Faces_1 = Import_1D2D.SourceFaces([ d_groups_1["Inlet_x"] ],0,0)
239 Import_1D2D_1 = Mesh_snap_on_group.UseExisting2DElements(geom=Inlet_z)
240 Source_Faces_2 = Import_1D2D_1.SourceFaces([ d_groups_1["Inlet_z"] ],0,0)
242 Import_1D2D_2 = Mesh_snap_on_group.UseExisting2DElements(geom=Outlet)
243 Source_Faces_3 = Import_1D2D_2.SourceFaces([ d_groups_1["Outlet"] ],0,0)
245 Import_1D2D_3 = Mesh_snap_on_group.UseExisting2DElements(geom=Wall)
246 Source_Faces_4 = Import_1D2D_3.SourceFaces([ d_groups_1["Wall"] ],0,0)
248 isDone = Mesh_snap_on_group.Compute()
252 for geom_group in geom_groups:
253 name = geom_group.GetName()
254 gr = Mesh_snap_on_group.Group(geom_group)
255 d_groups_3[name] = gr
257 # Compare whole mesh volume
258 mesh_3_volume = Mesh_snap_on_group.GetVolume()
259 volume_error_3 = abs(mesh_2_volume-mesh_3_volume)/mesh_2_volume
261 assert volume_error_3 < 1e-10
263 # Check viscous layers with snapping
264 for name in faces_imprinted:
265 gr_2 = d_groups_2[name]
266 gr_3 = d_groups_3[name]
267 assert gr_3.Size() > 0
268 # Nb of quadrangles is in 7th index of mesh info
269 assert gr_2.GetMeshInfo()[7] == gr_3.GetMeshInfo()[7]
271 # Compare mesh group with imprinting set by ids and mesh group with imprinting set by geom group
272 area_error_3 = abs(smesh.GetArea(gr_2)-smesh.GetArea(gr_3))/smesh.GetArea(gr_2)
273 assert area_error_3 < 1e-10