Minor Setback - Creating a Circular Annular Structured Grid With Salomé

I had a minor set back.

To write a report I usually solve the problem first and repeat the "solving" while making the report. Unfortunately, I had everything done (grid, chosen Openfoam solver, parameters, etc) but the problem was not getting stable. Courant numbers were reaching high values on the third iteration.

When you face something like that you have to attack the problem in two ways: Dt and Ds. In other words, pacing time and the grid. For this case, the problem was related to time. But when you reduce the pace so much, you'll have to wait a lot to see your results. I don't have a shitty computer and of course is not an academic cluster. It has 4 cores and I can divide into 8 threads, but even so, depending on the problem, can take from 3 minutes to 18 hours to solve it.

My Openfoam solver is not calibrated yet and it is quite boring to tweak every little detail.

I like Salomé ... a lot ... and to test the solver you must have some few grids to perform some sensitive test while simulating, because of this I decided to make a structured grid with this annular.

It is really not that complicated when you discover the way but it has some "magic" behind it.

Lets begin!!

Follow the base geometry script, I changed a little the previous one, but it has the same core:

######


import sys
import salome
import geompy
import math

salome.salome_init()
theStudy = salome.myStudy

import salome_notebook
notebook = salome_notebook.notebook

gg = salome.ImportComponentGUI("GEOM")

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

# BHA components parameters

Wellbore = 8.5 # in

OD_Component = 6.5 # in

Length_Component =2 # m

# Converting imperial units to SI

Wellbore = Wellbore*0.0254

OD_Component = OD_Component*0.0254

# Base vector for rotation and extrusion

Origin = geompy.MakeVertex (0,0,0)

EndVector = geompy.MakeVertex (0,0,1)

Vector = geompy.MakeVector(Origin ,EndVector)

# Base line for wellbore annular

pOD = geompy.MakeVertex(OD_Component,0,0)

pWell = geompy.MakeVertex(Wellbore,0,0)

vAnnular = geompy.MakeVector(pOD,pWell)

# Base annular disc

wAnnular = geompy.MakeWire([vAnnular])

Extrusion_1 = geompy.MakePrismVecH(wAnnular, OZ, 1)

[Edge_1,Seed,Edge_3,Edge_4] = geompy.ExtractShapes(Extrusion_1, geompy.ShapeType["EDGE"], True)

Annular_Simple= geompy.MakeRevolution(Extrusion_1, OZ, 360*math.pi/180.0)

Annular = geompy.MakeCompound([Annular_Simple, Extrusion_1])

[Inlet,Inner_Wall,Outer_Wall,Outlet,Face_Seed] = geompy.ExtractShapes(Annular, geompy.ShapeType["FACE"], True)

[Edge_2,Edge_Seed_Base,Edge_Seed_Top,Edge_7] = geompy.ExtractShapes(Face_Seed, geompy.ShapeType["EDGE"], True)

geompy.addToStudy( Annular, 'Annular' )

geompy.addToStudy( Extrusion_1, 'Extrusion_1' )

geompy.addToStudyInFather( Extrusion_1, Edge_1, 'Edge_1' )

geompy.addToStudyInFather( Extrusion_1, Seed, 'Seed' )

geompy.addToStudyInFather( Extrusion_1, Edge_3, 'Edge_3' )

geompy.addToStudyInFather( Extrusion_1, Edge_4, 'Edge_4' )

geompy.addToStudyInFather( Annular, Inlet, 'Inlet' )

geompy.addToStudyInFather( Annular, Inner_Wall, 'Inner_Wall' )

geompy.addToStudyInFather( Annular, Outer_Wall, 'Outer_Wall' )

geompy.addToStudyInFather( Annular, Outlet, 'Outlet' )

geompy.addToStudyInFather( Annular, Face_Seed, 'Face_Seed' )

geompy.addToStudyInFather( Annular, Edge_Seed_Base, 'Edge_Seed_Base' )

geompy.addToStudyInFather( Annular, Edge_Seed_Top, 'Edge_Seed_Top' )

