Invest in Openfoam

It has been a while since my last post. This was caused by a project inside Schlumberger that took me away 3 months from Openfoam.

This was good because was a quite challenging job and made me value a lot my field job. In the meantime, I bought 2 Xeon CPU's and some other computer parts. Now ... I have a "biffy" computer to run my simulations 20 threads of processing power. Compared to my previous computer, this is a lot. ahhahaahaha


The only issue regarding these two CPUs is because I had to buy a Supermicro motherboard ... it is an EATX form factor (quite big) and here in Brazil we don't have such type of EATX case ... so ... I need to build one.

Right now my Frankenstein looks like this:

Quite ugly. But it is powerful and well cooled. hahahahahahah

I hope I find a way to build this EATX case.

Cheers!!!

Drill Pipe Rotation Influence

8 1/2" Well with 5 1/2" FH DP and 7" TJ - 60 RPM

8 1/2" Well with 5 1/2" FH DP and 7" TJ - 100 RPM

8 1/2" Well with 5 1/2" FH DP and 7" TJ - 160 RPM

8 1/2" Well with 5 1/2" FH DP and 7" TJ - No RPM

Merging Meshes from Salomé to be Used on Openfoam

This post took me a while. The main problem in here is due to mesh conversion between Salomé and Openfoam. Depending on your geometry, when you convert Salomé's unv file to Openfoam it appear some inner parts that are pretty annoying to get rid off.

If the mesh is pure unstructured you can do it easily. Just create and export. But if you have some structured mesh within a mesh, these inner parts come to play. One of my goals is to observe the velocity profile behavior and if you have an unstructured grid this won't be good, mainly because you won't have an even spaced velocity vector.

So, because of this issue, I had to redesign my straight stabilizer in order to have a hybrid grid. I had divided him into three main meshes and he is going to follow basically the same divisions stated at http://chasingaftermystuff.blogspot.com/2013/04/creating-spiral-drilling-stabilizer_12.html.

Follows the mesh from inner wall:

Zooming close to SmSpa we can see that this part is made with triangles and tetrahedrons. This is the only part on this straight stabilizer that is unstructured.

I won't publish the script capable of doing this because, by now, you should be able to do the same, but if you want the full script I can send you with no charge. : ) . Just leave a comment.

Ok.

As we saw we have mainly 3 parts (Slick, SmSpa and Stb) to complete the the whole we must mirror SmSpa and Slick. So, in the end we have 5 parts in total. It is easier to merge them with Openfoam. For each part we need to create an Openfoam case.

It is going to be like this:

As you can see, we have the all 5 parts, 3 original and 2 mirrored (MSlick and MSmSpa). Every part has an interface. The Slick part has one entrance, the Inlet, and two interfaces between Slick and SmSpa, the interface that belongs to Slick and one that belongs to SmSpa,

Showing the interface.

After exporting all your meshes to unv format and putting everything on the respective cases we must merge and stitch them.

Let us take the case where we must connect the Slick part to SmSpa. Go to the Slick folder and type:

    $ mergeMeshes ../Slick/ ../SmSpa/

    $ stitchMesh Slick_Interface_SmSpa SmSpa_Interface_Slick

On the mergeMeshes you need to input the case folders. Mines are /Slick/ and /SmSpa/. On stitchMesh is going to be the name of your interfaces. I named mine as Slick_Interface_SmSpa (belongs to Slick) SmSpa_Interface_Slick (belongs to SmSpa). You will see you had created two new folders inside of the case. These folders are going to be created according to you controlDict write Interval. If you write interval is 0.001, your first folder is going to be 0.001 and the second, 0.002. Go to the second folder and check your boundary file. If your stitchMesh went smooth you are supposed to have 6 patches but 2 of them with 0 faces. If this is the fact. Perfect. Erase these patches and correct the first number that begins the ( ). It should be 6. Now, type 4. Like this:

Copy the content of the second folder into the folder constant/polyMesh. The idea is to substitute the previous content.

By now you merged the Slick and SmSpa. Erase the folders created through stitchMesh and mergMeshes. Now we need to merge the rest. 

I've found that Openfoam has some problems to stitch more meshes ... I like to be practical, so, a workaround for this is just change to another folder. We need to merge the mesh that we just created into Stb. Go to Stb/ folder and repeat the same process.

Cheers!!!

Getting back to business!!

Hi ...

It has been a while since my last post.

As you know, Brazil has been passing trough several "small" revolutions. People's minds are changing ... politicians are afraid of that ... well ... they were.

