quinta-feira, 6 de junho de 2013

Creating a Structured Annular Grid with Viscous Layers in Salomé

Hi, from the previous posts you should be able to read the following python code to create your structured annular grid.

What is different from previous posts? Right now the code is very simple and we have viscous layers on the walls.

I had learned how to do this by taking a look at Salomé's forums. But, the function that called my attention and made this possible was PROPAGATE (http://docs.salome-platform.org/salome_7_2_0/gui/GEOM/geompy_doc/group__l3__blocks__op.html#ga81185dc6c66632dff79c2f0a19f77537)

With this function I could divide my inlet/outlet in a way that I could use the quadrangle mesh.

Again, if you have any doubt, just send me an email.


Structured Annular Mesh with Viscous Layer


###########
import sys
import salome

salome.salome_init()
theStudy = salome.myStudy

import salome_notebook
notebook = salome_notebook.notebook

###
### GEOM component
###

import GEOM
import geompy
import math
import SALOMEDS

gg = salome.ImportComponentGUI("GEOM")
geompy.init_geom(theStudy)

O = geompy.MakeVertex(0, 0, 0)
OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)

Wellbore = 8.5 # in
OD_Component = 6.5 # in
Length_Component = 10 # m

Wellbore = geompy.MakeCylinder(O, OZ, 0.5*Wellbore, Length_Component)
OD_Component = geompy.MakeCylinder(O, OZ, 0.5*OD_Component, Length_Component)
Annular_1 = geompy.MakeCut(Wellbore, OD_Component)
Annular_2 = geompy.MakeScaleAlongAxes(Annular_1, O, 0.0254, 0.0254, 1) # converting OD to imperial units

# Preparing for mesh
Plane_1 = geompy.MakePlaneLCS(None, 2000, 3)
Annular = geompy.MakeHalfPartition(Annular_2, Plane_1)
[Compound_1, Compound_2, Compound_3, Compound_4] = geompy.Propagate(Annular)
[Face_1,Face_2,Face_3,Face_4,Face_5,Face_6,Face_7,Face_8,Face_9,Face_10] = geompy.ExtractShapes(Annular, geompy.ShapeType["FACE"], True)
listSubShapeIDs = geompy.SubShapeAllIDs(Annular, geompy.ShapeType["FACE"])
Inlet = geompy.CreateGroup(Annular, geompy.ShapeType["FACE"])
geompy.UnionIDs(Inlet, [42, 26])
Outlet = geompy.CreateGroup(Annular, geompy.ShapeType["FACE"])
geompy.UnionIDs(Outlet, [45, 14])
Outer_Wall = geompy.CreateGroup(Annular, geompy.ShapeType["FACE"])
geompy.UnionIDs(Outer_Wall, [4, 38])
Inner_Wall = geompy.CreateGroup(Annular, geompy.ShapeType["FACE"])
geompy.UnionIDs(Inner_Wall, [48, 34])
Inlet_Outlet_Wire = geompy.CreateGroup(Annular, geompy.ShapeType["COMPOUND"])
geompy.UnionList(Inlet_Outlet_Wire, [Compound_3, Compound_4])

geompy.addToStudy( Wellbore, 'Wellbore' )
geompy.addToStudy( OD_Component, 'OD_Component' )
geompy.addToStudy( Annular, 'Annular' )
geompy.addToStudyInFather( Annular, Compound_1, 'Compound_1' )
geompy.addToStudyInFather( Annular, Compound_2, 'Compound_2' )
geompy.addToStudyInFather( Annular, Compound_3, 'Compound_3' )
geompy.addToStudyInFather( Annular, Compound_4, 'Compound_4' )
geompy.addToStudyInFather( Annular, Inlet, 'Inlet' )
geompy.addToStudyInFather( Annular, Outlet, 'Outlet' )
geompy.addToStudyInFather( Annular, Outer_Wall, 'Outer_Wall' )
geompy.addToStudyInFather( Annular, Inner_Wall, 'Inner_Wall' )


###
### SMESH component
###