######

With the created geometry we can have our structured grid.

Activate the Salomé's Mesh Module.

Go to Mesh >> Create Mesh. Pay attention on this one. You'll choose <None> on everything in here. <None>!!!

Go to Mesh >> Create Sub-Mesh. The chosen geometry is a face. I call it Face_Seed because this guy helped us making the solid and we are going to use it to produce the grid as well. The algorithm is Quadrangle (Mapping) because is going to divide this face into structured small pieces. Regarding the Hypothesis ... just standard.

In the same box, choose "1D" and choose "Wire Discretisation". We will divide the z-axial length into equal pieces.

Click on Hypothesis and follow these parameters (or choose another one if you prefer). Apply. (Don't Apply and Close ... just Apply).

Select "Edge_Seed_Base" and choose "Wire Discretisation". Instead of 200 divisions I will choose 20 for this one on the Hypothesis. Apply.

Repeat the same procedure to "Edge_Seed_Top" but press Apply and Close.

Right Click over the Mesh_1 and Compute the mesh.

You can see you have some few warnings on your left just below the Mesh_1. Don't worry. Now, go to Modification >> Revolution. The angle 360 is because we will rotate 360 degrees and "Number of Steps" is 36 because it is going to divide 360 by 36 ... So we will divide the circle into 36 pieces and every piece corresponds to 10 degrees.

You should have something like this.

We have the structured grid already, but we can't export like this. We must define the boundaries and that is when comes the tricky part.

Go to Mesh >> Create Group.

Repeat the same procedure to "Outlet" and Apply. Remember, when setting the filter Threshold Value (third column) choose "Outlet" on the geometry.

Now we can proceed to "Inner_Wall". Apply.

Proceed the same way with "Outer_Wall". The same thing here ... when setting the filter Threshold Value choose "Outer_Wall" on the geometry. Apply and Close.

Now we have the mesh and the boundaries. Everything Ok to export to UNV format and then Openfoam?! NO!!!

If you export like this you won't be able to run your simulations. Unfortunately when you set the groups while using the filter you picked some nodes that belonged to another geometry. For example, when you created the Inlet you are picking up nodes that belong to Inner_Wall and Outer_Wall. We must erase the excess.

Go to Modification >> Transformation >> Merge Nodes.

Now we are done. We can export this grid!

Creating the mesh

It is well known that we have several types of grid generation methods and I want to be practical. Which means that I'll try to have the best method without taking care of the minor details. To have all my simulations running, I must choose from the beginning.

Of course, through the whole learning curve I should be adapting which kind of method best fits for the corresponding problem.

From now I have the following assumptions:

• Flow is going to happen within an annular space created by a circular well-bore and a circular tool ;
• We won't have any type of heat exchange;
• We are going to have a non-Newtonian fluid within annular space;
• Rheology values are going to be kept constant;
• Well-bore and the tool will be centralized;
• It is an incompressible flow;
• The grid will be a mix of structured and unstructured ones, but I will be trying to use structured grids as much as possible;

Creating the Grid for Simple Axial Annular Flow with Salomé

In order to have my simulations running I must create a grid and this geometry is quite simple.

• Length = 2 m
• Well Bore OD = 8.5 in
• Pipe OD = 6.5 in
• g = -9.81 m/s2

You must be captured we'll have our flow in this space created between 8.5 in and 6.5 in, 1 inch of annular space. Also, looking at this simple drawing is easy to see we must have some structured grid meshing, right? No.

Structured grids can be very time consuming. My focus right now is to have some few simulations running. You will see the base of this grid is unstructured but which makes it unstructured as a whole. You can use as reference the following link http://www.salome-platform.org/user-section/salome-tutorials/edf-exercise-4 to see how to create some structured grid with Salomé. But be aware, Salomé is not that good with structured grids either.

To create the geometry I'll use the following Salomé script:
###########


import sys
import salome
import geompy
import math

salome.salome_init()
theStudy = salome.myStudy

import salome_notebook
notebook = salome_notebook.notebook

gg = salome.ImportComponentGUI("GEOM")

# BHA components parameters
Wellbore = 8.5 # in
OD_Component = 6.5 # in
ID_Component = 1 # in
Length_Component = 2 # m


# Converting imperial units to SI
Wellbore = Wellbore*0.0254
OD_Component = OD_Component*0.0254
ID_Component = ID_Component*0.0254


# Base vector for rotation and extrusion
Origin = geompy.MakeVertex (0,0,0)
EndVector = geompy.MakeVertex (0,0,1)
Vector = geompy.MakeVector(Origin ,EndVector)


# Base line for wellbore annular
pOD = geompy.MakeVertex(OD_Component,0,0)
pWell = geompy.MakeVertex(Wellbore,0,0)
vAnnular = geompy.MakeVector(pOD,pWell)


# Base annular disc
wAnnular = geompy.MakeWire([vAnnular])
revolution = geompy.MakeRevolution(wAnnular, Vector, 2*math.pi)
Annular = geompy.MakePrismVecH(revolution, Vector, Length_Component)
id_Annular = geompy.addToStudy(Annular , "Annular")

[Inlet,Inner_Wall,Outer_Wall,Outlet] = geompy.ExtractShapes(Annular, geompy.ShapeType["FACE"], True)
geompy.addToStudyInFather( Annular, Inlet, 'Inlet' )
geompy.addToStudyInFather( Annular, Inner_Wall, 'Inner_Wall' )
geompy.addToStudyInFather( Annular, Outer_Wall, 'Outer_Wall' )
geompy.addToStudyInFather( Annular, Outlet, 'Outlet' )

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

###################

I will use the results from this script to work on my mesh.

Toggle the Mesh Environment or choose Mesh on the list.

Go to "Mesh" and choose "Create Mesh". You need to choose the geometry we will be working on.

Mesh >> Create Mesh

Now we must choose what method of grid generation we have to use. I had worked with several meshes regarding this annular problem and I had identified some few things. I'm choosing the extrusion of the face to generate the rest of the mesh, because it satisfies the resolution of our problem and is also faster. If we had generated the whole grid just using, lets say, NETGEN 3D-2D-1D, the grid would have a very dense resolution without being practical.

Choose 3D extrusion:

Don't "Apply and Close" we have some few things to do with this subject. While extruding in z axis we must define how the meshing is going to behave. I want quadrangle:

Let it in standard

We must define how many "cuts" we do want in z axis. I choose 200:

You can put more or less than 200

Now we can "Apply and Close"

We need to define this extrusion's base:

I want 0.01 of Max Element Area

Still at the extrusion base we will divide the circumference into 90 segments:

Choose 90 segments

Do the same for the "Outlet"!!

Now we must define in the grid what is Inlet, Outlet, Inner_Wall and Outer_Wall.

Without this definition we can't run Openfoam.

Repeat this procedure to Outlet, Inner_Wall and Outer Wall.

Now we can compute the mesh. Right click over Mesh_1 and choose Compute.

Setting the next moves

In order to keep advancing into the annular flow we need to identify the next immediate moves. I've been learning Openfoam (in my free time) for the last 5 months. I can quite say, Openfoam's learning curve is quite steep and that is why we must divide the task into small ones.

The next immediate moves:

1. Create a simple annular flow problem with its respective flow variables.
2. Create the mesh.
3. Pick one model from Openfoam's tutorial package and adapt to the problem.
4. Compare the result with analytical solution.

I believe I can finish steps one to three in two days. Step 4 is going to take a little bit longer because I need to solve Navier-Stokes equation regarding the simple problem and after that I must generate the graphics to compare the results

Being towed in Atlantic Ocean

Exercise - Moineau Rotor - Solution

As promised, we can observe the solution from the previous exercise.

I tried to record a video to show the solution but my recording gear is terrible. I'll search for a good microphone.

From the exercise we can extract all the important variables.

Follow:

• Rotor Lobes : 7
• Rotor Length : 5.8 m
• MAX OD : 6.5 in
• MIN OD : 5 in
• Pitch : 1.45 m/rev

One thing we might have to infer while drawing the rotor is the lobe offset. I will show you while performing the drawing.

To make our lives easier we will use all the above values at the end. I'll show you how to fix everything when we finished with the drawing.

• Rotor Lobes : 7
• Rotor Length : 2*
• MAX OD : 14 m
• MIN OD : 12 m
• Pitch : 2* m/rev

First, we know that hypocycloid is a set of parametric equations and can be taken at http://en.wikipedia.org/wiki/Hypocycloid. 

On the following picture you will observe I put the set of parametric equations and the rest of the values. We don't want a whole hypocycloid. We want only a piece of it and from this piece we can provide the rest. For "Min t" I'm using 0.15 /rad because if I put zero I'm going to have a spike (and we don't need that).  "Max t" is going to be  divided by the quantity of lobes. In this case, ~0.448799.