Due to this I couldn't avoid going to the streets also. I was part ... I am part of this ... Just got tired of been robbed by these guys. Brazil is a so good country, but unfortunately we have a very bad cast of people leading us.

Between these protests I had to work in a very demanding well, with a very demanding client.

Let's see what we must do.

1. We must create a good mesh for our stabilizer. I would like to create a structured mesh because we know we can create an unstructured one.
2. Attach this mesh to a simple BHA and run OF simulations.

Cheers!

How to set Icofoam to solve the Annular Flow problem

As you observed, I solved the problem already. This post is a step-by-step to help you to solve the same problem.

What do we need?

1. Prepare a set of folders to run Openfoam;
2. Export the created mesh to unv;
3. Convert the unv to Openfoam;
4. Prepare Icofoam to run.

1 - Prepare a set of folders to run Openfoam

I'm here not to reinvent the wheel. You can see a good Openfoam introduction in http://www.openfoam.org/docs/user/cavity.php. 

At this type of annular incompressible flow I choose to use Icofoam ... it is not the closest flow model to the problem ... it is the flow model. So, pick the cavity problem itself and copy to where you are going to run your simulation.

2 - Export the created mesh to unv

The best way, to me, is to export to a unv file. 

Right click over the mesh and choose to export to UNV file

When you export the file, please, export to the root folder of your problem. For example, if you copied "cavity" folder, this mesh should be on "cavity".

3 - Convert unv file to Openfoam

If you read some of my reports, you will see, we set the environment to succeed on this process. We have the Inlet, Outlet, Inner Wall and Outer Wall boundaries. So, this process should be straight forward. If you don't succeed you are maybe facing issues in your mesh or in your Salomé (6.5!) version (Take a look at http://openfoamwiki.net/index.php/IdeasUnvToFoam).

Is your unv file in the correct folder? If yes, just type:

  ideasUnvToFoam Annular.unv

My file is called "Annular.unv"

Ok. Type:

  

  checkMesh

And in your output you should be looking for this:

Checking patch topology for multiply connected surfaces ...

    Patch                      Faces    Points   Surface topology                  

    Inlet                       1296      1368     ok (non-closed singly connected)  

    Outlet                    1296      1368     ok (non-closed singly connected)  

    Outer_Wall           14400    14472    ok (non-closed singly connected)  

    Inner_Wall            14400    14472    ok (non-closed singly connected)  

As you can see, we have everything. Inlet, Outlet, Outer_Wall and Inner_Wall.

4 - Prepare Icofoam to run

Now comes the part that took me a while to understand. Take a look on the following link http://www.openfoam.org/docs/user/boundaries.php

Let's take care of the boundary conditions. Open the "p" file on the subfolder "0" (zero) and substitute for this one:


/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.1                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    Inlet
    {
        type            zeroGradient;
    }

    Outlet
    {
        type            fixedValue;
        value           uniform 0;
    }

    Inner_Wall
    {
        type            zeroGradient;
    }

    Outer_Wall
    {
        type            zeroGradient;
    }

}

// ************************************************************************* //

You can find an explanation of this file in here: http://www.openfoam.org/docs/user/basic-file-format.php.

But the there are some few things that I must call your attention.

On the line that says:

   dimensions      [0 2 -2 0 0 0 0];

Is saying that you have a dimension of m²/s². Well, it is a pressure file so, we should be expecting in SI units a dimension of kg/(m.s²). What is the problem? No problem at all. With Icofoam we will use kinematic viscosity and to keep the units coherent on the solver we must divide the pressure by density.

Type zeroGradient:

As you can see, I'm using zeroGradient in every patch except "Outlet". It is only because I'm not using the gravity for this problem and solving the Navier-Stokes Equations for this setting you will find a zero pressure gradient. On the Outlet I'm using fixed value because the flow runs is not constrained anymore and goes to the environment that has zero pressure.

On the same "0" subfolder let us see the "U" file and substitute for the following:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.1                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{
    Inlet
    {
        type            fixedValue;
        value           uniform (0 0 1.66);
    }

    Outlet
    {
        type            zeroGradient;
    }
 
    Inner_Wall
    {
type            fixedValue;
        value           uniform (0 0 0);
    }

    Outer_Wall
    {
type            fixedValue;
        value           uniform (0 0 0);
    }
}

// ************************************************************************* //


As you can see, is very straightforward. You have the no-slip condition on the walls and the velocity inlet.

