When developing your own code in OpenFOAM, it is almost certain that access to mesh information will be required in order to evaluate various parameters. In this post, I have provided a short description on how to access some of this information which could come in handy while programming. This post is comprised of the following topics:

Mesh connectivities


const labelListList& cellPoints = mesh.cellPoints();  // Cell to node
const labelListList& cellEdges = mesh.cellEdges();    // Cell to edge
const cellList& cells = mesh.cells();                 // Cell to face
const labelListList& cellCells = mesh.cellCells();    // Cell to cell


const faceList& faces = mesh.faces();                   // Face to node
const labelListList& faceEdges = mesh.faceEdges();      // Face to edge
const labelList& faceOwner = mesh.faceOwner();          // Face to owner cell
const labelList& faceNeighbour = mesh.faceNeighbour();  // Face to neighbour cell 


const labelListList& pointPoints = mesh.pointPoints();  // Node to node 
const labelListList& pointEdges = mesh.pointEdges();    // Point to edge 
const labelListList& pointFaces = mesh.pointFaces();    // Point to face 
const labelListList& pointCells = mesh.pointCells();    // Point to cell 


const edgeList& edges = mesh.edges();                   // Edge to node 
const labelListList& edgeFaces = mesh.edgeFaces();      // Edge to face 
const labelListList& edgeCells = mesh.edgeCells();      // Edge to cell

Mesh coordinates

Cell center, face center and node

const volVectorField& C = mesh.C();         // Cell center coordinates
const surfaceVectorField& Cf = mesh.Cf(); // Face center coordinates
const pointField& points = mesh.points(); // Node coordinates

Edge center

Since OpenFOAM is based on the conventional cell-centred Finite Volume Method, edge coordinates are not computed by the fvMesh class. However, if necessary they can be obtained very easily by averaging the 2 node coordinates associated to each edge.

// Include vectorList definition to store edge coordinates
#include "vectorList.H"

// Create edge list which contains edge to node connectivity
const edgeList& edges = mesh.edges();
// Create a vectorList to store edge centre coordinates 
// This list has size of edgeList 'edges' and is initialized to zero 
vectorList Ce(edges.size(), vector::zero);             

// Loop over all edges
forAll(edges, edge)
   const label& own = edges[edge][0];     // Index of edge owner node
   const label& nei = edges[edge][1];     // Index of edge neighbour node
   // Calculate edge center coordinates by averaging node coordinates
   Ce[edge] = 0.5*(points[own] + points[nei]);

Boundary mesh data

// Store boundary mesh information
const polyBoundaryMesh& boundaryMesh = mesh.boundaryMesh();

Distinguish between internal and boundary faces

// Loop over all faces (internal and boundary)
forAll(faces, face)
  if (mesh.isInternalFace(face)) // Internal face is found
    // Do your calculations for internal faces i.e.
// U[face] = vector::zero;
} else // Boundary face is found { const label& patch = boundaryMesh.whichPatch(face); // Boundary patch index const label& facei = boundaryMesh[patch].whichFace(face); // Local boundary face index
// Do your calculations for boundary faces i.e.
// U.boundaryField()[patch][facei] = vector::zero;


Boundary patches

// Loop over boundary patches
forAll(mesh.boundary(), patch)
  const word& patchName = mesh.boundary()[patch].name();            // Boundary patch name
  // Loop over all faces of boundary patch
  forAll(mesh.boundary()[patch], facei)
    const label& bCell = boundaryMesh[patch].faceCells()[facei];    // Boundary cell index
    const label& face = boundaryMesh[patch].start() + facei;        // Face index
// Do your calculations e.g. // U.boundaryField()[patch][facei] = vector::zero;
// Loop over all nodes of boundary patch forAll(boundaryMesh[patch].meshPoints(), pointi) {
const label& point = boundaryMesh[patch].meshPoints()[pointi]; // Node index
// Do your calculations e.g. // U[point] = vector::one;

Other useful mesh parameters

// Cell volumes
const scalarField& V = mesh.V();

// Face area normal vectors const surfaceVectorField& Sf = mesh.Sf();
// Face areas const surfaceScalarField& magSf = mesh.magSf();
// Face normals const surfaceVectorField& N = Sf/magSf;
// Total number of cells const label& nCells = mesh.nCells();
// Total number of nodes const label& nPoints = mesh.nPoints();
// Total number of internal faces const label& nInternalFaces = mesh.nInternalFaces();
// Total number of internal nodes const label& nInternalPoints = mesh.nInternalPoints();

OpenFOAM classes

For more information the following classes in OpenFOAM Foundation release version 7  can be referred:

Similar Posts


    1. Tobi, the aim of this post is to provide a brief description of the mesh information that I believe one is most likely to use while programming in OpenFOAM. I know that everything is there in the code but it is not always easy for novice OpenFOAM programmers to dig the relevant information. Nevertheless, based on your suggestion, I have added the link to the relevant class.

  1. Hello,
    This seems to be useful for my problem. I need to extract the face normal of a patch in Openfoam. Where can I put the codes you provided above in my folder case? Thanks

  2. Are there terminal commands that are equivalent to check? This would prevent from having to compile.

  3. Thanks for the concise summary. I was programming a solver that needs to check neighbor cells for a certain criteria, and found that cellCells doesn’t work in parallel if you try to address across processor boundaries, that is, the cells across decomposition boundaries aren’t included in cellCells() output. Hopefully using some of the other stuff in here will get this to work! 😀

    1. You are welcome and apologies for the delayed response! You are right because cellCells() is only aware of the cells who belong to that particular processor. If you want to access the complete mesh information then you could refer to globalMeshData class. It might have something useful for your needs.

  4. Thanks a lot for this!
    I am developing a solver which will generate boundary condition for solving Laplacian using an already calculated field. I calculated the boundary condition but I am struggling on how to use it for solving the Laplacian. Can you suggest something?

  5. Thanks a lot for your post. This was very helpful.
    Could you please help me calculate the face areas of cells in a boundary.
    scalar localArea = mesh.magSf().boundaryFieldRef()[patchi][facei];
    I have tried as above. But, am unsuccessful. [patchi] in my boundary patch.

  6. Hello,

    I found this incredibly helpful!. Thank you for putting the effort and clarifying not only the syntax but the detailed implementation. I am struggling with the idea of cell labels. It seems that all fields are saved at the center of any given cell, in the case of finite volume methods, this makes sense. However, in the PolyMesh folder there is no file such as cells. However, there is owner and neighbor information. How can I construct a connectivity diagram myself just to improve my understanding of creating meshes?.

    Thank you!

    1. Hi, I am glad that you found the post helpful. OpenFOAM is based on the cell-centred Finite Volume Method where the primary unknown variables are certainly stored at the cell centres. Then there are additional variables such as fluxes and boundary values which are generally stored at the face centres. In some cases you might have other variables which could be stored at the nodes.

      When it comes to mesh generation the information you generally require are the nodal coordinates and element connectivities. Therefore, cell centre coordinates are derived from the nodal coordinates by OpenFOAM and not part of the ‘polyMesh’ folder.

      Now the question is what are these files named ‘owner’ and ‘neighbour’? These essentially give you the connectivity of faces to cells. So every internal face will have an owner cell and a neighbour cell (generally speaking), however, boundary faces will only have one cell (owner cell) since it is at the boundary. So to cut it short, the size of owner file will correspond to the total number of faces in the mesh and each face will have an owner cell associated to it. On the other hand, the size of neighbour file will be equal to the internal faces of the mesh and each internal face will have a neighbour cell. I hope that it is clear now!

  7. Hi Jibran,

    First of all – thank you very much for this post, this helped me A LOT since I am just getting started with programming in OpenFOAM. I have some questions that I hope you could help me clarify. But first you should know – I am building my own functionObject to be able to output some desired data in a file in the end of the simulation and I used the command foamNewFunctionObject to get it going. So to my questions:

    1. I had to change all your “mesh.” statements to “mesh_” why is that? (I saw someone using “mesh_” so I tried it out and it worked)
    2. If I want to find a collection of all extensions that the boundaryMesh() has (e.g. .faceCells(), .start() etc,) where can I find that?
    3. You are writing U.boundaryField()[patch][facei] … I was trying to write that in the console because what I was expecting was that I would get an output showing me the velocity value for each face for each time step. However it gave me an error saying that “U was not declared in this scope”. How would I change it to make it work?

    I am super thankful for any help on this.

    Kind regards,

  8. Hi Haidar,

    Thanks for the post. I have a few questions for you.

    1. Look at the line below from your post:
      const edgeList& edges = mesh.edges();
      What is the first “edges”? Is it only a name you are giving to this data type? Can one use something else. I see this pattern in all of the data you created in the post.
      I do understand that the “edges” in the end of the line is a member function.
    2. I am currently studying one code in OpenFoam to better understand how I can write my own application.
      In one section of this sample application, there is this line:
      // Get the fvPatch
      const fvPatch& patch = U.boundaryField()[patchID].patch();

      I know it is trying to pull some mesh info out for a specific patch to use later, but my question is that what is the type of the data which is output from the above line? and how I can view that data for that particular patch?
  9. Hi Haidar,

    How to loop over cells in parallel run.
    I want to give ‘if’ condition for some cells. Please help me.


    1. Hi Bibin,
      OpenFOAM code is executed in exactly the same manner on all processors unless you explicitly specify it to do so. Therefore, if you loop over cells in you code, it will loop over cells in all processors separately. Now depending on what you want to implement in the if condition, additional stuff might be necessary (e.g. MPI reduce operations etc.).

      1. Hi Haider,

        I loop over all the cells inside my time loop and replace the field (e.g. T) in the cell of interest with the value I want. It works fine with the static mesh. But it doesn’t work with the adaptive mesh. Could you please tell me why and how can I make it work for the adaptive mesh? I have written the piece of code that works with the static mesh below.

        const pointField& ctrs = mesh.C();
        forAll(ctrs, cellI)
        scalar xc= magSqr(x0.value() – ctrs[cellI]);
        if (xc <= radius)
        T[cellI] = 2000;


  10. Thank you very much for your post above. It has been very useful for me. In my experiments with OpenFoam, I have a pointVectorField defined. Here in the following piece of code I am trying to access the boundary field of the defined pointVectorField. But I keep getting compilation errors. It would be very helpful for me to know how to access and modify the boundary field.
    auto& mappedPatchField =
    mappedPatchField[i_cell] = ....

    here mappedField is the pointVectorField and I want to access the field on the nodes of patch with patchId_ and modify it. When I try to compile it I get the following error message
    no match for ‘operator[]’ (operand types are ‘Foam::pointPatchField<Foam::Vector<double> >’ and ‘Foam::label {aka int}’)

    Thank you very much in advance.

  11. Hi Haidar,

    I’ve been trying to apply a nonuniform profile for my dynamic mesh with velocityComponentLaplacian but for some reason I can’t access the mesh. When I tried to construct the mesh point coordinates I got ‘mesh’ not declared in this scope, please let me know if I missed something and thanks for your time.

    Best regards,

    1. Actually please ignore my question, I realised I forgot to construct the polyMesh.

      Best regards,

  12. Hello,

    I am building my own finite volume solver but I want to use the mesh data from OpenFOAM.
    I make the .unv mesh in Salome and then convert it to an Foam mesh using the command ideasUnvToFoam.

    Now, I would like to extract the cell-to-node connectivity matrix. How do I do that?

    Where do I need to input the following statement?

    const labelListList& cellPoints = mesh.cellPoints(); // Cell to node

    I tried in controlDict but it does not work.



      1. Thanks for the reply.

        But I do not want to solve the equations, I only need to do use the above commands as a preprocessing step. How do I write the subroutine in such a way that mesh data is extracted without running the simulation?


  13. Hi Haidar
    Thanks for the post.
    I want to access volume of cells in volScalarField but in openFoam there is two ways to access to the mesh volume and none of them are volScalarField:
    const DimensionedField& V() const;
    const scalarField& cellVolumes() const;

    1. Hi Mojtaba,

      Glad that the post helped you. May I ask why would you want to store volume in a volScalarField?

      In OpenFOAM volFields have an internal field and a boundary field. Since cell volumes are only associated to the cells themselves they are not associated to a volField.
      Depending on your need, you could always work around and create a volScalarField whose internalField() is filled with cellVolumes().

      1. Thank you very much for your reply.

        I want to add a source term in my continuity equation. In this source term, there is a constant called “specific area” which has a value of non-zero on the boundary and zero value in other locations. My problem is 2D, and the value of the specific area is defined as the ratio of the area of the face, connected to the boundary on a cell boundary, to the volume of the boundary cell. I was wondering how could I program this definition to calculate the specific area, and which type of specific area should I add to the continuity equation?

        Thank you very much in advance!

        1. I think what you need is mentioned in the Boundary patches part of this post.
          You can loop over boundary patches and then faces of those patches.
          You will have access to the global face index of this patch face and also the boundary cell index.
          The code snippet is repeated below for convenience.

          // Loop over boundary patches
          forAll(mesh.boundary(), patch)
          const word& patchName = mesh.boundary()[patch].name(); // Boundary patch name

          // Loop over all faces of boundary patch
          forAll(mesh.boundary()[patch], facei)
          const label& bCell = boundaryMesh[patch].faceCells()[facei]; // Boundary cell index
          const label& face = boundaryMesh[patch].start() + facei; // Face index

          // Do your calculations e.g.
          // U.boundaryField()[patch][facei] = vector::zero;
          } // End loop over faces of boundary patch
          } // End loop over boundary patches

  14. The post is great, man!

    I have a boundary patch with a polygon shape. But with the fvPatch instance, I could not find a way to obtain the edges of this polygon. Could you please kindly share your thoughts on this issue? Thanks a lot!

  15. Hi,

    Thanks for the post!

    A question, should the the surface normal vector be a reference?
    Shouldn’t it be:
    const surfaceVectorField nf = Sf/magSf;

  16. Hi,
    Thank you for your post !! I have a query to ask you.
    My domain is divided into N elements (I= 1 to N). In my work, I need to duplicate the node as (IB, IC) for a purpose. Let’s say,
    for(I = 1 to N)
    if (T[I] == 100)
    IB = I;
    for(I=IB+1 to N)
    if (T[I] < 100)
    IC = I;
    (1) How to declare the for loop specifically in OpenFOAM?
    (2) How to address IB and IC ? Like if its scalar variable we can declare "scalar IB" but that's not gonna work here.
    Kindly someone share the ideas.
    Thank you

  17. Hello!

    My question is regarding ‘multiRegion’ in OpenFOAM. Suppose I have two regions – solid and liquid.
    The two regions share a common boundary, known as the coupled patch. From being in the solid region or vice versa, through which command can I find out the patchID of the coupled patch?

Leave a Reply

Your email address will not be published. Required fields are marked *