Limiting values in exponential formulations

[Arcadia197]

Hey Timofey,

I have noticed also already during earlier testing, that we observe some constraints on the stability of the solution when using the library, with crashes occurring due to Floating Point exceptions in the root-finder of the wall model (Spalding). We have observed this more often now when trying to start the WMLES simulation with coarser time-step to advance the initial transient more quickly. In fact, that's also referred to in your docs under troubleshooting.

My colleague Felix noticed, that in the original implementation of the root-finder in OpenFOAM a limiter for the term $\kappa * u / u_{\tau}$ was set (see current code for nutUSpaldingWallFunction - line 155). However, in your version this limiter is not present (see current code for SpaldingLawOfTheWall – around line 122 and 152). Is this intentionally left out? I suspect, that this might help with stabilizing the simulation when a few cells go haywire, but haven’t done much testing due to limited time resources.
A possible change could be, to change line 122 and 152 from

scalar uPlus = u/uTau;

to

scalar uPlus = min(u/uTau, scalar(100));

The limit of 100 was set arbitrarily and sets a somewhat similar factor to the OF implementation with standard kappa. I've experienced more stable runs using the limited version. If it's working well, this could be also be applied for the other wall-functions as well.

What are your thoughts on this? Do you agree that this could help? Are there any ways I can support testing here? I have some logs showing crashing without and smooth runs with the fix but I'm not sure if they really are of any assistance.

- Julius

[timofeymukha]

Hi Julius!
I have not intentionally removed the guard. I guess I just didn't see it. In principle I am not against this if it helps in practice. The only issue here is whether we are sure that the flow will recover from the errors that caused a u+>100. What I mean is that, for example, the ODE wall models are very stable, but as a downside I have had expensive simulations run to completion only to see that the I have nut equal to 1 million at the wall or whatever, but no crash. I guess in the scenario where people actually pay for compute time this could be a huge problem. So I kind of started appreciating the fact that the Newton-based models are honest and just die :-). But there is of course a balance here. Not sure if this can be "rigorously" tested, as this is in the end an engineering trick. I suggest you continue seeing how this works out for you and report back after some time.

[timofeymukha]

I've added limiters for the solution on develop. But yes, Spalding's law is particularly bad in this respect due to exponential terms which depend on u+. Additional stability is added by using derivative-free equation solver instead of Newton.