import smesh, SMESH, SALOMEDS

smesh.SetCurrentStudy(theStudy)
import StdMeshers

Mesh_1 = smesh.Mesh(Annular)

Z_axis_segments = 200
Radius_segments = 10
Radius_divisions = 36 # For each half of inlet/outlet - It is going to divide 180 degrees for Radius_divisions

# Number of segments in z axis
Regular_1D = Mesh_1.Segment()
Nb_Segments_1 = Regular_1D.NumberOfSegments(Z_axis_segments)
Nb_Segments_1.SetDistrType( 0 )

# Number of radius segments
Regular_1D_1 = Mesh_1.Segment(geom=Inlet)
Nb_Segments_2 = Regular_1D_1.NumberOfSegments(Radius_segments)
Nb_Segments_2.SetDistrType( 0 )
Regular_1D_2 = Mesh_1.Segment(geom=Outlet)
status = Mesh_1.AddHypothesis(Nb_Segments_2,Outlet)

# Number of radius divisions
Regular_1D_3 = Mesh_1.Segment(geom=Inlet_Outlet_Wire)
Nb_Segments_3 = Regular_1D_3.NumberOfSegments(Radius_divisions)
Nb_Segments_3.SetDistrType( 0 )

# Type of meshing - Quadrangle for this case
Quadrangle_2D = Mesh_1.Quadrangle(algo=smesh.QUADRANGLE)
Quadrangle_Parameters_1 = Quadrangle_2D.QuadrangleParameters(StdMeshers.QUAD_STANDARD)

Quadrangle_2D_1 = Mesh_1.Quadrangle(algo=smesh.QUADRANGLE,geom=Inlet)
status = Mesh_1.AddHypothesis(Quadrangle_Parameters_1,Inlet)

Quadrangle_2D_2 = Mesh_1.Quadrangle(algo=smesh.QUADRANGLE,geom=Outlet)
status = Mesh_1.AddHypothesis(Quadrangle_Parameters_1,Outlet)

# Creation of a hexahedron meshing - it is going to extrude the base along z axis
Hexa_3D = Mesh_1.Hexahedron(algo=smesh.Hexa)
Viscous_Layers_1 = Hexa_3D.ViscousLayers(0.001,4,1,[ 42, 26, 45, 14 ]) # Viscous layer

isDone = Mesh_1.Compute()
Inlet_1 = Mesh_1.GroupOnGeom(Inlet,'Inlet',SMESH.FACE)
Outlet_1 = Mesh_1.GroupOnGeom(Outlet,'Outlet',SMESH.FACE)
Outer_Wall_1 = Mesh_1.GroupOnGeom(Outer_Wall,'Outer_Wall',SMESH.FACE)
Inner_Wall_1 = Mesh_1.GroupOnGeom(Inner_Wall,'Inner_Wall',SMESH.FACE)

## set object names
smesh.SetName(Mesh_1.GetMesh(), 'Mesh_1')
smesh.SetName(Regular_1D.GetAlgorithm(), 'Regular_1D')
smesh.SetName(Nb_Segments_1, 'Nb. Segments_1')
smesh.SetName(Quadrangle_2D.GetAlgorithm(), 'Quadrangle_2D')
smesh.SetName(Quadrangle_Parameters_1, 'Quadrangle Parameters_1')
smesh.SetName(Hexa_3D.GetAlgorithm(), 'Hexa_3D')
smesh.SetName(Viscous_Layers_1, 'Viscous Layers_1')
smesh.SetName(Nb_Segments_2, 'Nb. Segments_2')
smesh.SetName(Nb_Segments_3, 'Nb. Segments_3')
smesh.SetName(Inlet_1, 'Inlet')
smesh.SetName(Outlet_1, 'Outlet')
smesh.SetName(Outer_Wall_1, 'Outer_Wall')
smesh.SetName(Inner_Wall_1, 'Inner_Wall')

if salome.sg.hasDesktop():
  salome.sg.updateObjBrowser(1)
##############

Nenhum comentário:

Postar um comentário