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

Documentation of trackToFace subroutine

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

Like Tree1Likes
  • 1 Post By boeleman

Reply
 
LinkBack Thread Tools Display Modes
Old   September 12, 2013, 19:56
Default Documentation of trackToFace subroutine
  #1
New Member
 
Join Date: May 2012
Posts: 17
Rep Power: 4
boeleman is on a distinguished road
Hello,

I am trying to write a Volume of Fluid (VoF) + Lagragian Particle Tracking (LPT) solver based on OF2.1.1 interDymFoam and the basicKinematicCloud template. However, my code hangs in the trackToFace subroutine in particleTemplates.H. Using gdb I figured out that tetLambda gives back a very small value for lam, lambdaMin is set to 0 and the correction algorithm kicks in, but my code still hangs. I would like to understand the the trackToFace and tetLambda code better, but there is not a lot of documentation online. Could anyone provide references to the algorithms used?

Thanks,

Arnout
boeleman is offline   Reply With Quote

Old   September 23, 2013, 19:05
Default
  #2
New Member
 
Join Date: May 2012
Posts: 17
Rep Power: 4
boeleman is on a distinguished road
Figured I'd post what I think is going on based on reading the code. Not sure this is correct:

1) A distance "dCorr" is picked based on either the parcel velocity or the maximum Courant number, whichever one is smaller.
2) findtris checks whether while going from the center of the tetrahedral, in which the parcel is located, to its current position + dCorr a face is crossed.
3) hitWallFaces checks whether a wall was hit.
4) if both findtris and hitWallFaces do not find any faces the trackToFace will return 1.
otherwise tetLambda will calculate the actual amount traveled before the face was hit.

What I am still confused about is why first an initial guess is calculated using findTris, and why not finding any faces with findTris guarantees that no faces where hit. It calculated from the center of a tet after all and not using the actual particle trajectory.
boeleman is offline   Reply With Quote

Old   September 28, 2013, 12:42
Default
  #3
New Member
 
Join Date: Apr 2013
Posts: 24
Rep Power: 3
GPesch is on a distinguished road
Hi Arnout,

sorry for just throwing that at you (as I didnt have time to look at it closely):

I am working with the DSMC Solver and sometimes also have the problem that one of my parcels is stuck in the trackToFace Routine. I could not figure out what exactly happenend and also didnt fully understand the code of trackToFace (as it looks rather complicated )...

Someone was giving me the hint of looking into

Graham Macpherson, Niklas Nordin and Henry G. Weller: "Particle tracking in unstructured, arbitrary polyhedral meshes for use in CFD and molecular dynamics" which was released in
Commun. Numer. Meth. Engng 2009, 55: 263-273.

I didnt look very closely into it but the paper seems to actually describe the trackToFace routine.

Please contact me if you are interested and if you need help gaining access to the publication.

Georg.
GPesch is offline   Reply With Quote

Old   September 28, 2013, 17:35
Default
  #4
New Member
 
Join Date: May 2012
Posts: 17
Rep Power: 4
boeleman is on a distinguished road
Hello Georg,

Thank you so much for the reference. I just downloaded it and will take a look at it. I did find a workaround for the hanging particles, which might be useful to you. In the trackToFace routine you twice find the line:

position_ += trackingCorrectionTol*(tet.centre() - position_);

in the particle position correction algorithm. trackingCorrectionTol is in the order of 10^-5. If you increase this value to, for example, 1 findTris and tetLambda will agree which other and there will be no more hanging particles (at least not for me). It is not a very elegant hack, and I am going to look at the paper to find a more elegant solution, but it does get the job done. Hopefully it will work for you too.

Cheers,

Arnout
boeleman is offline   Reply With Quote

Old   October 1, 2013, 17:31
Default
  #5
New Member
 
Join Date: May 2012
Posts: 17
Rep Power: 4
boeleman is on a distinguished road
Just read the article. Things make a lot more sense now. (The reason for having trackToFace is because it is computationally cheaper to store which cell a particle is in, than to find the cell it is in by calculating it every iteration. How lambda is part of the calculation which face was crossed). However, I feel there must be a follow up article, because this article does not mention tet decomposition as a workaround to particles hanging in concave cells. Would you happen to know which article this can be found in?
boeleman is offline   Reply With Quote

Old   November 6, 2013, 04:13
Default
  #6
New Member
 
Join Date: Apr 2013
Posts: 24
Rep Power: 3
GPesch is on a distinguished road
Hi Arnout,

sorry for not getting back to you earlier!

Unfortunately I do not know whether a follow-up article is existing! If you stumble over it at some point I would be very happy if you could let me know!

I just wanted to thank you for the little "hack" you posted: I haven't tried it yet, but I will definitely do in order to avoid the hanging particles! Thanks!

Georg
GPesch is offline   Reply With Quote

Old   November 13, 2013, 20:35
Default
  #7
New Member
 
Join Date: May 2012
Posts: 17
Rep Power: 4
boeleman is on a distinguished road
Hello Georg,

Figured I might as well share the modifications I made:

https://www.dropbox.com/s/9jc27rwjp9...ngianOf211.zip

I ended up writing an extra counter, which resets particles to the center of a grid cell whenever the number of iterations is larger than 100. As long as the counter is smaller than 100, the normal correction is used. Also, I modified the calculation at the wall. My particles were getting stuck because they were larger than grid cells at the wall. Last but not least I think I found a bug in particleTemplates.C where I replaced triI with tI in the trackToFace subroutine. You can do a diff to see the details. My code runs without hanging now, so hopefully it helps you too. The article was very helpful!

Cheers,

Arnout
GPesch likes this.
boeleman is offline   Reply With Quote

Old   February 27, 2014, 06:20
Default
  #8
New Member
 
Join Date: Apr 2013
Posts: 24
Rep Power: 3
GPesch is on a distinguished road
Finally had the time to test it out! I'm glad, that you made the modifications as I was to lazy to look into the code! Works well for me—you are my hero, thanks a bunch!!
GPesch is offline   Reply With Quote

Old   May 19, 2014, 09:40
Default
  #9
New Member
 
Join Date: Feb 2014
Posts: 2
Rep Power: 0
cbpf is on a distinguished road
Hello boeleman,

Can you say how can i use your files to overcome the lagrangian reconstructPar errors?

I'm very new in openfoam and i don't know how to recompile only that folders.

Thank you.
cbpf is offline   Reply With Quote

Old   May 19, 2014, 09:55
Default
  #10
New Member
 
Join Date: Apr 2013
Posts: 24
Rep Power: 3
GPesch is on a distinguished road
Hi,

what kind of reconstructPar errors are you getting?
GPesch is offline   Reply With Quote

Old   May 19, 2014, 10:45
Default
  #11
New Member
 
Join Date: Feb 2014
Posts: 2
Rep Power: 0
cbpf is on a distinguished road
Thank you for the fast answer.

The error is:

--> FOAM Warning :
From function void Foam:article::initCellFacePt()
in file particle/particleI.H at line 759
cell, tetFace and tetPt search failure at position (-0.0125212 -0.00987232 0.00017957)
for requested cell 0

--> FOAM FATAL ERROR:


FOAM aborting

#0 Foam::error:rintStack(Foam::Ostream&) in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#1 Foam::error::abort() in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#2 Foam:article::initCellFacePt() in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/lib/liblagrangian.so"
#3 Foam::Cloud<Foam:assiveParticle>::initCloud(bool ) in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/lib/liblagrangian.so"
#4 Foam::reconstructLagrangianPositions(Foam:olyMes h const&, Foam::word const&, Foam::PtrList<Foam::fvMesh> const&, Foam::PtrList<Foam::IOList<int> > const&, Foam::PtrList<Foam::IOList<int> > const&) in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/lib/libreconstruct.so"
#5
in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/bin/reconstructPar"
#6 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#7
in "/home/celio/OpenFOAM/OpenFOAM-2.3.0/platforms/linux64GccDPOpt/bin/reconstructPar"
Abortado (imagem do n˙cleo gravada)
cbpf is offline   Reply With Quote

Old   May 28, 2014, 11:35
Default
  #12
New Member
 
Craig White
Join Date: Nov 2013
Posts: 5
Rep Power: 3
DSMC123 is on a distinguished road
Quote:
Originally Posted by boeleman View Post
Hello Georg,

Figured I might as well share the modifications I made:

https://www.dropbox.com/s/9jc27rwjp9...ngianOf211.zip

I ended up writing an extra counter, which resets particles to the center of a grid cell whenever the number of iterations is larger than 100. As long as the counter is smaller than 100, the normal correction is used. Also, I modified the calculation at the wall. My particles were getting stuck because they were larger than grid cells at the wall. Last but not least I think I found a bug in particleTemplates.C where I replaced triI with tI in the trackToFace subroutine. You can do a diff to see the details. My code runs without hanging now, so hopefully it helps you too. The article was very helpful!

Cheers,

Arnout
Hi Arnout,

I have been seeing similar errors with hanging during the trackToFace function and would like to try your code to sort it, but the .zip file doesn't seem to be working? I'd be grateful to see what you have done.

For interest, I was having similar problems a few years ago with small (microscale) meshes and Graham Macpherson sent me a small loop to add to particleTemplates.C that fixed that issue, but does not solve this one. The loop is:

Code:
tetAreas[0] = tet.Sa();
        tetAreas[1] = tet.Sb();
        tetAreas[2] = tet.Sc();
        tetAreas[3] = tet.Sd();
        
            //******
            for (label i = 0; i < 4; i++)
            {
                tetAreas[i] /= (mag(tetAreas[i]) + VSMALL);
            }
            //******

        FixedList<label, 4> tetPlaneBasePtIs;

        tetPlaneBasePtIs[0] = basePtI;
        tetPlaneBasePtIs[1] = f[fPtAI];
        tetPlaneBasePtIs[2] = basePtI;
        tetPlaneBasePtIs[3] = basePtI;
Where the bold is the new code and the rest will show you the correct place in the code to place it. I suspect this stopped particles being larger than the grid cells at the wall as you discussed previously. This loop never seems to have made it into any OpenFOAM release, but Graham was still working for OpenCFD when he sent me this.

Thanks,

Craig

Edit: Downloaded the .zip file now. I hope my code snippet above is at least of interest, if not use.
DSMC123 is offline   Reply With Quote

Reply

Tags
kinematic cloud, lpt, tetlambda, tracktoface, vof

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
CoCoons Project - Community-driven Documentation on OpenFOAM« Technology holger_marschall OpenFOAM Announcements from Other Sources 3 December 13, 2013 10:14
Errorin subroutine appeared when applying cavitation model pitisrisuk CFX 1 July 2, 2012 04:36
Doxygen documentation Tanay OpenFOAM Installation 9 September 23, 2011 12:40


All times are GMT -4. The time now is 05:50.