Stability issue in dynamicLagrangian model

Hi,

There is a stability issue in the dynamicLagrangian model, which is caused by flm_/fmm_. Some cells may have 0 value in fmm_. If replacing flm_/fmm_ with flm_/(fmm_+fmm0_), the model will be stable, where fmm0_ is a VSMALL defined in the model.

Attached is a patch file for this model.V1806-Fix-stability-issue-in-dynamicLagrangian-model.patch

Thanks,

Ning Ren

Hi,

dynamicLagrangian model is not unconditionally stable. Yet, IMHO, if the denominator is heuristically conditioned as such, the derivation assumptions of flm and fmm , and the derivations themselves are likely violated.

When you do not let fmm->0, Eq. 3.5 of doi:10.1017/S0022112096007379 , that is Cs, will blow up anyways.

Have you played with the model's relaxation time scale, T?

Kind regards, Kuti

Hi Kutalmis,

Based on the Eq. 3.8 and 3.9 in that paper, the relaxation time scale may not solve the problem. The problem is caused by those zero velocity cells.

I am doing a fire modeling. Initially, the flow is stationary. Even when the fire is developed, very far away from the fire, the velocity is almost zero in some cells. In this case, Lij and Mij have zero values. From Eq. 3.8 and 3.9, both flm_ and fmm_ will be 0, regardless of the relaxation time scale. Now, we'll have the 0/0 problem and the code crashes. By adding, a VSMALL to fmm_, the code is very stable and will have negligible influence to the Cs value.

Best regards,

Ning

Hi Ning,

Many thanks for providing further information. Please see below for my informal remarks (if need be, further formal support can be requested through consulting services or software support services ):

• I presume that when you were mentioning 'the flow is stationary', you meant 'standstillness of the flow' (i.e. no fluid flow) rather than 'statistical stationarity of the flow'.

• dynamicLagrangian subfilter scale model is the traditional Smagorinsky subfilter scale model with dynamic Lagrangian-averaging procedure. In case of zero-valued velocity field -hence no turbulence, no fluid flow occurrence-, the Smagorinsky model's some of derivation assumptions are violated: e.g. Kolmogorov's hypothesis of local isotropy's sufficiently high Reynolds number presumption. Therefore, dynamicLagrangian may become prone to fail (though not necessarily failing) when a given case and its derivation scope are incompatible, seemingly such in this case. Actually, dynamicLagrangian is correctly and automatically going to zero when it encounters laminar flow/no flow, e.g. that's the reason why we suppress the Smagorinsky model's Cs towards walls where laminar shear flow occurs.

• I listed below some of the risks involved with your coding practice suggestion:

  • The spatiotemporal field of Cs is a function of the spatiotemporal fields of flm and fmm. When we force fmm locally very small, we inherently force the local Cs to be very large. Nevertheless, for laminar no-shear flows, Cs should go to zero (as dynamicLagrangian is trying to do herein). If we allow Cs to be very large locally, it may allow a computation numerically more stable; however, turbulent momentum fluxes are then allowed to exist/being large, where flow is not turbulent.
  • I can agree or disagree for a number of reasons but we would need formal quantifications/controlled studies for such heuristic tuning's effects on flow developments before we can give a definite answer. So far, to my knowledge, among peer-reviewed publications, such a scenario was not studied. And in my humble opinion, we need to be very cautious to make non-peer-reviewed/published changes to a given model.

• Having said that, I kindly offer two potential solutions:

  • Give a try to dynamicKEqn since it is not prone to 'Floating Point Exceptions' due to non-existence of a field in denominator.
  • The relaxation time scale T controls the memory length of the Lagrangian averaging, and several choices can be made. Also, the default T is the only heuristic part of the formal model. If you increase it, you may allow zero-velocity regions' adverse effects may blurred into non problematic regions.

Many thanks for your raising this question. This might be a good research question to improve this particular model.

Kind regards, Kuti

Disclaimer: These opinions are not officially approved or endorsed by ESI Group or OpenCFD Ltd.

Kuti,

Thanks for you in-depth discussion. I see this as simply a numerical issue where division by zero causes a floating point exception. This sort of fix that Ning proposes happens elsewhere all throughout the OpenFOAM source code where division by zero may occur. A quick search for VSMALL will yield thousands of instances where this operation is performed. Ning is not proposing to fundamentally change the dynamicLagrangian model. He is simply proposing to add VSMALL to the denominator in this fashion: flm_/(fmm_+VSMALL)

As for "formal support", we currently have such a support contract in place with OpenCFD and should be covered as such.

Karl

Hello Karl,

Many thanks for your highlights, and perspective.

I kindly then also suggest to Ning's attention to trap FP exceptions, and clip the local Cs to zero thereat. IMHO, this could prevent any violation on model derivation (i.e. avoiding high Cs across laminar regions), and ensure numerical stability (likely).

Kind regards,

Kuti

How would one go about "trapping FP exceptions"? I see that in etc/bashrc one can set FOAM_SIGFPE, but what happens when that is enabled?

assigned to @kuti

changed due date to January 31, 2019

changed due date to January 6, 2019

changed due date to January 13, 2019

added enhancement feature request models suggestion labels

mentioned in commit 8a63a8c87b1e1fd3ceab538b951ff373343961bb