CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Programming & Development

Convert coordinate system in OpenFOAM

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree5Likes

Reply
 
LinkBack Thread Tools Display Modes
Old   June 2, 2010, 05:03
Default Convert coordinate system in OpenFOAM
  #1
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Dear all,

I am working with a project where we would like to implement a surface capturing method that requires information of the coordinate system in the form x(i,j,k). I have understood that OpenFOAM uses a cell list of the form x(i), so you do not really know the ijk position of a cell. The grid we want to use is a simple orthogonal grid. I would like to know the following:
  1. Does anyone know of a method to retrieve the "ijk information" of each cell in OpenFOAM?
  2. How can the cell number (or cell position info) be accessed in OpenFOAM? Is it a "real cell number" or a e.g. a pointer to a location in memory?
Many thanks in advance!
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 2, 2010, 05:31
Default
  #2
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
Hi

Well, you can achieve the information from the following building blocks:

Code:
// Neighbour and owner cells to given face
const unallocLabelList & nei = mesh.neighbour();
const unallocLabelList & own = mesh.owner();

forAll( mesh.C(), celli)
{
    // This is a list of faces of which the cell is build
    const cell& cProp = mesh_.cells()[celli];

    // Loop over all faces
    forAll(cProp, facei)
    {
         // Using own and nei you can determine  x(i,j,k)
    }

}
What method is it you are considering to implement?

Best regards,

Niels
ngj is offline   Reply With Quote

Old   June 2, 2010, 05:45
Default
  #3
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Hi Niels,

Fantastic! Thank you for your quick reply!

The method I want to implement is called the Constrained Interpolation Profile (CIP) method. It is basically a subgrid reconstruction method which you can use to capture an interface. I am interested in modelling interaction between gas-liquid interfaces and particles (or drops of a third phase/species). By using CIP on an orthogonal grid I hope to be able to implement more coupled physics easier, like Marangoni effects, Eh-pH, electrophoretic forces etc.

Apologies for the noob-question , but would you have some hints on how to determine x(i,j,k) from own and nei? Thanks in advance!
linhan.ge likes this.
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 2, 2010, 05:52
Default
  #4
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
Hi Gijs

Okay, that is interesting. I did an implementation of a DDR scheme following Harvie & Fletcher earlier this year. It has some better properties with respect to keeper a clean interface, however my implementation unfortunately seems to be somewhat faulty.

With respect to neighbouring/owner, then
Code:
// Cell labels
label owner = own[cProp[facei]];
label neigh = nei[cProp[facei]];

if (owner == celli)
    // then the neighbour is actaully neigh
else // neigh == celli
  // then the neighbour is actually owner
Hope it makes sense.

Best regards,

Niels
mm.abdollahzadeh likes this.
ngj is offline   Reply With Quote

Old   June 3, 2010, 07:44
Default How to check alpha1 value on each side of a face and "grow"
  #5
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Hi Niels,

Thanks for your answer. Due to some other issues (related to the physico-chemistry we want to solve) we decided to use only part of the CIP method, so we don't really need the x(i,j,k) coordinates at this point. Sorry for the hassle about that .

To implement the (partial) CIP method I want to try transporting a solid using interFoam first, as a test. I intend to kind of freeze the solid and liwuid around it and then translate it according to the flow field. This is to prevent diffussion of the solid (which will otherwise happen with this VOF method). In order to do this I need to browse through all the faces and check whether there is on one cell center water and on the other side air. If I then find the faces that are "interface faces" I want to "grow" a (few) cell layer(s) around the solid object and then transport the "grown solid block". Then, I want to integrate the forces on the faces of the grown interface. Would you have some hints how to implement this? I have some pseudo-code below, but lots is still missing ...

Code:
// Neighbour and owner cells to given face
const unallocLabelList& nei = mesh.neighbour();
const unallocLabelList& own = mesh.owner();


// Loop through all cell faces
forAll(mesh.Sf(), facei)
{ 
    // Check for all faces whether it has on one side
    // a cell center with alpha1 = 0 and a cell center
    // with alpha1 = 1 on the other side.
    forAll( mesh.C(), celli)
    {
        // Grow a cell layer in the direction of alpha1 = 0
    }  
}
Thanks in advance!
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 7, 2010, 07:58
Default
  #6
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
Hi Gijs

Unfortunately I cannot help you with your problem, however if the object you are going to consider remain a simple geometry throughout the computation, then you might want to consider moving the boundary between the two fluids, e.g. [1].

