dispersion model with lagragian particle tracking model for incompressible flows
 

Dear All,

I started to explore the Lagrangian particle tracking capabilities of OF 2.1 last week. I would like to do particle tracking of non-reacting particles in an incompressible flow with turbulent dispersion. Standard forces like drag and buoyancy should be included. It seems that there is no standard model for this, as most models which come for OF are developed for reacting flows in compressible fluids. Therefore I tried to implement it myself. I have included the solver (pimpleCoupledKinematicParcelFOam ) + 2 test cases to this post. The solver works, except that I can not use the dispersionModels! I hope to find out what is wrong by opening this thread. First I will give some back ground of my approach and then I end with my question.

As a starting point I took the isoUncoupledKinematicParcelFoam solver in OF 2.1.0, which uses the kinematic particle foam class to solve the particle transport. As I understand, the particle tracking models in OF are all developed for compressible flow. This means that dt(rho,U) is solved in stead of dt(U), and moreover, the pressure is expressed in kg/m/s2 (Pa). The isoCoupledKinematicParcelFoam solver does not solve for U at all, so no problem here. However, I wanted to make a solver which also includes the U and p solution. If you formulate the solver as dt(rho,U), the turbulence modelling can not be included, because these require the pressure to have the unit m2/s2 (i.e. pressure normalized with density). Therefore I create a new pressure field p_rgh with units m2/s2 (just as the pimple solver), which I use in the solution of d(U)/dt and p_rgh, and calculate the p=rho*(rg+p_rgh) pressure which is the pressure in N/m2 including the hydrostatic pressure. This tricks works: you can include turbulence in the solution and still be able to call the kinematicParticle foam routines.

Well, here my version of the solver, which I call pimpleCoupledKinematicParcelFOam. It works and it allows you to use all functionality in the kinematic particles. It is actually very similar with the solver submitted michaelb in http://www.cfd-online.com/Forums/ope...-tracking.html, except that Michaels solver is only one-way coupled (it does not feed back the particle force on the fluid) and does not include the gravity term. So this solver is a bit more generic. I have taken Michaels test case of the 3D cavity with moving lid and show that the particle can be solved qualitatively correct. See the graph for a mass-less particle tracking and a particle with mass.

There is one problem which I can not solve though: if I want to use the dispesionModel in the kinematicParticleProperties dictionary, I get the following error:

lookup of RASProperties from objectRegistry region0 succesful, but is not RASModel, it is kEpsilon.

It seems to me that the dispersionModel again expects a turbulence model which is based on a compressible flow, whereas now in the solver the turbulence model is based on an incompressible flow.

My question is: can this be easily solved by tricking the dispersionModel, or do I need to rewrite the whole particle tracking class such that it is based on incompressible turbulence models?

Has somebody already implemented this Lagrangian particle tracing in transient incompressible flow which turbulent dispersion? Any hints or suggestions appreciated!

Regards
Eelco

Dear eelcovv,

I'm also dealing with lagrangian particle tracking.
I would like to use the stochasticDispersionRAS model with the icoUncoupledKinematicParcelFoam solver, and I get the same error :

--> FOAM FATAL ERROR:

lookup of RASProperties from objectRegistry region0 successful
but it is not a RASModel, it is a kEpsilon

So, unfortunatley I cannot help you, and I will be glad to read any advice.
And, let me know if you find a solution! :)

Furthermore, as I'm quite beginner with OpenFOAM, I've also some questions about lagrangian particle tracking :
As I said, I'm using the icoUncoupledKinematicParcelFoam solver, which does not solve U. It means that only one-way coupling is possible. However in the Properties file of the particles (in my case kinematicCloudProperties), one can precise if the solution is coupled or not. So, to what is refering this coupling ?

Thanks for the help.

Best,

Mat

Hi Mat,

Thanks for the reply. You got the same error message as I got. If you look into the lagrangian particle tracking libraries (e.g $FOAM_SRC/lagrangian/intermediate/submodels/Kinematic/DispersionModel/DispersionRASModel) you see that it is related to the fact that the dispersion term is calculate on a compressibel ras model indeed, so I suspect this is the root of the problem. Either a new library model should be implemented for incompressible flows, or a trick should be done at solver level to be able to you the compressible model.

Building the new library model is not that trivial because it requires to modefy the library upto the level at which the RAS model is not relevant anymore. I think the whole Kinematic directory and all below it should be duplicated to a directory KinematicIncompressible, which means that all the files below it should be modified to. Or am I wrong? Can I introduce an incompressble ras model at a lower level in the kinematic particle class ? Who can advice?

On you question. You are right that you have an option in the kinematicParticleProperties dictionary, but in the icoUncoupledKinematicParcelFoam solver this is not used because the velocity field is not solved at all. So the option does not work. If you take the pimpleCoupledKinematicParcelFoam solver I sent in the previous post, the U equation is solved and also the particle mometum source is put back in the U equation (line 9 in UEqn.H). So now it is possible to do a two way coupled simulation. This is only relevant if you suspect your particles to have an effect on the flow of course. Such as bubbles in water.

Anyway, hopefully somebody has already developped a transient two way coupled solver for non-reacting particles which includes stochastic dispersion

Regards
Eelco

Solver for particle tracking in incompressible fluid with turbulent dispersion
 

Dear All,

Today I found that turbulent dispersion can be included in particle tracking solvers based on incompressible fluid relatively easy. Just as I anticipated in my previous post, a new build of the particle tracking libary needs to be done in which all compressible turbulence models are replaced by incompressible turbulence models. If you then link your solver to this new library, it is possble to include turbulent dispersion in an incompressible solver during the run. The figure attached below gives an impression of the same cavity case as posted earlier, now with two-way coupled particle tracking in a gravity field with and without turbulent dispersion. As you can see, the particles indeed got dispersed by the stochasticDispersion model. Although I have validated the particle behaviour of individual particles in a stagnant fluid (both the terminal rising speed as well as the particle relaxation time are predicted correctly), I have not yet validated the dispersion models. Nevertheless, qualitatively the dispersion seems to be included correctly.

What I did is copying lagrangian/intermediate library to my own directory, and replace all compressible turbulence refereces with incompressible references, which actually only needs to be done in the files
DispersionRASModel (.C and .H) and BrownianMotionForce (C and H) files below the submodels directory. Also in the Make/options file all the compressible turbulent libraries need to be replaced by incompressible libraries. I have added the new library and the solver based on this library to this post. Here a short instruction to compile it:

1: make your own library directory under your home for the new particle tracking ilbary

mkdir -p $FOAM_USER_SRC/lagrangian/incompressible/
cd $FOAM_USER_SRC/lagrangian/incompressible/

2: copy the excisting Foam particle tracking library to your own directory (only the intermediate directory needs to be recompiled)

cp -r $FOAM_SRC/lagrangian/intermediate/ .
cd intermediate

3: Download the submodels.tar.bz2.gz attachment. Beware: the attachment is not really zipped with gzip, but the .gz extension is just to allow to upload the file to the forum (as gz extension is required); I actually zipped the archive with bzip2 in order to get the zipped archive below the required 100k. So put the attachment into the current directory and do the following

mv submodels.tar.bz2.gz submodels.tar.bz2
tar xjvf submodels.tar.bz2

This should unzip (using bzip2) the submodels directory plus the modified make files.

4: Compile the new library with

wmake libso

Now a new libary libincompressibleLagrangianIntermediate should have been put in your local library directory $FOAM_USER_LIBBIN

The new attachment of the pimpleCoupledKinematicParcelFoam solver can be compiled using the new libary and the dispersion models can this time be included in the kinemeticCloudProperties dict.

Finally, I have included my modified particleTrack utility, kinematicParticleTracksIncompressible, which adds functionality to the one which comes with OF. Next to writing the partlce track it also writes information on other properties such as the residence time and particle diameter. Besides that, some more options can be added in the particleTrackProperties dictonary, such as writing one file per track or writing the data in raw format. Just compile the file and try to run it.

Well, hopefully this is useful. Good luck!

Regards
Eelco

stochatischDispersionLES model added to incompressible pimple solver
 

For those interested: I added stochastic fluctuations in an incompressible LES context as well (the original dispersion model is only implemented for RANS). Just compile the library attached and add

