CFD Online Discussion Forums - Access raw field data (scalars, vectors) on mesh

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Programming & Development (https://www.cfd-online.com/Forums/openfoam-programming-development/)

- - Access raw field data (scalars, vectors) on mesh (https://www.cfd-online.com/Forums/openfoam-programming-development/107328-access-raw-field-data-scalars-vectors-mesh.html)

Access raw field data (scalars, vectors) on mesh

Hi,

I am trying to develop a writer module (.so shared library) to OpenFOAM to allow me to write the mesh and solution in a different format while running the analysis (in contrast to FIRST running the analysis and THEN converting to another format). To save disk space, I have a dictionary where I specify which fields I want to print:

Code:

objectNames (U p);

I have successfully written the mesh, and that works beautifully. Now I need to print volume field data, like scalars (pressure, turbulent quantities), and vectors (velocity field). To do that I need access to the raw field data, e.g. the velocity in all cells.

I tried to follow the code used in the sampling and probes libraries, since they do something not very far from what I need. I have successfully identified the fields I give in the dictionary, and made two functions that shall write the data:

Code:

void Foam::h5Write::fieldWriteScalar()

{

    forAll(scalarFields_, fieldI)

    {

        Info<< "fieldWriteScalar: " << scalarFields_[fieldI] << endl;

    }

}

And similarly for the vector data. This compiled nicely, and the example above give the following print:

Code:

fieldWriteScalar: p

fieldWriteVector: U

However, I cannot really figure out what to do from here. How do I access the data for U and p? I need in some way to loop through all values for U and p, and then print them in my own format. The following PSEUDOCODE is what I want to do with the data (in this case a vector):

Code:

//Initialize a plain continous array for the data

scalar vectorData[vectorField.size()][3];



// Loop through the field and construct the array

forAll(vectorField, iter)

{

    vectorData[iter][0] = vectorField[iter].x();

    vectorData[iter][1] = vectorField[iter].y();

    vectorData[iter][2] = vectorField[iter].z();

}



// Write the data to disk

writeFcn(vectorData);  // Imaginary write function

Does anybody have any tips in how to achieve this?

Thanks in advance.

I'm confused. Why does "<< scalarFields_[fieldI]" output the field name? It should output the field data in list form. This makes me wonder what scalarFields_ is.

Why are you copying the data from vectorField into vectorData? Shouldn't you make your write function work with the native OF data structures?

Also, using a forAll loop over all the cells is slow.

First I will try to explain a little more:

The scalarFields_ and vectorFields_ variables are lists containing references (names) to the different scalar- and vectorfields that are going to be written. The creation of these lists are a pure copy from src/sampling/probes/probesGrouping.C and associated files.

My final goal is to have the solution written in a HDF5 (possibly also pNetCDF or CGNS) archive, rather than the one-file-per-variable-per-process-per-timestep approach. To do this, I need the data to be stored continously in memory, and then point H5Dwrite to the first element of this continous data chunk. As far as I know I cannot give H5Dwrite a pointer to OF's datastructures. Then I am left with constructing a simple array, copy the data I need into this and then pass a pointer to the first element to H5Dwrite. Please correct me if I am mistaken.

The reason for doing this is to reduce the total number of files created by OpenFOAM during the analysis, and especially when involving lagrangian particles. The hopper tutorial produces 39 files per timestep per process! That is simply too much. Further work then also includes writing particle information to this archive.

However, I have found the solution on my own (after a week with frustrating trial-and-error):

Code:

    forAll(vectorFields_, fieldI)

        Info<< "fieldWriteVector: " << vectorFields_[fieldI] << endl;

        

        const volVectorField& field = obr_.lookupObject<volVectorField>

            (

                vectorFields_[fieldI]

            );

        

        //Initialize a plain continous array for the data

        scalar vectorData[field.size()][3];

        

        // Loop through the field and construct the array

        forAll(field, iter)

        {

            vectorData[iter][0] = field[iter].x();

            vectorData[iter][1] = field[iter].y();

            vectorData[iter][2] = field[iter].z();

        }



        //

        // Initialize and create HDF5 dataset here

        // (cut away - not intersting)

        //



        // Do the actual write

        H5Dwrite

            (

                dsetID, 

                H5T_SCALAR,

                H5S_ALL,

                H5S_ALL,

                plistID,

                vectorData

            );

    }

If anyone can come up with a better solution than first using a forAll loop to create vectorData and then passing this array to H5Dwrite, please let me know.

I will of course release this work to the public after I have done some initial testing and verification.

If you replace the forAll loop with a pointer loop, it will be about 10x faster. OpenFOAM implements these in src\OpenFOAM\fields\Fields\Field\FieldM.H... but this will make your code ugly and difficult to maintain.