View Issue Details

IDProjectCategoryView StatusLast Update
0002292OpenFOAM[All Projects] Bugpublic2016-10-20 09:11
ReporterqiqiAssigned Tohenry 
Status resolvedResolutionfixed 
PlatformGNU/LinuxOSUbuntuOS Version15.04
Product Version 
Fixed in Versiondev 
Summary0002292: SpalartAllmaras in pisoFoam not restarting correctly
Descriptionturbulence->validate(); in pisoFoam.C changes nut, thereby making the flow different from the case without save-and-restart.
Steps To ReproduceRun a simulation for 10 steps, save it, and continue for another 10 steps. Compare the outcome of the second 10 steps with the outcome of running 20 steps from the beginning. The difference will be in the 5th digit. But for chaotic turbulent flows, this difference will amplify over time, making simulations with and without restart totally different.
TagsNo tags attached.



2016-10-14 14:42

manager   ~0007008

Chaotic turbulent flows are indeed chaotic and only the averaged statistics are reproducible. If you average the results for a long time, restart and restart the averaging do the averaged statistics agree between the two runs?


2016-10-14 14:49

reporter   ~0007009

That's correct. I found this bug, which does not occur for other LES models, when experimenting with numerical methods for sensitivity analysis. Being able to restart and continue on exactly the same trajectory is important for me, and likely for other numerical methods, e.g., for computing Lyapunov exponents of the flow dynamics.


2016-10-14 15:00

manager   ~0007010

All turbulence->validate() does is call correctNut which updates the nut field in the same manner as it does after the calculation of nuTilda. This call is necessary to ensure the nut field is physical and consistent with nuTilda at the beginning of the run, particularly when starting for user-input fields in which nut may be specified as 0. There is no difference in this procedure between SpalartAllmaras and any of the other models.


2016-10-15 02:23

reporter   ~0007011

You are right. turbulence->validate() should not have changed the nut that was computed in the last time step. But somehow it does when using SpalartAllmarasDES.


2016-10-15 07:29

manager   ~0007012

Here is the correctNut code from the SpalartAllmarasDES model:

template<class BasicTurbulenceModel>
void SpalartAllmarasDES<BasicTurbulenceModel>::correctNut
    const volScalarField& fv1
    this->nut_ = nuTilda_*fv1;


template<class BasicTurbulenceModel>
void SpalartAllmarasDES<BasicTurbulenceModel>::correctNut()

I don't see anything wrong with this.


2016-10-15 07:44

manager   ~0007013

Incidentally SpalartAllmarasDDES and SpalartAllmarasIDDES are derived from SpalartAllmarasDES and inherit the same correctNut code.


2016-10-18 13:47

updater   ~0007027

I suspect that the original diagnosis isn't fully correct. The re-calculation of the "nut" field is very likely only a symptom of the real problem.
Furthermore, the original report refers to "SpalartAllmaras" and not "SpalartAllmarasDES", which increases my suspicion.

@qiqi: Can you please provide a test case that replicates the problem that you're seeing?
I ask this because the most likely suspect is actually the boundary conditions at the inlet(s), given that if it's using a random inlet source for flow speed and/or pressure values, that would justify the problem.

If this is the case, either the inlet values need to be pre-planned for each time step or a re-programmable random number generator would have to be implemented. The random generators in C++11 do have this feature, namely to get a seed for at a particular time step, to be re-use in the next run.


2016-10-18 15:11

reporter   ~0007028

The testcase is encapsulated in

The case is in

Note that all the boundary conditions are "clean" -- fixedValue, zeroGradient.

The test passes only with a modified version of pisoFoam, the only difference between which and the shipped version is the validate(); line.


2016-10-19 00:27

updater   ~0007031

Many thanks for the test case! With it I've been able to reproduce the same problem, but I haven't fully figured out yet the real source of the problem. The narrowest I've gotten when trailing back was that this expression seems to be at the genesis of the problem:

   this->nut_ = nuTilda_*fv1;

I was able to do a quick test after running pisoFoam once for 1s and then using the following hack in pisoFoam:

      return 0;

The "nut" field was considerably different, although still in similar ranges of values.

OK, I believe I've got it. I see what's wrong. In "SpalartAllmarasDES<BasicTurbulenceModel>::correct()", "fv1" is calculated with the "nuTilda" from before "nuTildaEqn" was solved. "correctNut(fv1)" at the end of the method uses the outdated fv1 value, resulting in the discrepancy. Removing "fv1" from the call should re-establish consistency.

Attached is the file "SpalartAllmarasDES.C" that replaces "src/TurbulenceModels/turbulenceModels/LES/SpalartAllmarasDES/SpalartAllmarasDES.C".

@qiqi: I haven't enough time today to test this, so any chance you can test this fix? It should work as intended, but it's better to make sure ;)


2016-10-19 00:28


SpalartAllmarasDES.C (10,873 bytes)


2016-10-19 00:28


SpalartAllmarasDES-2.C (10,873 bytes)


2016-10-19 00:30

updater   ~0007032

OK... I'm not use to yet to the possibility to attach files along with the comment... please ignore the file "SpalartAllmarasDES-2.C", which is a duplicate of the previous upload "SpalartAllmarasDES.C".


2016-10-19 10:47

updater   ~0007036

@qiqi: I forgot to mention that after replacing the file, please run the "Allwmake" script in the folder "src/TurbulenceModels", because it needs to update at least a couple of libraries that depend on this template class.


2016-10-19 23:20

updater   ~0007049

@Henry: I've confirmed that the provided test case works as intended when using the attached file "SpalartAllmarasDES.C", for replacing "src/TurbulenceModels/turbulenceModels/LES/SpalartAllmarasDES/SpalartAllmarasDES.C"

The change was simply this:

  - correctNut(fv1);
  + correctNut();

so that the final calculated "nut" field within the method "correct()" is based on the latest "nuTilda" and not based on the one from before solving the equation.

I've also done a quick look through the other sibling template classes within LES and the only ones that use "correctNut()" with additional arguments are "dynamicKEqn" and "dynamicLagrangian", but those two rely on calculations done based on "U", which doesn't change within the "correct()" method.


2016-10-20 02:10

reporter   ~0007051

Thanks for fixing this issue. Great job Henry.


2016-10-20 02:15

reporter   ~0007052

I should be thanking @wyldckat. Thanks wyldckat. Great job.


2016-10-20 09:11

manager   ~0007056

Resolved in OpenFOAM-dev by commit 667dc10af61341690286ad2420b996937f94464f

Issue History

Date Modified Username Field Change
2016-10-14 14:31 qiqi New Issue
2016-10-14 14:42 henry Note Added: 0007008
2016-10-14 14:49 qiqi Note Added: 0007009
2016-10-14 15:00 henry Note Added: 0007010
2016-10-15 02:23 qiqi Note Added: 0007011
2016-10-15 07:29 henry Note Added: 0007012
2016-10-15 07:44 henry Note Added: 0007013
2016-10-18 13:47 wyldckat Note Added: 0007027
2016-10-18 15:11 qiqi Note Added: 0007028
2016-10-19 00:27 wyldckat Note Added: 0007031
2016-10-19 00:28 wyldckat File Added: SpalartAllmarasDES.C
2016-10-19 00:28 wyldckat File Added: SpalartAllmarasDES-2.C
2016-10-19 00:30 wyldckat Note Added: 0007032
2016-10-19 10:47 wyldckat Note Added: 0007036
2016-10-19 23:13 wyldckat Assigned To => henry
2016-10-19 23:13 wyldckat Status new => assigned
2016-10-19 23:20 wyldckat Note Added: 0007049
2016-10-20 02:10 qiqi Note Added: 0007051
2016-10-20 02:15 qiqi Note Added: 0007052
2016-10-20 09:11 henry Status assigned => resolved
2016-10-20 09:11 henry Resolution open => fixed
2016-10-20 09:11 henry Fixed in Version => dev
2016-10-20 09:11 henry Note Added: 0007056