In free surface flows, the specific energy is a fundamental quantity governing the physics of relevant phenomena, e.g, Chow (2009). Recently, Computational Community (e.g Le Floch & Duc Thanh (2011), Castro et al. (2013), Bernetti et al. (2008)) has re-discovered its crucial role in solving some Shallow Water problems on discontinuous bottom. In the work of Valiani & Caleffi (2008) the analytical inversion of the depth - specific energy relationship is given.
The reference:
Valiani, A. & Caleffi, V., Depth-energy and depth-force relationships in open channel flows: Analytical findings (2008) Advances in Water Resources, 31 (3), pp. 447-454.
# preparation of the numerical environment
%matplotlib inline
from numpy import sqrt, cos, pi, arctan, linspace
from matplotlib import pyplot as plt
For a fixed value of the specific discharge $q$, the classical expressions of the specific energy $E$ as functions of the depth $Y$ is:
$$ E = Y + \frac{{q^2 }}{{2g\,Y^2 }}; $$where $g$ the is gravity acceleration, and can be interpreted as a third-order equation in $Y$.
The problem to solve is to analytically found the two pysically meaningful depths $Y_{sb}$ and $Y_{sp}$ (subcritical and supercritical, respectivelly) for given values of $E=E_0$ and $q=q_0$.
The analytical solution described in Valiani & Caleffi (2008) is given by the following steps:
where $\alpha = \arctan \left( {\sqrt {\Gamma _0^3 - 1} } \right)$;
The following python code is proposed to highlight the simplicity of the method.
def radix(q_0,E_0,g):
"""
Compute analitically the roots (Y_sb, Y_sp) of the equation:
E_0 - Y - q_0^2/(2 g Y^2) = 0
Parameters
----------
Input: | Output:
q_0 : float | Y_sb : float
[m^2/2] specific discarge | [m] subcritical depth
E_0 : float | Y_sp : float
[m] specific energy | [m] supercritical depth
g : float |
[m/s^2] gravity |
"""
Y_c = (q_0**2/g)**(1.0/3.0) # critical depth
E_c = 3.0/2.0*Y_c # crtitical energy
# check the esistence of physical solutions!
if E_0 > E_c:
Gamma_0 = E_0/E_c
else:
print('no physical solutions!')
return 0.0, 0.0
alpha = arctan(sqrt(Gamma_0**3 - 1.0))
# non-dimensional subcritical depth
eta_sb = Gamma_0 / 2.0 * (1.0 + 2.0*cos(( pi - 2.0*alpha)/3.0))
# non-dimensional supercritical depth
eta_sp = Gamma_0 / 2.0 * (1.0 + 2.0*cos(( pi + 2.0*alpha)/3.0))
Y_sb = eta_sb * Y_c # subcritical depth
Y_sp = eta_sp * Y_c # supercritical depth
return Y_sb, Y_sp
To test the solution we find the two physical roots of the specific energy function for $q$ = 2.00 [m$^2$/2] and $E$ = 2.50 [m]:
gg = 9.81 # [m/s^2] gravity
qq = 2.00 # [m^2/2] specific discarge
EE = 2.50 # [m] specific energy
Y_ss = radix(qq,EE,gg)
Y = linspace(0.2,3,100)
E = Y + qq**2 / (2.0*gg*Y**2)
plt.figure(2,figsize=(10, 7))
plt.title('Speficic Energy vs. Depth')
plt.plot(E,Y,label='E = E(Y)')
plt.plot((EE,EE),Y_ss,'or',label='solutions')
plt.xlabel('E [m]')
plt.ylabel('Y [m]')
plt.axis((1.0,4.0,0.0,3.0))
plt.grid('on')
plt.legend()
plt.show()
Bibliography: