|
[Sponsors] |
SimpleFoam inconsistency in results over time |
|
LinkBack | Thread Tools | Search this Thread | Display Modes |
January 27, 2022, 11:34 |
SimpleFoam inconsistency in results over time
|
#1 |
Member
Chris
Join Date: Dec 2020
Posts: 45
Rep Power: 5 |
Hey, I would like to implement a inner loop for simpleFoam but if I do, the results differ from each other:
I added just a while loop to compute UEqn and pEqn: Code:
while (simple.loop()) { Info<< "Time = " << runTime.timeName() << nl << endl; int i = 0; // new line while(i<5) // new line { // Pressure-velocity SIMPLE corrector { #include "UEqn.H" #include "pEqn.H" } laminarTransport.correct(); turbulence->correct(); i++; // new line } runTime.write(); runTime.printExecutionTime(Info); } Another question: How can I print/set relaxationFactors (which are defined in fvSolution) inside the code to change it in my loop? I know that relaxation is applied by e.g. UEqn.relax() but how can I reset it? Thanks in advance for any help. Pyro Last edited by Pyrokrates; February 9, 2022 at 07:22. |
|
January 29, 2022, 07:27 |
|
#2 |
Senior Member
Join Date: Apr 2020
Location: UK
Posts: 740
Rep Power: 14 |
Chris
remember that the SIMPLE algorithm is an iterative algorithm - it does not get to a converged solution in one iteration. Indeed, it can't do this if you have under-relaxation set. What you are seeing in your successive inner loops is just convergence of the solution. I guess you are trying to implement an alternative to the PIMPLE scheme? Have you looked at the PIMPLE scheme? As for changing the RFs during a run - these are read from the fvSolution dictionary .. you can change these manually "on the fly" during a run, by editing the fvSolution file. I dare say that you could probably code up something using a coded function object to update the contents of the dictionary file ... I don't really recommend it though. My approach is to control the RFs from outside the run. Again, you could write a script to monitor the residuals and adjust the RFs ... you may spend more time getting the logic right, though, than just doing it manually! Good luck! |
|
January 31, 2022, 07:04 |
|
#3 |
Member
Chris
Join Date: Dec 2020
Posts: 45
Rep Power: 5 |
Hey,
thank you for the answer. What I´m trying to do is changing the relaxation on the fly... I´m using porous media where I get oscillations inside my velocity right next to the porous media. If I lower the relaxation, I can reduce this oscillations... Now I have the problem, that my velocity changes locally over time so I need "high relaxation" for letting U change and "low relaxation" for the correct caluclation of my U field. I know that SIMPLE loop is interative... thats why I would like to make one iterative loop to lets say more iterative loops where the relaxation changes after convergence. If I do it manually e.g. for my initial solution I run normal simpleFoam with Urelax = 0.1... after 500 time steps I change my Urelax to 0.001 and do 1000 more steps ... until I get the correct Ufield... That is exactly what I would like to do for my time simulation (I don´t want to simulate 3600 seconds with dt=0.1 by manually change of relaxation for each step). Thanks in advance Pyro |
|
February 1, 2022, 18:24 |
|
#4 |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
In my opinion your for-loop does not make sense at all. The only difference is, that you will get a faster steady-state solution during the calculation == speedup factor of 6 in your case, as you do 6 more iterations (but its not a speed-up in terms of run-time; simply you do 6 times less simple-loops; the simple-loop is simply a for-loop).
The simple-loop itself, is just an iterative algorithm around UEqn.H and pEqn.H. So if you do 6 simple loops, you should get the same result compared to 1 simple loop and your additional 6 loops. Both do make 6 loops around UEqn.H and pEqn.H. As we don´t have any time-information here, its simply as follows:
With the new p-U fields, you start again on the first point and construct the momentum matrix. Here, and this is the main point, it does not matter if you go one SIMPLE iteration or one of your iterations. Hence, each iteration in your for-loop will differ except if you reached the steady state solution. Your for-loop will not gain anything in your case. Furthermore, the relaxation factors can be set by using the <Function1> Class. So you have more freedom to set a time-dependent or what ever-dependent relaxation factor. Despite of that you can use Doxygen to search for the function you are searching: https://www.openfoam.com/documentati...1fvMatrix.html https://www.openfoam.com/documentati...ebf646011c212e As you can see, the relax() function is defined in the fvMatrix class. We also have the <fvMatrix>.relax(const scalar) function. You simply can add any number here which is used as relaxation factor rather than the one provided and stored in the class (read from fvSolutions).
__________________
Keep foaming, Tobias Holzmann Last edited by Tobi; February 2, 2022 at 11:47. |
|
February 2, 2022, 04:32 |
|
#5 |
Senior Member
Join Date: Apr 2020
Location: UK
Posts: 740
Rep Power: 14 |
I have to confess that I agree with Tobi - I am confused by your discussion of "time" when you are using simpleFoam ... it's a steady state solver. Perhaps you need to tell us a little more about the problem you are trying to solve.
|
|
February 2, 2022, 09:47 |
|
#6 | |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
Private message from Chris
Quote:
If so, simply use simpleFoam and converge your solution (only p/U). Use the U or phi field of the converged solution for a second run using e.g., scalarTransportFoam or use simpleFoam -postProcess xyFO to calculate the concentration based on a passive transport equation on the fly (FO); here you could also add source terms if necessary. If you are interested in time-based data, you need to change the ddt scheme for the passive scalar. If you want to have an in-time velocity field and an in-time concentration field, you cannot use simpleFoam (maybe SIMPLEC + almost no relaxation + time-derivative activated; but a quick look into the source code of simpleFoam gives us the hint, that there is no time derivative implemented, so its not working and not worth to investigate). Therefore, you need to use pimpleFoam. Here, you can use the PIMPLE algorthm to converge between each time step the p-U correlation much more (but this is only needed if the Courent number is way much larger than 1; see my book for the usage of pimple algorithm for example)... If you want to do anything else, you need to clarify your statements and what you want to do. Out of the box, I have no other idea what you would like to do.
__________________
Keep foaming, Tobias Holzmann |
||
February 9, 2022, 12:39 |
Update on inner loop results
|
#7 | |
Member
Chris
Join Date: Dec 2020
Posts: 45
Rep Power: 5 |
Hey,
thank you for the interesting and helpful input. Quote:
Normally I would totally agree with you saying that e.g 5 simple loops are the same as 1 simple loop with an inner loop doing 5 times the same but the results differ. I have uploaded two simulations where I compare my results. What I`m doing is simulating the velocity field over a porous media (left hand side of both images). On the right hand side you can see the results (turquoise line represent the normal simpleLoop and the brown line the inner looping method). I have done 2 simulations https://we.tl/t-N3pm5wCFzd: 1) init5000in1.png for initial flow field with 3 steps: 500 time steps with Urelax of 0.1, than 2000 time steps with Urelax of 0.001, than 2500 time steps with Urelax of 0.00001 For normal simple loop I´m doing this manually by changing Urealx after the specific number of time steps until I get the final result after 5000 steps. For the modified simple loop with inner for (or while) loop I have automated this with two arrays: One with the time steps in it rUTS(500,2000,2500) and one with the relaxation factors rU(0.1,0.001,0.00001). This way I can build two loops around UEqn and pEqn to get the 5000 steps in one single simple loop. As you can see, the results differ. 2) With the manually generated initial solution of 1) I start calculating my "final solution". Therefore I only have one time step number (5 here) and one relaxation factor (0.001 here). As I expect for both results, with the increase of relaxation from 0.00001 (final init relaxation) to 0.001 (relaxation for final solution), the velocity shows more oscillations (which can be seen here). But you can also see the same behavior like in 1)... 5 simple loops are not the same as 1 simple loop with 5 inner for loops. I´m not sure whats the reason for it. Maybe the algorithm uses some "backup" or "presaved" data as another inital guess or something like this. That would explain the difference for me because for 1 loop with 5 inner iterations I would not update this "guess" or something like this. For the point of setting relaxation... I did this all (thanks for your help) with UEqn.realx(rU[rUi]). At the moment I get the rU[rUi] by a new defined dictionary in fvSolution and set it over the inner loops afterwrds. This works quite good for what I would like to do but is there also a way to get the normal defined relaxation value of fvOptions and overwrite it? Than I could leave my UEqn.realx() as default with the advantage, that I can also leave the standard relaxation factor unchanged and change it only if needed. Here is some code to clarify the inner loop part and the way, I read my fvSolution: simpleFoam: Code:
while (simple.loop()) { Info<< "Time = " << runTime.timeName() << nl << endl; forAll(rU,rUi) { for(int i = 0; i < rUTS[rUi]; i++) { // Pressure-velocity SIMPLE corrector { #include "UEqn.H" #include "pEqn.H" } laminarTransport.correct(); turbulence->correct(); } } runTime.write(); runTime.printExecutionTime(Info); } Code:
IOdictionary fvSolution // Reading different relaxation factors ( IOobject ( "fvSolution", runTime.system(), runTime, IOobject::MUST_READ ) ); dictionary& fvSol = fvSolution.subDict("relaxationFactors").subDict("equations"); label NRelaxU(readLabel(fvSol.lookup("NRelaxU"))); List<scalar> rU(NRelaxU); // List of relaxation parameter List<label> rUTS(NRelaxU); // List of relaxation time steps forAll(rU,rUi) { rU[rUi] = fvSol.subDict("relaxU").lookupOrDefault<scalar>(("rU" + name(rUi)),0.1); rUTS[rUi] = fvSol.subDict("relaxU").lookupOrDefault<label>(("rUTS" + name(rUi)),1); } Code:
relaxationFactors { fields { p 0.1; } equations { U 0.1; //(0.1, 0.001, 0.000001) for init flow field k 0.1; epsilon 0.1; NRelaxU 3; relaxU { rU0 0.1; rU1 0.001; rU2 0.000001; rUTS0 500; rUTS1 2000; rUTS2 2500; } } } Chris |
||
February 10, 2022, 14:15 |
|
#8 |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
No idea why you does make it so complicated:
Code:
relaxationFactors { equations { U { type table; values ( (0 0.1) (10 0.2) (50 0.3) ); } And, I don't expect that there are any differences between the normal SIMPLE loop and the SIMPLE loop with additional for-loop. Its the same. If you need a proof, I will give it to you.
__________________
Keep foaming, Tobias Holzmann |
|
February 11, 2022, 05:30 |
|
#9 |
Member
Chris
Join Date: Dec 2020
Posts: 45
Rep Power: 5 |
Hey Tobi,
thanks you for the answer. Using this kind of linear blending result in following error: Code:
Attempt to return dictionary entry as a primitive Problem is that I can only define the changing relaxation once per run... What I need is the change of relaxtion in each time step (I know for simpleFoam there is no real time, its more about iterations). But for the part I have programmed for my concentration I need time steps. This is why I would like to change the relaxation every calculation step (increase Urealx to let U change local in my mesh but reduce it in the same step again to get correct pressure drop and lower fluctuations inside my velocity for the concentration calculation). A proof for the same results of simpleFoam with and without for loop would be great because I´m confused about my results. I thought the results should be the same but as you can see in my uploaded images, it differs. Maybe its the porousMedia... I will do some testing without porous media and let you know if there is a change... Thanks for all the time and help until now Chris |
|
February 11, 2022, 07:05 |
|
#10 |
Member
Chris
Join Date: Dec 2020
Posts: 45
Rep Power: 5 |
Hey,
I have made some comparisons and need to correct/clarify something I wrote not correctly in my previous posts... The first thing is, for my calculations I`m using porousSimpleFoam which I have modified with my concentration models, not the simpleFoam... I thought the two solvers should not differ much only for the part of using porous media but now it looks like they do. I have compared four simulations doing 100 timesteps with Urelax 0.1 in a curved pipe: (Link to the plots: https://wetransfer.com/downloads/113...1113935/8e36d6) 1) simpleFoam [red line] 2) simpleFoam (with implemented inner loop, 10 time steps with each 10 inner loops) [pink line] 3) porousSimpleFoam (without porous media) [blue line] 4) porousSimpleFoam (without porous media, with implemented inner loop, 10 time steps with each 10 inner loops) [pink line] The results are following: 1) and 2) differ from each other. 1) differ from 3) ,which confuses me a little bit. 3) and 4) also differ. Maybe this helps a little to solve the problem. Thanks Chris Last edited by Pyrokrates; February 11, 2022 at 09:31. |
|
February 11, 2022, 14:55 |
|
#11 |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
So here is the proof that both SIMPLE and SIMPLE + Inner-Loop are identical:
I simply made a copy of simpleFoam and changed the source file accordingly: Code:
Info<< "\nStarting time loop\n" << endl; label iter = 0; while (simple.loop()) { Info<< "Time = " << runTime.timeName() << nl << endl; for (int i = 0; i<4; ++i) { // --- Pressure-velocity SIMPLE corrector { #include "UEqn.H" #include "pEqn.H" } laminarTransport.correct(); turbulence->correct(); ++iter; } runTime.write(); runTime.printExecutionTime(Info); if (iter == 100) { Info<< "We made 100 Iterations" << endl; runTime.writeAndEnd(); } } Info<< "End\n" << endl; return 0; SimpleFoam and porousSimpleFoam might give you slightly different results. I never checked the seceond one but there is definitely some other implementation which you have to take into consideration. You analysis seems to modify the relaxation factors and hence, it is obvious that you don't get the same guys. For me it is still unclear what you want to achieve as it does not make any sense to me. Nevertheless, as I showed and Tobermony also backed,
Cheers.
__________________
Keep foaming, Tobias Holzmann |
|
February 14, 2022, 06:19 |
|
#12 |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
Hey all,
I got a private message from Chris in which he clarifies his ideas and what he want to do. I will summarize and give feedback now:
The last question Chris mentioned:
__________________
Keep foaming, Tobias Holzmann |
|
February 14, 2022, 06:39 |
|
#13 |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
For comparison reasons and transparency:
https://holzmann-cfd.com/forum/simpl...eration.tar.gz The solver is not attached as you simply can built it by copy simpleFoam and change the source file accordingly (as given two posts above).
__________________
Keep foaming, Tobias Holzmann |
|
February 14, 2022, 07:35 |
|
#14 | |
Senior Member
Join Date: Apr 2020
Location: UK
Posts: 740
Rep Power: 14 |
Quote:
For your case, if you are calculating a filter blockage/deposition rate in your simulation, and if the rate of deposition is slow compared to the flow timescales, then the above pseudo-steady approach for the flow field will possibly give you a reasonable estimate. Otherwise, you need a fully coupled transient as per Tobi's suggestion. PS - I forgot to add - Tobi is absolutely right - a part converged SIMPLE solution could contain absolute junk ... and the rate of convergence of the field is totally case dependent, so there's no real way of taking advantage of a SIMPLE solution after only a handful of iterations; you probably need 100's, in which case just go to full convergence, and you are heading towards the 1st approach above. |
||
February 14, 2022, 08:36 |
|
#15 |
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 52 |
Totally agree Tobermory, thank you for bringing up another field in which we use a similar approach.
However, I want to update this talk in terms of differences, of the given implementation I presented in prior posts, if we use under-relaxation for p field. So lets reconsider the code snippet I presented above: Code:
while (simple.loop()) { Info<< "Time = " << runTime.timeName() << nl << endl; for (int i = 0; i<4; ++i) { // --- Pressure-velocity SIMPLE corrector { #include "UEqn.H" #include "pEqn.H" } laminarTransport.correct(); turbulence->correct(); ++iter; } runTime.write(); runTime.printExecutionTime(Info); if (iter == 100) { Info<< "We made 100 Iterations" << endl; runTime.writeAndEnd(); } } However, Chris mentioned, that using under-relaxation for the p field, we obtain different results (at least for experienced user it is obvious why; hence, I want to clarify it). The field relaxation is simply: Therefore, we update the new iteration (n+1) with the old iteration and the new value. Having no relaxation means , and hence, the new evaluated value is used completely. Despite of that, any value of less than 1, takes the old p field into account. And now the differences:
Code:
while (simple.loop()) { Info<< "Time = " << runTime.timeName() << nl << endl; for (int i = 0; i<4; ++i) { // --- Pressure-velocity SIMPLE corrector { #include "UEqn.H" // We need to store the previous iteration (overwriting the p-old) p.storePrevIter(); #include "pEqn.H" } laminarTransport.correct(); turbulence->correct(); ++iter; } runTime.write(); runTime.printExecutionTime(Info); if (iter == 100) { Info<< "We made 100 Iterations" << endl; runTime.writeAndEnd(); } }
__________________
Keep foaming, Tobias Holzmann |
|
February 14, 2022, 08:59 |
|
#16 |
Member
Chris
Join Date: Dec 2020
Posts: 45
Rep Power: 5 |
Hey Tobi,
thank you very much for the help. The part with the storing of p was exactly what I supposed when I said, that it looks like inner looping does not use some pre saved data correctly. Thanks and have a nice day Chris |
|
Tags |
innerloop, relaxation factor, simplefoam, timestep |
Thread Tools | Search this Thread |
Display Modes | |
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
[General] Extracting ParaView Data into Python Arrays | Jeffzda | ParaView | 30 | November 6, 2023 22:00 |
laplacianFoam with source term | Herwig | OpenFOAM Running, Solving & CFD | 17 | November 19, 2019 14:47 |
Setting up Lid driven Cavity Benchmark with 1M cells for multiple cores | puneet336 | OpenFOAM Running, Solving & CFD | 11 | April 7, 2019 01:58 |
Inconsistencies in reading .dat file during run time in new injection model | Scram_1 | OpenFOAM | 0 | March 23, 2018 23:29 |
High Courant Number @ icoFoam | Artex85 | OpenFOAM Running, Solving & CFD | 11 | February 16, 2017 14:40 |