New Entity >> Basic >> Curve

Right know comes the "engineering" assumption. A Moineau pump, of course, doesn't have all those cuspids and we need to round them up. To create the lobes we must have a rolling small inner circle inside of a big one. For this case the big circle has a 7 m radius and to have a seven lobes we must have the a small circle with 1 m  radius. So ... we need to create an "insider" tangent circle.

Create the point:

New Entity >> Basic >> Point

Create the circle based on the last point:

New Entity >> Basic >> Circle

Build a wire:

New Entity >> Build >> Wire

Build an edge (you do want to do this because we can extract some values from this edged curve):

New Entity >> Build >> Edge

By now, we should observe the hypocycloid curve is intersecting the small circle. That is when we have the offset. We want this curve to be  as much tangent as possible. I'm using an offset of 1.17 ... it is going to give you a nice tangent figure.

Operations >> Transformation >> Scale Transform

Fuse the results from the scaled up hypocycloid curve and the small circle.

Operations >> Boolean >> Fuse

You must be asking yourself "What are you doing?!!". But this is the idea, to make your life easy. If you have symmetry in your subject you must use it in your favor, otherwise, you would be wasting time. 

We must clean the resultant fuse. Explode the edges and create a curve for only what we need to have. Remember, we are aiming on a lobe design. 

New Entity >> Explode

Choose which edges we need:

In my case I'm going to need edges 3,5 and 6.

Create a wire from the chosen edges:

New Entity >> Build >> Wire

Mirror the resultant curve and fuse both curves. We can see something close to a lobe. 

Operations >> Transformation >> Mirror Image

Explode it into vertexes.

New Entity >> Explode

In my case, I will need the beginning and the of the curves. Vertex 2 and Vertex 3. By having these two vertexes we can create a reverse arc.

New Entity >> Basic >> Arc

We can see a wire from a lobe. We are almost there. Just create a a wire with this arc and the connected curve.

New Entity >> Build >> Wire

Build a face:

New Entity >> Build >> Face

Do you remember that I stated in a Moineau pump the axial section would look like a helix? We need to make two helical curves. One curve is going to give me the helical path and the other curve is going to give the binormal vector. Without this binormal vector the lobe would get flat while rounding up. We don't want that!!! 

First curve must have a radius of 6 m (small circle's center) and the second, 7 m (hypocycloid radius). We want only a full circle that is why we must have 2 at "Max t".

Build the first curve and I will let you build the second:

New Entity >> Basic >> Curve

Create an extrusion along a path. Remember we will use a binormal!!

New Entity >> Generation >> Extrusion Along a Path

We have one lobe! Time to make 6 more.

Perform a multi-rotation:

Operations >> Transformation >> Multi-Rotation

We have a nice screw. Time to put everything in corrected values. To do this we need to scale down every axis by the corrected amount. It takes into account unit and size conversion only.

Operations >> Transformation >> Scale Transform

1.45 m Moineau Rotor , if you want the complete 5.8 m just translate a copy and fuse them.

Answer

Any doubt, just ask.

Cheers!!!!