tests/layers_imprinting.py

   1 # -*- coding: utf-8 -*-
   2
   3 import sys
   4 import salome
   5
   6 salome.salome_init()
   7 theStudy = salome.myStudy
   8
   9 ###
  10 ### GEOM component
  11 ###
  12
  13 import GEOM
  14 from salome.geom import geomBuilder
  15 import math
  16 import SALOMEDS
  17
  18
  19 geompy = geomBuilder.New(theStudy)
  20
  21 # first cylinder
  22 r1 = 0.5
  23 h1 = 5
  24
  25 # second cylinder
  26 r2 = 0.3
  27 h2 = 3
  28
  29 length_piquage = 1.5
  30
  31 O = geompy.MakeVertex(0, 0, 0)
  32 OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
  33 OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
  34 OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
  35 geompy.addToStudy( O, 'O' )
  36 geompy.addToStudy( OX, 'OX' )
  37 geompy.addToStudy( OY, 'OY' )
  38 geompy.addToStudy( OZ, 'OZ' )
  39
  40 Cylinder_1 = geompy.MakeCylinderRH(r1, h1)
  41 Cylinder_2 = geompy.MakeCylinderRH(r2, h2)
  42 Rotation_1 = geompy.MakeRotation(Cylinder_2, OY, -90*math.pi/180.0)
  43 Translation_1 = geompy.MakeTranslation(Rotation_1, 0, 0, length_piquage)
  44
  45 piquage = geompy.MakeFuseList([Cylinder_1, Translation_1], True, True)
  46 geompy.addToStudy( piquage, 'piquage' )
  47
  48 Inlet_z = geompy.GetFaceNearPoint(piquage, O)
  49 geompy.addToStudyInFather( piquage, Inlet_z, 'Inlet_z' )
  50
  51 p_inlet_x = geompy.MakeVertex(-h2, 0, length_piquage)
  52 Inlet_x = geompy.GetFaceNearPoint(piquage, p_inlet_x)
  53 geompy.addToStudyInFather( piquage, Inlet_x, 'Inlet_x' )
  54
  55 p_outlet = geompy.MakeVertex(0, 0, h1)
  56 Outlet = geompy.GetFaceNearPoint(piquage, p_outlet)
  57 geompy.addToStudyInFather( piquage, Outlet, 'Outlet' )
  58
  59 Wall = geompy.CreateGroup(piquage, geompy.ShapeType["FACE"])
  60 faces = geompy.SubShapeAll(piquage, geompy.ShapeType["FACE"])
  61 geompy.UnionList(Wall, faces)
  62 geompy.DifferenceList(Wall, [Inlet_x, Inlet_z, Outlet])
  63 geompy.addToStudyInFather( piquage, Wall, 'Wall' )
  64
  65 p_corner = geompy.MakeVertex(-r2, 0, length_piquage+r2)
  66 corner = geompy.GetVertexNearPoint(piquage, p_corner)
  67 geompy.addToStudyInFather( piquage, corner, 'corner' )
  68
  69 geom_groups = [Inlet_x, Inlet_z, Outlet, Wall]
  70 d_geom_groups = {}
  71
  72 ###
  73 ### SMESH component
  74 ###
  75
  76 import  SMESH, SALOMEDS
  77 from salome.smesh import smeshBuilder
  78
  79 from salome.StdMeshers import StdMeshersBuilder
  80 from salome.HYBRIDPlugin import HYBRIDPluginBuilder
  81
  82 smesh = smeshBuilder.New(theStudy)
  83
  84
  85 # Viscous layers without imprinting
  86 # =================================
  87
  88 Mesh_1 = smesh.Mesh(piquage, "Mesh_without_imprinting")
  89 NETGEN_2D = Mesh_1.Triangle(algo=smeshBuilder.NETGEN_1D2D)
  90 NETGEN_2D_Parameters = NETGEN_2D.Parameters()
  91
  92 NETGEN_2D_Parameters.SetMinSize( 0.01 )
  93 NETGEN_2D_Parameters.SetLocalSizeOnShape(corner, 0.01)
  94 NETGEN_2D_Parameters.SetFineness( 5 )
  95 NETGEN_2D_Parameters.SetGrowthRate( 0.1 )
  96 NETGEN_2D_Parameters.SetNbSegPerEdge( 2 )
  97 NETGEN_2D_Parameters.SetNbSegPerRadius( 3 )
  98
  99 wall_faces = geompy.SubShapeAll(Wall, geompy.ShapeType["FACE"])
 100 wall_ids = geompy.GetSubShapesIDs(piquage, wall_faces)
 101
 102 MG_Hybrid = Mesh_1.Tetrahedron(algo=smeshBuilder.HYBRID)
 103 MG_Hybrid_Parameters_1 = MG_Hybrid.Parameters()
 104 MG_Hybrid_Parameters_1.SetLayersOnAllWrap( 0 )
 105 MG_Hybrid_Parameters_1.SetElementGeneration( HYBRIDPluginBuilder.Generation_Tetra_Dominant )
 106 MG_Hybrid_Parameters_1.SetHeightFirstLayer( 0.01 )
 107 MG_Hybrid_Parameters_1.SetBoundaryLayersProgression( 1.1 )
 108 MG_Hybrid_Parameters_1.SetNbOfBoundaryLayers( 3 )
 109 MG_Hybrid_Parameters_1.SetFacesWithLayers( wall_ids )
 110
 111
 112 isDone = Mesh_1.Compute()
 113
 114 if not isDone:
 115   raise Exception("Error when computing Mesh_without_imprinting")
 116
 117 d_groups_1 = {}
 118
 119 for geom_group in geom_groups:
 120   name = geom_group.GetName()
 121   gr = Mesh_1.Group(geom_group)
 122   d_groups_1[name] = gr
 123   d_geom_groups[name] = geom_group
 124
 125 assert Mesh_1.NbQuadrangles() == 0
 126
 127 shape_volume = geompy.BasicProperties(piquage)[2]
 128
 129 # Compare whole mesh volume
 130 volume_error_1 = abs(shape_volume-Mesh_1.GetVolume())/shape_volume
 131
 132 assert volume_error_1 < 0.02
 133
 134 # Viscous layers with imprinting
 135 # ==============================
 136
 137 Mesh_2 = smesh.Mesh(piquage, "Mesh_with_imprinting")
 138
 139 NETGEN_2D_1_1 = Mesh_2.Triangle(algo=smeshBuilder.NETGEN_1D2D)
 140 status = Mesh_2.AddHypothesis(NETGEN_2D_Parameters)
 141
 142 Inlet_x_id, Inlet_z_id, Outlet_id = geompy.GetSubShapesIDs(piquage, [Inlet_x, Inlet_z, Outlet])
 143
 144 MG_Hybrid_2 = Mesh_2.Tetrahedron(algo=smeshBuilder.HYBRID)
 145 MG_Hybrid_Parameters_2 = MG_Hybrid_2.Parameters()
 146 MG_Hybrid_Parameters_2.SetLayersOnAllWrap( 0 )
 147 MG_Hybrid_Parameters_2.SetElementGeneration( HYBRIDPluginBuilder.Generation_Tetra_Dominant )
 148 MG_Hybrid_Parameters_2.SetHeightFirstLayer( 0.01 )
 149 MG_Hybrid_Parameters_2.SetBoundaryLayersProgression( 1.1 )
 150 MG_Hybrid_Parameters_2.SetNbOfBoundaryLayers( 3 )
 151 MG_Hybrid_Parameters_2.SetFacesWithLayers( wall_ids )
 152 MG_Hybrid_Parameters_2.SetFacesWithImprinting( [ Inlet_x_id, Inlet_z_id, Outlet_id ] )
 153
 154 isDone = Mesh_2.Compute()
 155
 156 if not isDone:
 157   raise Exception("Error when computing Mesh_with_imprinting")
 158
 159 assert Mesh_2.NbQuadrangles() > 0
 160
 161 d_groups_2 = {}
 162
 163 for geom_group in geom_groups:
 164   name = geom_group.GetName()
 165   gr = Mesh_2.Group(geom_group)
 166   d_groups_2[name] = gr
 167
 168 faces_imprinted = ["Inlet_x", "Inlet_z", "Outlet"]
 169
 170 # Check viscous layers with imprinting
 171 for name in faces_imprinted:
 172   geom_area = geompy.BasicProperties(d_geom_groups[name])[1]
 173   gr_1 = d_groups_1[name]
 174   gr_2 = d_groups_2[name]
 175   #assert gr_1.Size() > 0
 176   assert gr_2.Size() > 0
 177   # Nb of quadrangles is in 7th index of mesh info
 178   assert gr_2.GetMeshInfo()[7] > 0
 179
 180   # Compare mesh group and geom group
 181   area_error_1 = abs(geom_area-smesh.GetArea(gr_2))/geom_area
 182   assert area_error_1 < 0.025
 183   # Compare mesh group with imprinting and mesh group without imprinting
 184   area_error_2 = abs(smesh.GetArea(gr_1)-smesh.GetArea(gr_2))/smesh.GetArea(gr_1)
 185   assert area_error_2 < 1e-10
 186
 187 # Compare whole mesh volume
 188 mesh_2_volume = Mesh_2.GetVolume()
 189 volume_error_2 = abs(shape_volume-mesh_2_volume)/shape_volume
 190
 191 assert volume_error_2 < 0.02
 192
 193 # Viscous layers with imprinting set by groups
 194 # ============================================
 195
 196 Mesh_3 = smesh.Mesh(piquage, "Mesh_with_imprinting_set_by_groups")
 197
 198 NETGEN_2D_1_1 = Mesh_3.Triangle(algo=smeshBuilder.NETGEN_1D2D)
 199 status = Mesh_3.AddHypothesis(NETGEN_2D_Parameters)
 200
 201 MG_Hybrid_3 = Mesh_3.Tetrahedron(algo=smeshBuilder.HYBRID)
 202 MG_Hybrid_Parameters_3 = MG_Hybrid_3.Parameters()
 203 MG_Hybrid_Parameters_3.SetLayersOnAllWrap( 0 )
 204 MG_Hybrid_Parameters_3.SetElementGeneration( HYBRIDPluginBuilder.Generation_Tetra_Dominant )
 205 MG_Hybrid_Parameters_3.SetHeightFirstLayer( 0.01 )
 206 MG_Hybrid_Parameters_3.SetBoundaryLayersProgression( 1.1 )
 207 MG_Hybrid_Parameters_3.SetNbOfBoundaryLayers( 3 )
 208 MG_Hybrid_3.SetFacesWithLayers( Wall )
 209 MG_Hybrid_3.SetFacesWithImprinting( [ Inlet_x, Inlet_z, Outlet ] )
 210
 211 isDone = Mesh_3.Compute()
 212
 213 if not isDone:
 214   raise Exception("Error when computing Mesh_with_imprinting_set_by_groups")
 215
 216 d_groups_3 = {}
 217
 218 for geom_group in geom_groups:
 219   name = geom_group.GetName()
 220   gr = Mesh_3.Group(geom_group)
 221   d_groups_3[name] = gr
 222
 223 # Compare whole mesh volume
 224 mesh_3_volume = Mesh_3.GetVolume()
 225 volume_error_3 = abs(mesh_2_volume-mesh_3_volume)/mesh_2_volume
 226
 227 assert abs(mesh_3_volume-mesh_2_volume)/mesh_2_volume < 1e-10
 228
 229 # Check viscous layers with imprinting
 230 for name in faces_imprinted:
 231   gr_2 = d_groups_2[name]
 232   gr_3 = d_groups_3[name]
 233   assert gr_3.Size() > 0
 234   # Nb of quadrangles is in 7th index of mesh info
 235   assert gr_2.GetMeshInfo()[7] == gr_3.GetMeshInfo()[7]
 236
 237   # Compare mesh group with imprinting set by ids and mesh group with imprinting set by geom group
 238   area_error_3 = abs(smesh.GetArea(gr_2)-smesh.GetArea(gr_3))/smesh.GetArea(gr_2)
 239   assert area_error_3 < 1e-10
 240
 241 if salome.sg.hasDesktop():
 242   salome.sg.updateObjBrowser(1)