dispersionModel stochasticDispersionLES;

or

dispersionModel gradientDispersionLES;

to the kinematicCloudPropertries dictionary. Have a look at the 3D cavity example.

Note: The same trick has been just to the intermediate.tar.bz2.gz archieve: remove gz extension and unpack it with tar xjvf (using bzip2).

I also update the kinematicParticleTracks utility compiled against the new library. Also some more info per particle track is written (UTubr and tTurb for instance).

One question from my side again: does anybody now how to write the particles independently from the fields ? I would like to write the particle at a higher frequency than I write the field files. Hopefully any knows how to do this. Thanks!

Regards
Eelco

Dear eelcovv,

I'm really interested by your dispersion model for incompressible fluid. Thanks a lot for that.
I've tried to compile it, and I get the following errors on the RAS dispersion model :

I followed your instructions, Do you have an idea about this error ?
I don't understand, because the reference of the incompressible class is also made in BrownianMotionForce.H, for which there is no error.

Hi Mat,

It seems that some reference to the incompressible turbulence libraries can not be found. First of all: I never mentioned that I developped everything in OpenFOAM 2.1.0. So if you are using a diffent version that is likely the cause of the problem. If you are indeed using 2.1, then you should check the Make/options file to see if all the reference to the incompressible RAS turbulence models are correct.

Regards
Eelco

Can't this be done with the partialWrite functionObject? An example of the usage can be found in /opt/openfoam2xx/src/postProcessing/functionObjects/IO/controlDict

Dear eelcovv,

It was indeed a problem of version. I'm looking forward to postprocess the results.
Thanks for this solver :)
I think, I will add also a Saffman lift force for the particles, which is significant in most of the cases.

Best,

Mathieu

That would be an excellent idea! The inclusion of the Saffman lift force and the stochastic dispersion model of Eelco would bring the particle solvers more on par with the dpm solvers in Ansys.

Please share it to the world when it is implemented :)

Lift and dispersion
 

I actually did implement lift, but I am still kind of testing it. It seems to be ok, though.

As for the dispersion model which I modified for incompressible flows: I submitted this to OpenCFD and they have improved the way the dispersion model is implemented in the 2.1.x version: see http://www.openfoam.org/mantisbt/view.php?id=440. Now you can use the original library and still call the dispersion model from my solver.

My lift implementation is still done in 2.1.0 (so on top of my own modified library) and I actually create a whole new particle class for that as I didn't know how to elegantly incorporate that on top of the kinematicParticle class. Probably can also be done much more elegant, as soon as I think everything is ok I will post my solution. In the hope somebody picks it up and improves it in the way how it should be done of course:)

Dear all,

I have a doubt about the pressure gradient force (maybe I should create a new thread?). Usually, this force is evaluated as it is written in the doc file in attachment.

However, in OpenFOAM, the pressure gradient force acting on a particle is computed as following :
value.Su() = mass*p.rhoc()/p.rho()*(p.U() & gradU[p.cell()]);
where p.U() is the velocity of the particle. And this is wrong, because it should be p.Uc(), which is the velocity of the carrier phase. Isn'it ?

Best,

Mathieu

Problems with escaping particles
 

At first: Thanks for this nice and useful tool :)

I'm currently trying to postprocess a simulation with kinematicParticleTracks and it seems that there is a problem with escaping particles at least when kinematicParticleTracks is started at a time at which not all particles are active anymore. In the header files for output (writeRawComposite.H etc.) it tries to write also data for escaped particles and this gives an error (segmentation fault).

Is there any possibility to check if the particle ID is valid before writing out the data for this particle?

I think there is a field called 'active' in the particle properties which is probably used to indicate wether a particle is active or not. You could retrieve this and store it in het list just as I did for the other quantities and used it to test the particle before writing. I don't have time for it right now my self. Good luck!
Eelco

I just found out that partialWrite is rather inefficient: it writes all the files at the interval you specify in the controlDict (which represents the high frequently written fields) and then removes all but the specified fields. It only keeps all the fields at and interval you specify in the partialWriteDict. So you write everything and trow away what you don't need, which is really really inefficient if you for instance only want to keep a small particle tracking file every time step.

Does anybody know a more efficient way to do this ?

I think there has to be more elegant way to do it, but I did it by changing in the solver

runTime.write();

to

kinematicCloud.write();

Now it writes only particle tracking data.

Hello Eelco!

I hope that your incompressible particle tracking library that OpenCFD included in 2.1.X is also included in the official 2.1.1 release! This wil make your solver already a lot cleaner.

Did you already manage to optimize and to incorporate your lift implementation into the kinematicParticleclass? I am very interested in your solver and would like to help you by comparing it to the dpm implementation in Ansys 14.

I am looking forward to an updated version!

New version particle tracking with lift
 

Hi Cedric and others,

Yes, I have lift included to the particle library. I actually created a new class call bubbleCloud as I didn't want to modify the original kinematicParticle class.

Well, I can give my code, as I hope it can be improved. The code contains a library which you should unpack in you own local library location (in my case
/home/eelco/OpenFOAM/eelco-2.1.x/src). Build it with wmake libso

Then there is a solver call pimpleLPTBubbleFoam which uses the library to include the bubble including lift.

Finally, I have added a sample case which calculates the benchmark with bubbels in a square bubble columm performed by Deen et al (2001, Chemical engineering science 56). The bubble tracking information can be found in the constant directory. I have added and extra field called writeEvery in this directory to solve the problem I mentioned earlier in this thread: it allows you to write the bubbles a higher frequency than the data fields.

Note that the files are zipped with bzip2 (as the resulting files are smaller) but that extension is not allowed so I added .gz to the file name. To unpack the files do

mv filename.tar.bz2.gz filename.tar.bz2
tar xjvf filename.tar.bz2

The file must be compiled using OpenFOAM 2.1.x (probably also 2.1.1) because I make use of the bug fixed released on my earlier post.

The pictures show the velocity profiles at three heights in the square bubble column compared to experiments, one without and one with lift.
You can see that the lift improves the results as more spread is established, however, it is not completely correct yet.

One of the thinks I am looking at is the turbulent dispersion, which was on for this case, but did not have an effect. The reason is that I used the SpalartAllmarasDDES turbulence model, which does not internally calculates the k-value. Therefore, the Uturb contribution is 0, where as it should have a value. I also run the case with the oneEqEddy turbulence model and then turbulence fluctuations are included because k is calculated in this turbulence model. My question is: how can I force the k-values to be calculated from the library , even this k is not used by the turbulence model?

Another thing I can not figure out is how to set the value of the force in the new fields Flift, Fdrag and Fbuoy I have created. I was trying to set the fields of the pointer to the particle p in the calcNonCoupled routines, but since it is defined as const it is not allowed to change the values. Removing the const from the type definitation does not work as well, even if done consistently throughout the code. What would be an appropriate way to access those field so that the force can be stored on the particle properties such you can write them to file to analyse later?

Well, so far. Hopefully this is enough information to check the code. Please let me know if you can make any improvements. Hopefully something like this code can be included in the official OF release as it allows to use Lagrangian particle tracking of bubble flows in OF.

Regards
Eelco

And in case you would like to use a utility to create the particle tracks:
here you are. Information can be set in the constant/lptBubbleTrackProperties file.

Regards
Eelco

just a question...
 

Dear Eleco
thanks for your files.
I used them but, when i want to compile kinematicparticletracks, with using "wmake" command, I encounter with forward error:

morsal@morsal-laptop:/opt/openfoam201/applications/solvers/kinematicParticleTracks$ wmake
Making dependency list for source file kinematicParticleTracks.C
SOURCE=kinematicParticleTracks.C ; g++ -m64 -Dlinux64 -DWM_DP -Wall -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -O3 -DNoRepository -ftemplate-depth-100 -I/opt/openfoam211/src/lagrangian/basic/lnInclude -I/opt/openfoam211/src/lagrangian/intermediate/lnInclude -I/opt/openfoam211/src/thermophysicalModels/specie/lnInclude -I/opt/openfoam211/src/thermophysicalModels/basic/lnInclude -I/opt/openfoam211/src/thermophysicalModels/reactionThermo/lnInclude -I/opt/openfoam211/src/thermophysicalModels/radiationModels/lnInclude -I/opt/openfoam211/src/turbulenceModels/incompressible/turbulenceModel -I/opt/openfoam211/src/turbulenceModels/incompressible/RAS/lnInclude -I/opt/openfoam211/src/turbulenceModels/incompressible/LES/lnInclude -I/opt/openfoam211/src/turbulenceModels/incompressible/turbulenceModel -I/opt/openfoam211/src/transportModels -I/opt/openfoam211/src/transportModels/incompressible/singlePhaseTransportModel -I/opt/openfoam211/src/finiteVolume/lnInclude -I/opt/openfoam211/src/meshTools/lnInclude -I/opt/openfoam211/src/regionModels/regionModel/lnInclude -I/opt/openfoam211/src/regionModels/surfaceFilmModels/lnInclude -I/opt/openfoam211/src/sampling/lnInclude -IlnInclude -I. -I/opt/openfoam211/src/OpenFOAM/lnInclude -I/opt/openfoam211/src/OSspecific/POSIX/lnInclude -fPIC -c $SOURCE -o Make/linux64GccDPOpt/kinematicParticleTracks.o
g++ -m64 -Dlinux64 -DWM_DP -Wall -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -O3 -DNoRepository -ftemplate-depth-100 -I/opt/openfoam211/src/lagrangian/basic/lnInclude -I/opt/openfoam211/src/lagrangian/intermediate/lnInclude -I/opt/openfoam211/src/thermophysicalModels/specie/lnInclude -I/opt/openfoam211/src/thermophysicalModels/basic/lnInclude -I/opt/openfoam211/src/thermophysicalModels/reactionThermo/lnInclude -I/opt/openfoam211/src/thermophysicalModels/radiationModels/lnInclude -I/opt/openfoam211/src/turbulenceModels/incompressible/turbulenceModel -I/opt/openfoam211/src/turbulenceModels/incompressible/RAS/lnInclude -I/opt/openfoam211/src/turbulenceModels/incompressible/LES/lnInclude -I/opt/openfoam211/src/turbulenceModels/incompressible/turbulenceModel -I/opt/openfoam211/src/transportModels -I/opt/openfoam211/src/transportModels/incompressible/singlePhaseTransportModel -I/opt/openfoam211/src/finiteVolume/lnInclude -I/opt/openfoam211/src/meshTools/lnInclude -I/opt/openfoam211/src/regionModels/regionModel/lnInclude -I/opt/openfoam211/src/regionModels/surfaceFilmModels/lnInclude -I/opt/openfoam211/src/sampling/lnInclude -IlnInclude -I. -I/opt/openfoam211/src/OpenFOAM/lnInclude -I/opt/openfoam211/src/OSspecific/POSIX/lnInclude -fPIC -Xlinker --add-needed -Xlinker --no-as-needed Make/linux64GccDPOpt/kinematicParticleTracks.o -L/opt/openfoam211/platforms/linux64GccDPOpt/lib \
-llagrangian -L/home/morsal/OpenFOAM/morsal-2.1.1/platforms/linux64GccDPOpt/lib -lincompressibleLagrangianIntermediate -lthermophysicalFunctions -lbasicThermophysicalModels -lspecie -lradiationModels -lincompressibleRASModels -lincompressibleLESModels -lincompressibleTurbulenceModel -lincompressibleTransportModels -lfiniteVolume -lmeshTools -lregionModels -lsurfaceFilmModels -lsampling -lOpenFOAM -ldl -lm -o /home/morsal/OpenFOAM/morsal-2.1.1/platforms/linux64GccDPOpt/bin/kinemeticParticleTracks
/usr/bin/ld: cannot find -lincompressibleLagrangianIntermediate
collect2: ld returned 1 exit status
make: *** [/home/morsal/OpenFOAM/morsal-2.1.1/platforms/linux64GccDPOpt/bin/kinemeticParticleTracks] Error 1
morsal@morsal-laptop:/opt/openfoam201/applications/solvers/kinematicParticleTracks$

I apologize you to help me...
thanks for your help...