tests/enforced_mesh.py

   1 # -*- coding: utf-8 -*-
   2 # Copyright (C) 2017-2022  CEA/DEN, EDF R&D
   3 #
   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.
   8 #
   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.
  13 #
  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
  17 #
  18 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  19 #
  20
  21 import sys
  22 import salome
  23
  24 salome.salome_init()
  25
  26 ###
  27 ### GEOM component
  28 ###
  29
  30 import GEOM
  31 from salome.geom import geomBuilder
  32 import math
  33 import SALOMEDS
  34
  35
  36 geompy = geomBuilder.New()
  37
  38 # first cylinder
  39 r1 = 0.5
  40 h1 = 5
  41
  42 # second cylinder
  43 r2 = 0.3
  44 h2 = 3
  45
  46 length_piquage = 1.5
  47
  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' )
  56
  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)
  61
  62 piquage = geompy.MakeFuseList([Cylinder_1, Translation_1], True, True)
  63 geompy.addToStudy( piquage, 'piquage' )
  64
  65 Inlet_z = geompy.GetFaceNearPoint(piquage, O)
  66 geompy.addToStudyInFather( piquage, Inlet_z, 'Inlet_z' )
  67
  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' )
  71
  72 p_outlet = geompy.MakeVertex(0, 0, h1)
  73 Outlet = geompy.GetFaceNearPoint(piquage, p_outlet)
  74 geompy.addToStudyInFather( piquage, Outlet, 'Outlet' )
  75
  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' )
  81
  82 p_corner = geompy.MakeVertex(-r2, 0, length_piquage+r2)
  83 corner = geompy.GetVertexNearPoint(piquage, p_corner)
  84 geompy.addToStudyInFather( piquage, corner, 'corner' )
  85
  86 p1 = geompy.MakeVertex(0, -0.25, 1)
  87 p2 = geompy.MakeVertex(0, 0.25, 1)
  88 p3 = geompy.MakeVertex(0, 0.25, 3)
  89 p4 = geompy.MakeVertex(0, -0.25, 3)
  90
  91 wire = geompy.MakePolyline([p1, p2, p3, p4, p1])
  92 face = geompy.MakeFace(wire, 1, theName="face")
  93
  94 ###
  95 ### SMESH component
  96 ###
  97
  98 import  SMESH, SALOMEDS
  99 from salome.smesh import smeshBuilder
 100
 101 from salome.StdMeshers import StdMeshersBuilder
 102
 103 smesh = smeshBuilder.New()
 104
 105 Mesh_faces = smesh.Mesh(face)
 106 Mesh_faces.Triangle(algo=smeshBuilder.MG_CADSurf)
 107 Mesh_faces.Compute()
 108
 109
 110 # Viscous layers with Netgen additional hypothesis
 111 # ================================================
 112
 113 Mesh_1 = smesh.Mesh(piquage)
 114 NETGEN_2D = Mesh_1.Triangle(algo=smeshBuilder.NETGEN_1D2D)
 115 NETGEN_2D_Parameters = NETGEN_2D.Parameters()
 116
 117 NETGEN_2D_Parameters.SetMinSize( 0.01 )
 118 NETGEN_2D_Parameters.SetLocalSizeOnShape(corner, 0.01)
 119 NETGEN_2D_Parameters.SetFineness( 5 )
 120 NETGEN_2D_Parameters.SetGrowthRate( 0.1 )
 121 NETGEN_2D_Parameters.SetNbSegPerEdge( 2 )
 122 NETGEN_2D_Parameters.SetNbSegPerRadius( 3 )
 123
 124 wall_faces = geompy.SubShapeAll(Wall, geompy.ShapeType["FACE"])
 125 wall_ids = geompy.GetSubShapesIDs(piquage, wall_faces)
 126
 127 NETGEN_3D = Mesh_1.Tetrahedron()
 128 Viscous_Layers_1 = NETGEN_3D.ViscousLayers(0.05,3,1.1,[],1,StdMeshersBuilder.SURF_OFFSET_SMOOTH)
 129 Viscous_Layers_1.SetTotalThickness( 0.05 )
 130 Viscous_Layers_1.SetNumberLayers( 3 )
 131 Viscous_Layers_1.SetStretchFactor( 1.1 )
 132 Viscous_Layers_1.SetMethod( StdMeshersBuilder.SURF_OFFSET_SMOOTH )
 133 Viscous_Layers_1.SetFaces( wall_ids, 0 )
 134
 135
 136 #isDone = Mesh_1.Compute()
 137 #Mesh_1.SplitVolumesIntoTetra( Mesh_1, 1 )
 138
 139 Outlet_1 = Mesh_1.GroupOnGeom(Outlet,'Outlet',SMESH.FACE)
 140 Wall_1 = Mesh_1.GroupOnGeom(Wall,'Wall',SMESH.FACE)
 141 Inlet_z_1 = Mesh_1.GroupOnGeom(Inlet_z,'Inlet_z',SMESH.FACE)
 142 Inlet_x_1 = Mesh_1.GroupOnGeom(Inlet_x,'Inlet_x',SMESH.FACE)
 143
 144 # Enforced mesh without layer
 145 # ===========================
 146
 147 Mesh_2 = smesh.Mesh(piquage, "Mesh_without_layer")
 148
 149 NETGEN_2D_1_1 = Mesh_2.Triangle(algo=smeshBuilder.NETGEN_1D2D)
 150 status = Mesh_2.AddHypothesis(NETGEN_2D_Parameters)
 151
 152 MG_Hybrid_2 = Mesh_2.Tetrahedron(algo=smeshBuilder.HYBRID)
 153 MG_Hybrid_Parameters_2 = MG_Hybrid_2.Parameters()
 154 MG_Hybrid_Parameters_2.SetLayersOnAllWrap( 0 )
 155 MG_Hybrid_Parameters_2.SetElementGeneration( 0 )
 156 #MG_Hybrid_Parameters_2.SetHeightFirstLayer( 0.01 )
 157 #MG_Hybrid_Parameters_2.SetBoundaryLayersProgression( 1.1 )
 158 #MG_Hybrid_Parameters_2.SetNbOfBoundaryLayers( 3 )
 159 MG_Hybrid_Parameters_2.SetEnforcedMeshWithGroup( Mesh_faces.GetMesh(), SMESH.FACE, "LayersGroup" )
 160 MG_Hybrid_Parameters_2.SetKeepFiles(1)
 161
 162 isDone = Mesh_2.Compute()
 163
 164 if not isDone:
 165   raise Exception("Error when computing Mesh_without_layer")
 166
 167 # Check that a group has been created with the enforced mesh
 168 assert len(Mesh_2.GetGroups()) == 1
 169 assert Mesh_2.GetGroups()[0].GetName() == 'LayersGroup'
 170 assert Mesh_2.GetGroups()[0].Size() > 0
 171
 172 # Enforced mesh with layers
 173 # =========================
 174
 175 Mesh_3 = smesh.Mesh(piquage, "Mesh_with_layer")
 176
 177 NETGEN_2D_1_1 = Mesh_3.Triangle(algo=smeshBuilder.NETGEN_1D2D)
 178 status = Mesh_3.AddHypothesis(NETGEN_2D_Parameters)
 179
 180 MG_Hybrid_3 = Mesh_3.Tetrahedron(algo=smeshBuilder.HYBRID)
 181 MG_Hybrid_Parameters_3 = MG_Hybrid_3.Parameters()
 182 MG_Hybrid_Parameters_3.SetLayersOnAllWrap( 0 )
 183 MG_Hybrid_Parameters_3.SetElementGeneration( 0 )
 184 MG_Hybrid_Parameters_3.SetHeightFirstLayer( 0.01 )
 185 MG_Hybrid_Parameters_3.SetBoundaryLayersProgression( 1.1 )
 186 MG_Hybrid_Parameters_3.SetNbOfBoundaryLayers( 3 )
 187 MG_Hybrid_Parameters_3.SetEnforcedMeshWithGroup( Mesh_faces.GetMesh(), SMESH.FACE, "LayersGroup" )
 188 MG_Hybrid_Parameters_3.SetKeepFiles(1)
 189 MG_Hybrid_3.SetFacesWithLayers( Wall )
 190
 191
 192 isDone = Mesh_3.Compute()
 193
 194 if not isDone:
 195   raise Exception("Error when computing Mesh_with_layer")
 196
 197 # Check that a group has been created with the enforced mesh
 198 assert len(Mesh_3.GetGroups()) == 1
 199 assert Mesh_3.GetGroups()[0].GetName() == 'LayersGroup'
 200 assert Mesh_3.GetGroups()[0].Size() > 0
 201