Z. Tukovic has further developed a FAM-routine, where you solve given differential equations on a surface. As you have chemistry in your problem it might be relevant so forth some of the reactions are limited to take place on the interface.

Good luck,

Niels

[1]: Tukovic, Z. and Jasak, H. (2007). Updated Lagrangian finite volume solver
for large deformation dynamic response of elastic body. Transactions of
FAMENA, 31(1), 5570.
ngj is offline   Reply With Quote

Old   June 7, 2010, 08:42
Default
  #7
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Hi Niels,

Thanks for your help.

Quote:
Z. Tukovic has further developed a FAM-routine, where you solve given differential equations on a surface. As you have chemistry in your problem it might be relevant so forth some of the reactions are limited to take place on the interface.
Indeed, good idea. I will keep digging .
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 7, 2010, 08:51
Default
  #8
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
You are welcome. I remembered having a link to the mentioned article, see [1]. Especially section 5 might interest you.

Further, regarding FAM-framework it is already implemented in 1.5-dev and there is a tutorial where the dispersion of a species on a curved surface is modelled.

Best regards,

Niels


[1] http://powerlab.fsb.hr/ped/kturbo/Op...tionFAMENA.pdf
ngj is offline   Reply With Quote

Old   June 17, 2010, 03:37
Default
  #9
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Hi Niels,

