Resolve instabilities under supercritcal flow

[RegMeasures]

Possibly by relaxing part of the momentum equation at high froude numbers or by implementing an adaptive timestep?

[RegMeasures]

One option to improve stability in transcritical/supercritical flows which is used in Flood Modeller Pro (previously ISIS) is to phase out part of the momentum term of the St-Venant equation when the froude number exceeds a threshold:

Supercritical flow can be modelled approximately in Flood Modeller. This is achieved by neglecting the dA/dx part of the convective momentum term in the momentum equation when the Froude number exceeds a specified upper value. Between this upper value and a specified lower value, the term is gradually phased out so that a smooth transition is achieved.

For steady supercritical flow in a uniform channel, the method should be acceptably accurate but will become more approximate as the channel becomes more non-uniform.

This approach is adequate for problems where supercritical flow prevails locally in isolated areas of a model and when low flows are required as initial conditions for an unsteady run.

Since Flood Modeller solves the differential form of the momentum equation, the solution at a hydraulic jump or bore can never be accurate. Instead of a sharp change in stage, the change will be spread over several nodes.

http://help.floodmodeller.com/floodmodeller/#t=Technical_Reference%2F1D_Modelling_Theory%2FSupercritical_Flows.htm

The default thresholds applied in Flood Modeller are:

• 0.75 = Value of Froude number below which there is no simplification of the St Venant momentum equation
• 0.9 = Value of Froude number above which part of the convective acceleration term is removed so that supercritical flow can be modelled.

http://help.floodmodeller.com/floodmodeller/#t=Technical_Reference%2F1D_Nodes_Reference%2FGeneral_System_Parameters.htm

[RegMeasures]

Response from Konrad Adams re how flood modeller pro implements this:

Indeed the phasing out of the term is as you say. It’s calculate the Fr, then If Fr<=0.75, term is fully present If Fr>=0.9, term is fully absent Else the term is linearly phased out between the two ( = (0.9-Fr)/(0.9-0.75) * “theoretical term” ) (The 0.75 and 0.9 are defaults, but can be changed; obviously(?) the upper bound cannot be 1 or greater since the thing goes mad for Fr=1

As for a reference, I’m afraid I couldn’t find anything in the usual suspects, and tbh don’t know its provenance. A quick Google didn’t reveal much instantly – the closest I could find is this, which analyses the MIKE method, which is similar method by the looks (uses a different weighting term, though I think we only phase out part of the term that is mentioned here https://www.mdpi.com/2073-4441/8/12/562/pdf ).

[RegMeasures]

Some reading on the subject following from Konrad's suggested reference turns up several useful papers regarding solving 1D transcritical flow with a Preissman scheme:

Abebe Y.A., Seyoum S.D., Vojinovic Z., Price R.K. (2016) Effects of reducing convective acceleration terms in modelling supercritical and transcritical flow conditions. Water 8(12).

As suggested by Konrad. It is a benchmarking type study which compares various models and how they perform for different test cases. Potentially I could replicate some of the tests for my code? It contains contains a clear description of the mathematics behind different options for simplifying the St-Venant equations as well as useful references which led me to the following papers.

Kutija, V. (1993) On the numerical modelling of supercritical flow. Journal of Hydraulic Research 31, 841–858.

This appears to be the original reference suggesting relaxation of the convective acceleration terms in the momentum equations. I don't have access to it so have requested an interloan from the library.

Djordjević S., Prodanović D., Walters G.A. (2004) Simulation of Transcritical Flow in Pipe/Channel Networks. Journal of Hydraulic Engineering 130(12):1167–1178.

This paper investigating using different Froude no. thresholds for triggering relaxation of the convective acceleration term, but it relaxes the whole term rather than only part of it.

Interestingly it also uses a simplified Preissman discretisation which would simplify the maths somewhat - see eq(8).

Sart C., Baume J.P., Malaterre P.O., Guinot V. (2010) Adaptation of Preissmann’s scheme for transcritical open channel flows. Journal of Hydraulic Research 48(4):428–440.

This well explained paper presents an alternative more accurate (and more complex) method of maintaining stability while solving the full 1D St-Venant equations using a Preissmann scheme. The approach involves identifying trans-critical reaches and replacing the standard st-venant equations for these reaches with equations describing internal boundary conditions. The paper has a good intro summarising the problem and other solutions (i.e. referencing some of the above papers). E.g.

A classical procedure allows the treatment of transcritical flows with a numerical method created for subcritical flows (Kutija 1994). The solution is to apply the method to a degenerated system of equations, reducing the influence of the convective momentum term as the flow becomes supercritical. But this trick does not conserve volumes.

[RegMeasures]

Based on the above reading I realise that it is possible to split the convective acceleration term of the momentum equation by applying the product rule. Based on the flood modeller pro manual it sounds like they do something like split the convective acceleration term:

Then relax the dA/dx part of it as the froude number increases.

[RegMeasures]

But Kutija splits/reduces slightly differently - I think like this: relaxing the first (du/dx) term.

Also they relax it always, rather than just for trans- and super-critical conditions.

[RegMeasures]

Had another idea to test a quasi-steady outlet channel (see quasi-steady branch) but it didn't seem to work.

[RegMeasures]

Setting "Beta" to zero has the effect of fully relaxing the convective acceleration terms of the momentum equation resulting in a universally stable model... I'm reasonably happy with this solution for now. I need to write this up properly and then I can probably close this issue.

[RegMeasures]

Fully relaxing the convective acceleration term is unfortunately not working as it massively underpredicts water level in the lagoon. Hence back to investigating better solutions...