Now we should move to the subfolder "constant". In there you should find the polyMesh folder (don't need to touch it. It came from the unv conversion) and "transportProperties" file. Open it and substitute for:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.1                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "constant";
    object      transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu              nu [ 0 2 -1 0 0 0 0 ] 7.27e-5;


// ************************************************************************* //

As I said, IcoFoam uses the kinematic viscosity. You have the correct dimensions and the attribute value that came from http://chasingaftermystuff.blogspot.com.br/2013/05/using-icofoam-to-solve-annular-flow.html.

To the folder "System" I would advise you to wait a little bit more. I don't understand it as a whole yet and would suggest you to look into Openfoam's website for further information.

But if you have any doubt, just contact me ...

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)
##############

Back to Business

Ufa ... this last embark was quite something. We had two influxes: one from water and another from gas. It has been some time that I didn't participate of a kick. But, thanks for the people we had on the rig (Toolpusher, CoMan, Drillers and us, DD's) we had identified everything from the beginning. A real team work. No harms were done to the people nor to the environment.

There are some few things I observed that I would like to test using CFD:

• How the pump pressure is going to behave when we have a lighter mud on the annular?
• How the pump pressure is going to behave when we have a heavier mud on the annular?
• How the gas is going to behave inside of the annular while migrating?
• How mud weight acts on gas migration? 

I didn't have time to search on SPE papers these items, and they do exist in there ... but ... I want to try to solve them using CFD.

While on the rig I couldn't run any simulation but I had studied Salomé a little bit more. Now I'm able to create a quite simple annular structured mesh with viscous layers. 

So, to a next post I'll use this method and solve again the annular problem using Openfoam. Expect to find on the next two posts:

• Creating a Structured Annular Grid with Viscous Layers in Salomé;
• How to set Icofoam to solve the Annular Flow problem.

Visit from the sky.

Dolphins trying to hijack the rig. We had to use water hoses to disperse them. :P

Comparing Velocity Plot from Simulation's Result Against Analytical Equation

We have several ways of comparing the data. So ... my comparison is only based on mean values. If you want something deep you should look for some statistical book regarding data analysis.

Taking the last graphic and rearranging the data to be fitted on the radius together with the analytical solution we can see the following result:

Graph generated using QtiPlot

Using the "eye" we can see that the simulation was very good, but taking into account velocity's mean value we have a small difference of 0.45% between analytical solution and the simulated one.

Remember this is a weak data analysis but we are paving the road to a broader area.

Right know I must prepare myself to embark, so, the next post regarding of how to get the simulation done and the graphics is going to take a little longer.

Usually, my embark rotation takes 28 days and will try to post something from the rig, but don't expect too much ... maybe some sunset pictures.

Cheers and see you in some few days!!!!

Using Icofoam to Solve Annular Flow

It took me a while to solve this problem. More than expected in fact.

I was able to run the problem since the beginning but the values were too away from what I expected and I thought the problem was related to something behind Icofoam ... not the code, far from that ...  but, my problem setting.

The values were wrong due to two things a coarse mesh and , I believe, floating number precision. I had to produce a finer mesh to see the actual results and increase the tool's length. I came with the following conclusion: we need to evaluate how the tool's length is going to affect the results previously. Took me a while to simulate everything.

This steady, laminar flow with Newtonian fluid has an analytical solution and can be easily derived from Navier-Stokes equations. Obs: For this time I will not consider the gravity.

Velocity profile is given by the following equation:

And the pressure drop on an annular section of length l would be:

Where:

p  → pressure

ρ  → density

ν  → kinematic viscosity

r_o → outer wall radius 

r_i  → inner wall radius

v_m → average annular velocity  

All units are in SI.

From the dimensions in http://chasingaftermystuff.blogspot.com.br/2013/04/creating-mesh.html we will change the tool's length to 10m , ν = 7.27e-5 m²/s, ρ = 1.10e+3 kg/m³ and v_m = 1.66 m/s we will have a Reynolds number of 1159 and we can safely assume that the flow is laminar. So, Δp/ρ is assumed to be 22.43 m²/s²,  Δp = 24.67 kPa.

We have the pressure drop and the velocity profile. Taking the results from the simulation we can see the values agree with minor error (we still have space to increase the mesh definition). You can see that in the graphics the "no-slip" condition on the walls are obeyed.

Pressure divided by density

Velocity magnitude contour

Velocity vector over velocity contour

Velocity plot against cells. 

For a next post I will compare simulated velocity plots against its analytical solution and solve another exercise.

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