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

Distance of a 'North' and 'South' surface of a cell for 2D mesh

Register Blogs Community New Posts Updated Threads Search

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   October 6, 2020, 10:22
Default Distance of a 'North' and 'South' surface of a cell for 2D mesh
  #1
Senior Member
 
René Thibault
Join Date: Dec 2019
Location: Canada
Posts: 114
Rep Power: 6
Tibo99 is on a distinguished road
Hi everyone!

My problem could be easy for some people, but I personally don't have any idea on how to code this thing here.

I would like to calculate, for a 2D mesh, the distance 'z' from the wall to the 'north' (Zn) and 'South' (Zs) faces of a cell (see image in attachment).

Im able to obtain the distance for a center of a cell (ZN, ZP, ZS) with 'mesh_.C' like so:
Code:
const vectorField& CellC = mesh_.C();
But I don't have any clue on how to get the distance 'z' for the 'North' and 'South' faces of the same cell. What I have understand so far, and I'm maybe wrong, is that I need to find from a particular cell the indices what are the 'North' and 'South' surface are compose of and calculate the normal vector in order to to calculate the height 'z' after.

But I don't have any idea whatsoever how to code this thing.

Can I have some example and some help?

Thank you everyone!
Attached Images
File Type: png Cells.png (2.7 KB, 15 views)

Last edited by Tibo99; October 6, 2020 at 14:52.
Tibo99 is offline   Reply With Quote

Old   November 27, 2020, 05:58
Default
  #2
New Member
 
Join Date: Feb 2019
Posts: 23
Rep Power: 7
dolfus is on a distinguished road
Hi,

You can loop over all cell faces and take into account only the faces which have the z component of its normal vector different from 0. Then, see which face is the N and which is S. Finally, just take the distance from the wall (assuming z=0) to the face centers.
The code could be something like (for cell ID = 10, for example):
Code:
const cellList& cellFaces = mesh.cells();
const surfaceVectorField& Sf = mesh.Sf();
const surfaceVectorField& Cf = mesh.Cf();
const surfaceVectorField& C = mesh.C();
scalar zWall = 0.0;
label cellI = 10;
scalar zn = 0.0;
scalar zs = 0.0;

forAll(cellFaces[cellI],faceI)
{
   label faceIG = cellFaces[cellI][faceI];
   scalar zSf = Sf[faceIG][2];
   if (zSf != 0.0)
   {
      if (Cf[faceIG][2] > C[cellI][2]) 
      {  
          zn = Cf[faceIG][2] - zWall;
      }
      if (Cf[faceIG][2] < C[cellI][2]) 
      {  
          zs = Cf[faceIG][2] - zWall;
      }
   }
 }
Notice that I haven't taken into account boundary faces. If so, the code would be slightly different.
dolfus is offline   Reply With Quote

Old   November 27, 2020, 09:06
Default
  #3
Senior Member
 
René Thibault
Join Date: Dec 2019
Location: Canada
Posts: 114
Rep Power: 6
Tibo99 is on a distinguished road
Thank you very much for replying!

The fact that the code don't take into account any boundary faces in the domain is not that a big deal for me since I use it, for the moment, for a 2D channel.

My 1st question is, because I need to apply this from the 2nd cell to the 4th cells from the wall in the entire domain, how I could modify the loop to take this into account?

The 2nd one is, is this option work with parallel computing?

I was thinking using 'topoSet' to get the ID for these cells and then, I would use your loop inside a new loop in order to apply it to all of the cells that 'topoSet' has selected.

If this approach is right, I don't know at this point how to get all of these cells ID, stored them in a constant ''IDList'' (for instance) and using it after in the new loop.

I keep working on it.

Even though, my main issue is pretty much solved with the loop you suggested.

Thank you very much again and be safe.

Regards,

Last edited by Tibo99; November 29, 2020 at 19:32.
Tibo99 is offline   Reply With Quote

Old   November 28, 2020, 06:04
Default
  #4
Senior Member
 
Mark Olesen
Join Date: Mar 2009
Location: https://olesenm.github.io/
Posts: 1,685
Rep Power: 40
olesen has a spectacular aura aboutolesen has a spectacular aura about
Your use case is fairly specific, but could still be worth generalizing. Starting with a simple struct with four label members (north, south, east, west) that will hold the cellId of the respective neighbours. You could also just use a FixedList of 4 labels for this purpose.
Create a List of these structs that is nCells long, with each element initialized to -1.
Code:
enum compass { NORTH = 0, SOUTH, EAST, WEST };
typedef FixedList<label, 4> compassType;

List<compassType> compassAddr(mesh.nCells(), compassType(-1));
Next you want to loop across each mesh face. Internal faces have own/neighbour cells, boundary faces only have an owner cell. From the face you get its normal. Since you only care about the general direction, just use the area normal (possibly slightly cheaper). You don't need the scalar (dot) product, but just find the cmptMax (the max component) of the face normal. With this you know if this is east-west or north-south. The combination of the sign of that component and the owner/neighbour relationship tells you which sub elements of the respective compassAddr entries to set.
Now that the addressing is complete, you can now handle everything else that you need.
A simple example (untested)
Code:
// cellCentres() for polyMesh
const pointField& cc = mesh.C();

scalar dist = 0;
label celli = someCellID;

// walk four to the east
for (label i=0; i < 4; ++i)
{
    const label next = compassAddr[celli][compass::EAST]:
    if (next == -1) break;
    dist += mag(cc[next] - cc[celli]);
    celli = next;
}
I think this is enough for you to work with.
olesen is offline   Reply With Quote

Old   November 28, 2020, 06:16
Default
  #5
Senior Member
 
Mark Olesen
Join Date: Mar 2009
Location: https://olesenm.github.io/
Posts: 1,685
Rep Power: 40
olesen has a spectacular aura aboutolesen has a spectacular aura about
Just re-read your first post and noticed that you actually wanted the centre to face distances. In this case you need to store the face id instead of the cell id in the north/south/east/west slots. This will give you an immediate calculation of centre to face, but when walking to the next cell you will need to use the mesh face owner/neighbour information for the connectivity and to walk in the correct direction.
This is probably more flexible than what I originally posted.
olesen is offline   Reply With Quote

Old   November 29, 2020, 13:29
Default
  #6
Senior Member
 
René Thibault
Join Date: Dec 2019
Location: Canada
Posts: 114
Rep Power: 6
Tibo99 is on a distinguished road
Hi Olesen,

I maybe should have mentioned in the first post that I still considered myself at a beginner level with OF.

Regarding the fact that the solution you suggested is generalized and probably the best one for my application, I have to admit, I don't exactly understand the code. I feel that this is an another level for me, which I'm not there yet. Sorry about that.

But, here is a thread that I posted and it's related to. That should help everyone to understand why I asked how to code this.

Applying correction to the k-e transport equation

Thank you very much for your time and be safe!

Last edited by Tibo99; November 30, 2020 at 10:57.
Tibo99 is offline   Reply With Quote

Reply


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 Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
[snappyHexMesh] Edge refinement ashghan OpenFOAM Meshing & Mesh Conversion 4 May 13, 2014 05:45
lid-driven cavity in matlab using BiCGStab Don456 Main CFD Forum 1 January 19, 2012 15:00
Mass flux for west, south, bottom cell face marcus Siemens 0 April 19, 2006 03:43
hex vs. tet grids Tom Plikas Main CFD Forum 3 June 9, 1999 02:46


All times are GMT -4. The time now is 06:40.