Is there an easy way to loop through all (internal) faces in a mesh? Something like the following perhaps (although this doesn't compile)?

Code:
const unallocLabelList& listOfFaces = mesh.faces();
forAll(listOfFaces, faceI)
{
    // do stuff
}
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 17, 2010, 06:02
Default
  #10
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
Hi Gijs

I have looped over the internal faces of a surfaceScalarField, e.g. phi, when doing this. I guess you need to flux anyway for your two phase solver, so that might be a way to go.

Best regards,

Niels
ngj is offline   Reply With Quote

Old   June 17, 2010, 09:01
Default
  #11
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Hi Niels,

Right, so you're suggesting to interpolate the volScalarField onto the faces and then loop through them? I mean, phi is already the face flux, isn't it? I am dealing with cell centered values (alpha). I guess it can be done, but is it efficient?
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 17, 2010, 09:07
Default
  #12
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
What I mean is that when I have needed what you are doing, I was considering the direction of the flux over the cell faces, hence I needed to consider the sign of phi on every face relative to the orientation of the hence. This made it a obvious choice.

Could you elaborate on what it is you need?

Best regards,

Niels
ngj is offline   Reply With Quote

Old   June 17, 2010, 09:22
Default
  #13
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Sure, sorry for the fuzzyness .

I have a phase field alpha1, which is set to 1 inside a small blob and 0 in the main fluid. What I want to do is "pad" the blob with another phase, which is actually the same as the main fluid, just name differently for convenience. Let's call that phase beta. So, I want to loop through all faces and compare alpha1 on each side. When alpha1 is significantly different on one side than the other, it is an interface and I want to "pad" the adjacent main fluid as beta.

Now I have this:

Code:
forAll(alpha1, cellI)
{
    beta[cellI] = alpha1[cellI];

    label ownerCell     = mesh.faceOwner()[cellI];
    label neighbourCell = mesh.faceNeighbour()[cellI];

    // then, compare owner and neighbour
}
beta is also a volScalarField and simply a working copy of alpha1. This kind of works, but for some reason only cells below the blob are padded, not on top. So my square blob of 2x2 cells becomes a strawberry . I like strawberries, but not like this . Would you have any suggestions? Thanks in advance!
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 17, 2010, 11:33
Lightbulb Jep, face flux did the trick!
  #14
Senior Member
 
Gijsbert Wierink
Join Date: Mar 2009
Posts: 383
Rep Power: 10
gwierink is on a distinguished road
Hi Niels,

I got rid of my strawberry!

Quote:
flux over the cell faces
As you suggested, I used the face fluxes to loop over all the internal faces! phic is already there, so I used that to find a list of faces. Then, for all those face I compared alpha on each side. Now the blob is also padded "on the top"! A big thanks for the hint!
__________________
Regards, Gijs
gwierink is offline   Reply With Quote

Old   June 18, 2010, 08:16
Default
  #15
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,698
Rep Power: 27
ngj will become famous soon enoughngj will become famous soon enough
You are welcome.

Have a good weekend,

Niels
ngj is offline   Reply With Quote

Old   November 19, 2010, 05:32
Default
  #16
Member
 
Erik Arlemark
Join Date: Mar 2009
Location: Eindhoven, Netherlands
Posts: 47
Rep Power: 9
Erik is on a distinguished road
Hi all and Gijsber,

I would like to find the minimum and maximum x, y and z coordinates in my geometry while running my solver. I can see these values from the case's blockMeshDict file but I would like to avoid manually entering these scalar values by hand in the solver code.

I find that when NOT running in parallel the following lines works fine:

//values below somewhere in the middle of the domain
scalar minX=0, minY=0, minZ=0, maxX=0, maxY=0, maxZ=0;

forAll( mesh.cellCentres(), cellPrev)
{
const cell& cProp = mesh.cells()[cellPrev];
forAll(cProp, facei)
{
if (mesh.faceCentres()[cProp[facei]].x() < minX)
minX=mesh.faceCentres()[cProp[facei]].x();
if (mesh.faceCentres()[cProp[facei]].y() < minY)
minY=mesh.faceCentres()[cProp[facei]].y();
if (mesh.faceCentres()[cProp[facei]].z() < minZ)
minZ=mesh.faceCentres()[cProp[facei]].z();

if (mesh.faceCentres()[cProp[facei]].x() > maxX)
maxX=mesh.faceCentres()[cProp[facei]].x();
if (mesh.faceCentres()[cProp[facei]].y() > maxY)
maxY=mesh.faceCentres()[cProp[facei]].y();
if (mesh.faceCentres()[cProp[facei]].z() > maxZ)
maxZ=mesh.faceCentres()[cProp[facei]].z();

}
}

However, as I run these lines in a parallel simulation they do not yield the min and max coordinate values for the respective processor domain. I take it that this is because some processors do not own their boundary face, but a cell of the neighboring processor cell does, and therefor these bordering vertices are not seen.

Any suggestion for how to retrieve the min and max coordinate values is appreciated!

Cheers,
Erik Arlemark
mm.abdollahzadeh likes this.
Erik is offline   Reply With Quote

Old   November 19, 2010, 06:09
Default
  #17
Senior Member
 
Laurence R. McGlashan
Join Date: Mar 2009
Posts: 370
Rep Power: 15
l_r_mcglashan will become famous soon enough
If the mesh doesn't change during your simulation then the easiest thing to do is to have a separate 'preprocessing' application that works things like this out for you in serial. Otherwise you have to do something more complicated.
__________________
Laurence R. McGlashan :: Website
l_r_mcglashan is offline   Reply With Quote

Old   November 20, 2010, 05:20
Default
  #18
Member
 
Erik Arlemark
Join Date: Mar 2009
Location: Eindhoven, Netherlands
Posts: 47
Rep Power: 9
Erik is on a distinguished road
Thank you Laurence R. McGlashan,

I chose to get the coordinates by creating an extra dictionary in my case, where I put in all the user inputs being specific for my case.

Cheers,
Erik
Erik is offline   Reply With Quote

Old   September 12, 2011, 10:08
Default
  #19
Member
 
Santiago
Join Date: Dec 2009
Posts: 85
Rep Power: 8
gascortado is on a distinguished road
Hi all.

I'm working in developing a fluid structure interaction solver for solid-fluid. I was looking into obtaining the coordinates of each of the cells at the interface so this post helps.

While searching the forums I notice how people seem to know all kind of classes in openfoam and what they actually do. My big question is: Is there a list or guide or something where I can learn about the different classes that I can use in OF? At this point I would like to know about class that operate on the mesh.

I will appreciate any help. Thanks
gascortado is offline   Reply With Quote

Old   September 12, 2011, 14:45
Default
  #20
Senior Member
 
Mirko Vukovic
Join Date: Mar 2009
Posts: 159
Rep Power: 9
mirko is on a distinguished road
The source code documentation (see here) has a list of all the classes.

Mirko
mirko is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
spherical coordinate system in Fluent Amir FLUENT 1 March 8, 2015 13:25
Defining local coordinate system bastil OpenFOAM Pre-Processing 11 August 8, 2011 08:26
Modified OpenFOAM Forum Structure and New Mailing-List pete Site News & Announcements 0 June 29, 2009 05:56
Adventure of fisrst openfoam installation on Ubuntu 710 jussi OpenFOAM Installation 0 April 24, 2008 14:25
OpenFOAM Debian packaging current status problems and TODOs oseen OpenFOAM Installation 9 August 26, 2007 13:50


All times are GMT -4. The time now is 02:04.