CFD Online URL
[Sponsors]
Home > Forums > OpenFOAM Programming & Development

Access raw field data (scalars, vectors) on mesh

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

Reply
 
LinkBack Thread Tools Display Modes
Old   September 24, 2012, 09:05
Default Access raw field data (scalars, vectors) on mesh
  #1
Senior Member
 
Håkon Strandenes
Join Date: Dec 2011
Location: Norway
Posts: 110
Rep Power: 10
haakon will become famous soon enough
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.

Last edited by haakon; September 24, 2012 at 09:37.
haakon is offline   Reply With Quote

Old   September 24, 2012, 17:29
Default
  #2
Senior Member
 
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11
marupio is on a distinguished road
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.
__________________
~~~
Follow me on twitter @DavidGaden
marupio is offline   Reply With Quote

Old   September 25, 2012, 10:12
Default
  #3
Senior Member
 
Håkon Strandenes
Join Date: Dec 2011
Location: Norway
Posts: 110
Rep Power: 10
haakon will become famous soon enough
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.
haakon is offline   Reply With Quote

Old   September 25, 2012, 11:01
Default
  #4
Senior Member
 
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11
marupio is on a distinguished road
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.
__________________
~~~
Follow me on twitter @DavidGaden
marupio is offline   Reply With Quote

Reply

Tags
access, data, field, raw

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
Add Mesh Layers doesnt work on the whole surface Kryo OpenFOAM Native Meshers: snappyHexMesh and Others 8 September 13, 2012 10:28
How to access to field point data zxj160 OpenFOAM 5 August 14, 2012 04:32
The Field Scalars in Tascflow Darcy CFX 0 November 19, 2003 23:14
Gambit problems Althea FLUENT 21 February 6, 2001 08:05
Mesh for 3 dim Geometry Phil FLUENT 9 July 12, 2000 05:39


All times are GMT -4. The time now is 23:22.