CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Post-Processing (https://www.cfd-online.com/Forums/openfoam-post-processing/)
-   -   Write cells and data intersecting a plane cuttingPlane (https://www.cfd-online.com/Forums/openfoam-post-processing/61321-write-cells-data-intersecting-plane-cuttingplane.html)

dmoroian January 28, 2008 10:57

Hi, I would like to save durin
 
Hi, I would like to save during the run time, data (and cells) that intersect a plane. I saw the cuttingPlane class, but is there a code snippet on the web showing how to save those cells in a top level solver?
I want later on to be able to read those files with paraFoam or paraview.

Dragos

ville January 28, 2008 11:48

Hi Dragos, yes, check this th
 
Hi Dragos,
yes, check this thread. In case you need further
advice please continue the thread and I'll help
you. The cutting plane is very practical if you
like to take your data to matlab, gnuplot etc for
further analysis.

http://www.cfd-online.com/OpenFOAM_D...tml?1197913395

best regards,
Ville

dmoroian January 29, 2008 02:14

Hi Ville, Actually your answe
 
Hi Ville,
Actually your answer was the source of my information in the first place. There you show how to access the cells and data, but I am interested in saving them in openfoam format. Can you show me how to do that?
If you help me, I promise to put toghether a wiki page so other people would benefit from it.

Dragos

ville January 29, 2008 04:41

Hi, unfortunately I doubt if
 
Hi,
unfortunately I doubt if that is possible since
the cutting plane option is needed to extract data from the mesh so that connectivity information (which the post processor should need) is lost..
-Ville

mattijs January 29, 2008 04:57

Get the cells from the cutting
 
Get the cells from the cutting plane and construct a cellSet from them. Something like

const labelList& cutCells = cutPlane.cells();

cellSet someCells(mesh, "someCells", cutCells);

someCells.write();

Then use foamToVTK with -cellSet option.

dmoroian January 29, 2008 14:58

Thanks both of you! This is w
 
Thanks both of you!
This is what I was after.

dmoroian January 30, 2008 01:40

...me again. The following co
 
...me again.
The following code works very well:

point pnt(0.5,0.25,-0.25);
vector direction(1,0,0);
plane pl1(pnt,direction);
cuttingPlane cutPlane(mesh,pl1);
const labelList& cutCells = cutPlane.cells();
cellSet someCells(mesh, "someCells", cutCells);
someCells.write();


Here is a slice obtained with it:
http://www.cfd-online.com/OpenFOAM_D...ges/1/6528.png
... but I would like more. Is it possible to save together with the above cell set also the data associated to it at a certain time?

Dragos

jaswi January 30, 2008 10:17

Hi Dragos Thanks for the co
 
Hi Dragos

Thanks for the code.

Could you also please post the list of include files so that this code compiles successfully. I put the code in my post processing routine exactly as you have psoted


point pnt(0.0,0.0,0.04);
vector direction(1,0,0);
plane pl1(pnt,direction);
cuttingPlane cutPlane(mesh,pl1);
const labelList& cutCells = cutPlane.cells();
cellSet someCells(mesh, "someCells", cutCells);
someCells.write();

and get the following error:

postProcessing.C:654: error: 'plane' was not declared in this scope
postProcessing.C:654: error: expected `;' before 'pl1'
postProcessing.C:655: error: 'cuttingPlane' was not declared in this scope
postProcessing.C:655: error: expected `;' before 'cutPlane'
postProcessing.C:656: error: 'cutPlane' was not declared in this scope
postProcessing.C:657: error: 'cellSet' was not declared in this scope
postProcessing.C:657: error: expected `;' before 'someCells'
postProcessing.C:658: error: 'someCells' was not declared in this scope
postProcessing.C:656: warning: unused variable 'cutCells'


The list of files I have already included are:

#include "fvCFD.H"
#include "wallFvPatch.H"
#include "interfaceProperties.H"
#include "twoPhaseMixture.H"
#include "incompressible/turbulenceModel/turbulenceModel.H"
#include "nearWallDist.H"
#include "fixedValuePointPatchFields.H"


Thanks alot

With Kind Regards
Jaswinder

jaswi January 30, 2008 10:37

Hi Dragos I just checked o
 
Hi Dragos

I just checked out the cuttingPlane code and it has a member function called


//- Sample the cell field
template<class>
tmp<field<type> > sample(const Field<type>&) const;

It seems that if you add the following to your code snippet

point pnt(0.0,0.0,0.04);
vector direction(1,0,0);
plane pl1(pnt,direction);
cuttingPlane cutPlane(mesh,pl1);
const labelList& cutCells = cutPlane.cells();

----------------
volScalarField slicedDesiredField = cutPlane.sample(desiredField);
-----------------

cellSet someCells(mesh, "someCells", cutCells);
someCells.write();

You might be able to extract the field you want.

I could not test whether if it works as i am unable to compile the original code offered by you.

With Kind Regards
Jaswinder

jaswi January 30, 2008 11:51

Hi Dragos I have solved the
 
Hi Dragos

I have solved the compilation problem but now I have problem visualizing it.

I give the following command:

foamToVTK <root> <case> -cellSet someCells

and it writes


Exec : foamToVTK . 7_flask_10ml_320rpm -cellSet someCells
Date : Jan 30 2008
Time : 17:48:12
Host : taifun
PID : 5343
Root : <root>
Case : <case>
Nprocs : 1
Create time

Number of cells in new mesh : 6034
Number of faces in new mesh : 24264
Number of points in new mesh: 12342


--> FOAM FATAL IO ERROR : cannot open file

file: /scratch1/power_consumption_sims/1st_10ml/7_flask_10ml_320rpm/system/subsetSubse t/fvSchemes at line 0.

From function regIOobject::readStream(const word&)
in file db/regIOobject/regIOobjectRead.C at line 66.

FOAM exiting


Could you please give me a clue what is missing

With kind Regards
Jaswinder

dmoroian January 30, 2008 15:42

Hi Jaswinder, Thanks for the
 
Hi Jaswinder,
Thanks for the information, I'll definitely try it (I have to check if volScalarField class can write down the data). On the other hand I have no idea why you get that error.

Dragos

mattijs January 30, 2008 16:09

foamToVTK .. cavity -cellSet c
 
foamToVTK .. cavity -cellSet c0

works for me (1.4.1):

Time constant

Internal : "../cavity/VTK/c0_0.vtk"
Original cells:5 points:24 Additional cells:0 additional points:0

Patch : "../cavity/VTK/movingWall/c0_0.vtk"
Patch : "../cavity/VTK/fixedWalls/c0_0.vtk"
Patch : "../cavity/VTK/frontAndBack/c0_0.vtk"
Patch : "../cavity/VTK/oldInternalFaces/c0_0.vtk"

dmoroian January 31, 2008 01:46

Hi Mattijs, Will you give me
 
Hi Mattijs,
Will you give me a clue on how to get only the slice (both geometry and data)? I don't want to save the entire domain and then to get the slice.

Dragos

mattijs January 31, 2008 05:04

pull out the bits of code from
 
pull out the bits of code from subsetMesh. It takes a cellSet and creates a whole new fvMesh from it and subsets the cells as well.

dmoroian February 1, 2008 10:40

Ok, new questions... In the s
 
Ok, new questions...
In the subsetMesh utility the following lines read the mesh and the volume scalars from a specific time directory:

IOobjectList objects(mesh, runTime.timeName());
wordList scalarNames(objects.names(volScalarField::typeName ));
PtrList<volscalarfield> scalarFlds(scalarNames.size());


then a subset of it is constructed by

subsetVolFields(subsetter, scalarNames, scalarFlds);

Now my question is: how do I construct the subset from the variables I have in memory (during the calculation) and not read them from a dictionary?

Dragos

dmoroian February 5, 2008 05:41

wiki: TurbPlaneCutFoam
 
wiki: TurbPlaneCutFoam

dmoroian February 5, 2008 05:59

Hmm, I need an idea: how to ma
 
Hmm, I need an idea: how to make it run in parallel?
Anybody can give me a clue?
This is the log file with the error I get when I run the solver in parallel:
http://www.cfd-online.com/OpenFOAM_D...hment_icon.gif channel_par.log

Dragos

dmoroian February 6, 2008 09:25

Hello again, Running in paral
 
Hello again,
Running in parallel and getting the slices saved, seems to be ok. But how should I reconstruct the plane back from each processor? For instance, two processors will give me two pieces like in the following image:
http://www.cfd-online.com/OpenFOAM_D...ges/1/6593.png
How should I put them together?
reconstructPar doesn't seem to work for the cells saved from the cuttingPlane.
I can merge the mesh using mergeMeshes, but how can I merge the data?

Dragos

braennstroem February 26, 2008 13:06

Hi Dragos, I tried your cod
 
Hi Dragos,

I tried your code snip, but have some problems compiling it:

-lincompressibleLESmodels -lincompressibleTransportModels -lfiniteVolume -lmeshTools -lOpenFOAM -liberty -ldl -lm -o /home/gcae504/OpenFOAM/gcae504-1.4.1/applications/bin/linux64GccDPOpt/no_avg_ood les
Make/linux64GccDPOpt/no_avg_oodles.o(.text+0x31b5): In function `main':
: undefined reference to `Foam::cuttingPlane::cuttingPlane(Foam::primitiveM esh const&, Foam::plane const&)'
Make/linux64GccDPOpt/no_avg_oodles.o(.text+0x31c2): In function `main':
: undefined reference to `Foam::cuttingPlane::cells() const'
collect2: ld returned 1 exit status
make: *** [/home/gcae504/OpenFOAM/gcae504-1.4.1/applications/bin/linux64GccDPOpt/no_avg_oo dles] Error 1

I added these libraries to the oodles solver:

#include "plane.H"
#include "cuttingPlane.H"
#include "cellSet.H"
#include "fvMeshSubset.H"

Would be nice, if you have a hint, how to compile it!?

Regards!
Fabian

dmoroian February 27, 2008 05:09

Hi Fabian, It seems that you
 
Hi Fabian,
It seems that you have all the correct "#include" directives, since the error is from the link stage and not from compilation. What you basically need is the "-lsampling" option. Here is how my options file looks like:

EXE_INC = \
-I$(LIB_SRC)/LESmodels \
-I$(LIB_SRC)/LESmodels/LESdeltas/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude

EXE_LIBS = \
-lincompressibleLESmodels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lmeshTools \
-lsampling


dmoroian February 27, 2008 06:27

...and by the way, to close th
 
...and by the way, to close the problem: in order to put back results from a parallel computation, I used mergeMeshes for the mesh and mapFields for the fields:
http://www.cfd-online.com/OpenFOAM_D...ges/1/6828.png

braennstroem February 28, 2008 14:18

Hi Dragos, thank! But how c
 
Hi Dragos,

thank! But how can one know this from the error message?

Fabian

dmoroian February 29, 2008 01:44

Hi Fabian, Well the error com
 
Hi Fabian,
Well the error comes from the linker because it belongs to an object file not to a source file:
Make/linux64GccDPOpt/no_avg_oodles.o(.text+0x31c2): In function `main':
: undefined reference to `Foam::cuttingPlane::cells() const'

Then, a little bit trickier is to find the library that contains the missing function (in this case cuttingPlane::cells()). For that I first go to the source directory of OpenFOAM, and then search for the class cuttingPlane:

home2/dragosm> src
OpenFOAM-1.4.1/src> find ./ -name cuttingPlane.C
./sampling/cuttingPlane/cuttingPlane.C
./sampling/lnInclude/cuttingPlane.C

...so it is in the sampling directory. Go there and list it:

OpenFOAM-1.4.1/src> cd sampling/
src/sampling> ls
cuttingPlane
include
Make
probes
graphField
lnInclude
meshToMeshInterpolation
samplingLine

As seen above, a Make directory exists. Let see what does it produce:

src/sampling> cat Make/files
probes/probes.C
probes/probesFunctionObject/probesFunctionObject.C
probes/IOprobes.C

samplingLine/pointAddressing.C
samplingLine/samplingLine.C

cuttingPlane/cuttingPlane.C

graphField/writePatchGraph.C
graphField/writeCellGraph.C
graphField/makeGraph.C

meshToMesh = meshToMeshInterpolation/meshToMesh
$(meshToMesh)/meshToMesh.C
$(meshToMesh)/calculateMeshToMeshAddressing.C
$(meshToMesh)/calculateMeshToMeshWeights.C

LIB = $(FOAM_LIBBIN)/libsampling


So all it generates is a library called libsampling. In this library you will find the requested function, which means that you need to include it in your own options file, if you want the compilation to be completed.

I hope this is usefull,
Dragos

braennstroem February 29, 2008 08:06

Hi, great, thanks for the e
 
Hi,

great, thanks for the explanation!

Fabian

fabian_korn March 31, 2008 06:23

Hi Dragos, i am writing a
 
Hi Dragos,

i am writing a Bachlor thesis at Lth, i have a probelm with the code you posted.
I am really new with openFOAM and i am running a LES over a 2 cylinder flow.

My problem is, in which file should i ionclude your code?
And did i understand it correctly, that it creates an output at every timestep for one defined slice of the field for all the variables

Thanks for help

Fabian

dmoroian March 31, 2008 06:48

Hi Fabian, The code goes in t
 
Hi Fabian,
The code goes in the main top solver. For simplicity you can put it entirely in the main time loop. That means it will construct and save the slice every time step.
It is true, it saves every time step but only the variables that you tell it to save. For instance, if you want pressure you have to use something like:
scalarNames[0] = "p"

Dragos

braennstroem March 31, 2008 07:23

Hi Dragos, do you have a hi
 
Hi Dragos,

do you have a hint, how to write out more than one plane at a time? I tried using the same code twice with different names, which does not work. I think 'subsetter' makes some problems!?

Regards!
Fabian

dmoroian March 31, 2008 08:07

Hm, it is strange! I would exp
 
Hm, it is strange! I would expect it to work. You have 2 subsetters, right? One for each slice, as well as with the rest of the variables.

braennstroem March 31, 2008 10:03

Yes, I have 2 subsetters like
 
Yes, I have 2 subsetters like this:

#--------------------------------------------
// Plane 1
point pnt(0,25,0);
vector direction(0,1,1);
plane pl1(pnt,direction);
cuttingPlane cutPlane1(mesh,pl1);
const labelList& cutCells1 = cutPlane1.cells();
word setName("someCells");
cellSet currentSet1(mesh, setName, cutCells1);

// Create mesh subsetting engine
fvMeshSubset subsetter1(mesh);
label patchI = -1;
subsetter1.setLargeCellSubset(currentSet1, patchI, true);


wordList scalarNames(1);
scalarNames[0] = "p";
PtrList<volscalarfield> scalarFlds1(scalarNames.size());

wordList vectorNames(1);
vectorNames[0] = "U";
PtrList<volvectorfield> vectorFlds1(vectorNames.size());

// Plane 2
point pnt2(0,25,25);
vector direction2(0,1,1);
plane pl2(pnt2,direction2);
cuttingPlane cutPlane2(mesh,pl2);
const labelList& cutCells2 = cutPlane2.cells();
word setName2("someCells2");
cellSet currentSet2(mesh, setName2, cutCells2);

// Create mesh subsetting engine
fvMeshSubset subsetter2(mesh);
label patchI2 = -1;
subsetter2.setLargeCellSubset(currentSet2, patchI2, true);

//wordList scalarNames(1);
scalarNames[0] = "p";
PtrList<volscalarfield> scalarFlds2(scalarNames.size());

//wordList vektorNames(1);
vectorNames[0] = "U";
PtrList<volvectorfield> vectorFlds2(vectorNames.size());



...


scalarFlds1.set(0, subsetter1.interpolate(p));
vectorFlds1.set(0, subsetter1.interpolate(U));
scalarFlds2.set(0, subsetter2.interpolate(p));
vectorFlds2.set(0, subsetter2.interpolate(U));

Info<< "Writing subsetted mesh and fields to time " << runTime.value()
<< endl;
subsetter1.subMesh().write();
subsetter2.subMesh().write();
forAll(scalarFlds1, i)
{
scalarFlds1[i].write();
}
forAll(vectorFlds1, i)
{
vectorFlds1[i].write();
}

forAll(vectorFlds2, i)
{
vectorFlds2[i].write();
}

#--------------------------------------------


but I just get one subsetP and one subsetU in my time directory. Strange...

dmoroian March 31, 2008 10:57

Both slices have the same name
 
Both slices have the same name: subsetp and subsetU, and they get overwritten.
You have to rename them for one of them (let say subset1p and subset1U) before you save them with scalarFlds1[i].write() and vectorFlds1[i].write().
For renaming check the ~/OpenFOAM/OpenFOAM-1.4.1/applications/utilities/mesh/manipulation/subsetMesh/su bsetMesh.C from line 242 and below.

Dragos

braennstroem March 31, 2008 11:36

Thanks for the file hint, but
 
Thanks for the file hint, but it does not work due to the different number of cells in the two planes, because it writes just one 'polyMesh' for both planes. Maybe, there is a chance to write it in vtk format, just like in sampleSurf!? Or even 'better' one somehow integrates sampleSurf; then one could use the sampleSurfDict as well...

Fabian

fabian_korn March 31, 2008 11:51

Hi Dragos, thanks for the fa
 
Hi Dragos,
thanks for the fast answer,
Now a second question how did you mamage it to reconnect the meshes after a parrallel run with 4 processors. I saw you used the mergeMesh command, but how does it work?
Hopefully the question is not to stupid.

Thanks Fabian

dmoroian April 1, 2008 02:05

Hello Fabian B. Indeed you ar
 
Hello Fabian B.
Indeed you are right, the polyMesh gets rewritten too. Beside your excellent ideas, another one would be to collect both planes in one set and write that down, maybe write sets information to be able to extract later the planes separately.

Hi Fabian K.
To put together two pieces of mesh you can use the mergeMesh. Check this answer on how to do it: mergeMesh. When you will be ready with that you have to "merge" the data too. For that you have to use mapFields... but first fix the merging part.

Dragos

fabian_korn April 2, 2008 02:37

Hi, hopefully the last que
 
Hi,

hopefully the last question, during compiling i got the follwing massage, i gues caused by PtrList. Here is the the entry of oodles.C

wordList scalarNames(1);
scalarNames[0] = "c";
PtrList<volscalarfield>

and the Error massage

oodles.C: In function 'int main(int, char**)':
oodles.C:98: error: 'volscalarfield' was not declared in this scope
oodles.C:98: error: template argument 1 is invalid
oodles.C:98: error: invalid type in declaration before '(' token
oodles.C:102: error: 'volvectorfield' was not declared in this scope
oodles.C:102: error: template argument 1 is invalid
oodles.C:102: error: invalid type in declaration before '(' token
oodles.C:98: warning: unused variable 'scalarFlds1'
oodles.C:102: warning: unused variable 'vectorFlds1'
make: *** [Make/linuxGccDPOpt/oodles.o] Error 1


Thank you again

dmoroian April 2, 2008 03:07

Hi Fabian, On short: volscala
 
Hi Fabian,
On short: volscalarfield -> volScalarField
I don't know why is written in the wrong way all over the thread, because is wrong, and when I put the above example I've just extracted it from a working code.

Dragos

fabian_korn April 2, 2008 06:10

Hi again, it is possible n
 
Hi again,

it is possible not to write every timestep, but let us say every 100 timestep?

fabian_korn April 2, 2008 08:13

Sorry for asking again, i h
 
Sorry for asking again,

i hope it is not too much. When compiling int i get

oodles.C: In function 'int main(int, char**)':
oodles.C:111: error: 'samplingCount' was not declared in this scope
make: *** [Make/linuxGccDPOpt/oodles.o] Error 1

I tried to find a library to find an idea how to fix it, but i do not have any guess.

Thanks again

dmoroian April 2, 2008 08:53

Don't appologise, if I'm getti
 
Don't appologise, if I'm getting disturbed I won't write http://www.cfd-online.com/OpenFOAM_D...part/happy.gif
The solution: you have to declare this variable somewhere; one place is imediately after the main function.
Quote:

main(...){
label samplingCount = 0;
...

Dragos

fabian_korn April 2, 2008 11:21

Hi Dragos, first hte good t
 
Hi Dragos,

first hte good thing, it compiles.
But when i run my case this happens.


Courant Number mean: 0.00138039 max: 0.0242767
DILUPBiCG: Solving for Ux, Initial residual = 0.00148907, Final residual = 3.64592e-09, No Iterations 1
DILUPBiCG: Solving for Uy, Initial residual = 0.00176973, Final residual = 7.19327e-09, No Iterations 1
--> FOAM Warning :
From function Foam::cuttingPlane::walkCell
in file cuttingPlane/cuttingPlane.C at line 189
Did not find closed walk along surface of cell 1646 starting from edge 9894 in 0 iterations.
Collected cutPoints so far:1(88)
--> FOAM Warning :
From function Foam::cuttingPlane::walkCell
in file cuttingPlane/cuttingPlane.C at line 189
Did not find closed walk along surface of cell 1647 starting from edge 9911 in 1 iterations.
Collected cutPoints so far:2(89 38)
DICPCG: Solving for p, Initial residual = 0.0217737, Final residual = 0.00108258, No Iterations 9
time step continuity errors : sum local = 1.70955e-10, global = -7.82466e-21, cumulative = 4.23193e-20
DICPCG: Solving for p, Initial residual = 0.00648299, Final residual = 8.2598e-07, No Iterations 93
time step continuity errors : sum local = 1.28292e-13, global = -1.27384e-21, cumulative = 4.10454e-20
ExecutionTime = 4.37 s ClockTime = 5 s

it generates a solution, but it slowes down the simulation and if i run my real big case, than the list of Collected cutPoints so far:2(89 38) is really long.

Any guess?

dmoroian April 3, 2008 01:32

It happens the same thing in m
 
It happens the same thing in my case. Just a wild guess is that some cells have a face exactly in that plane and the solver doesn't know which of the two cells sharing the face should be picked. To speed up the process, I separate the cell set construction from the interpolation and writing. I put the construction part. For instance, if you take the code from cuttingSlice I put everything above scalarFlds.set(0, subsetter.interpolate(c)) outside the time loop, and kept the rest inside. In this way I get only once those warnings.

I hope this is helpful,
Dragos


All times are GMT -4. The time now is 10:19.