各分割断面で流下方向に等流が成立すると仮定すると、
各分割断面で次式が成り立つ。
ここで、
分割断面間の流れ同士の干渉によるせん断力
また、
樹木群と流れの干渉によるせん断力
ここに、
河床摩擦損失はマニング則を用いて次式とする。
具体例として、次の断面のうち、n=3と2の分割断面について運動方程式を立ててみる。
分割断面n=3の運動方程式
分割断面n=2の運動方程式
すべての分割断面に対応できるように、一般化すると次式となる。
よって、分割断面の数だけ方程式が成立する。
なお、分割断面n=5は流速0のため、これを除いた5分割断面の計算を行なう。
さらに、以下の連続式を連立させることにより、各分割断面の流速と損失勾配
ここで得られた流速から求めるせん断力を流下方向の運動方程式に代入することにより縦断水位を計算する。
なお、流水間のせん断力は前述のとおり分割断面間でキャンセルされるため、計算しないことに留意する。
さらに、以下のように有効断面の外側に死水域を設定する場合は下式のように境界部分のせん断力を考慮する。
境界混合係数の標準値は以下のように定義されている。
私の感覚では厳密にこの図に則って計算している事例は少ないように思う。
高水敷・低水路間だけしっかり決めてあげれば良い気がする。
一般断面の基礎式を用いて、限界水深、等流水深を矩形断面と同様に定義で設定することは難しいため、便宜的に以下のとおりに設定する。
全ての損失による水頭の勾配が河床勾配
なお、一般断面では水深を用いないため、正確には等流時の水位となる。
となり、これを満足する水位
フルード数が1となる水位を限界流時の水位とする。
フルード数はエネルギー保存則の分母が0より、
となるが、平方根の中に水深
そのため、分割断面法では水深
横断方向流速分布の式は以下となり、分割断面数分の式が成り立つ。
以下の連続式を連立させることにより、各分割断面の流速と損失勾配
縦断水位の離散化はこれまでと同様に次式となる。
なお、
常流の場合、下流から逐次計算を行なうため未知数は
縦断水位はこれまで同様のため説明を省略し、横断方向流速分布の計算方法をのみを示す。この計算にはかなり高度な数値計算手法が必要なため、理解できなくても問題ない。
計算するための前知識として以下の2手法を紹介する。
高次のニュートン法は、通常のニュートン法と同様に以下のように示される。 なお、ボールド体は行列を示す。
上式の右辺第二項を
よって、下式で解の更新を行なう。
TDMA法(トーマス法)とは係数行列が三重対角行列の連立一次元方程式の効率的な数値解法である。
ガウスの消去法(直接法)の計算量は分割数Nの3乗であるが、TDMA法はNの数倍程度の計算量となる。
導出はwikipediaを参照ください。https://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm
def TDMA(a,b,c,d):
x, P, Q = np.empty_like(a), np.empty_like(a), np.empty_like(a)
n = len(x)
P[0], Q[0] = -c[0]/b[0], d[0]/b[0]
for i in range(1, n):
P[i] = - c[i] /(a[i]*P[i-1] + b[i])
Q[i] = (-a[i]*Q[i-1] + d[i])/(a[i]*P[i-1] + b[i])
x[-1] = Q[-1]
for i in range(n-2, -1, -1):
x[i] = P[i]*x[i+1] + Q[i]
return x
参考: https://nekodamashi-math.blog.ss-blog.jp/2017-12-03-2
def lmTDMA(a,b,c,d):
n = len(a)
#! 前進消去
for i in range(1, n):
ratio=a[i]/b[i-1]
b[i]=b[i] - ratio*c[i-1]
d[i]=d[i] - ratio*d[i-1]
#! 後退代入
d[-1]=d[-1]/b[-1]
for i in range(n-2, -1, -1):
d[i]=(d[i]-c[i]*d[i+1])/b[i]
return d
分割断面の運動方程式を
以下のとおりに式変形できる。
ここに、
全分割断面の方程式をベクトル形式で記述すると以下となる。なお、分割断面5は、水深が0となるため、省いて行列を構成する。
そのため、高次のニュートン法を用いて反復計算を行なう。
なお、上式は三重対角行列となるため、TDMAにより容易に逆行列が計算可能である。
ニュートン法による計算結果の
よって、以下の式により
上記からわかるとおり、本手法では水位が高めに計算されがちである。
本手法には多くの問題点あるが、より正確に水位を計算したいという考えは正しいと思われる。
本手法:平均流速レベル3による一次元(準二次元)不等流計算方法について説明を行った。
上式の右辺を次式と定義する。 $$ \begin{align} \dfrac{T}{\rho g A} \equiv I \end{align} $$ ここに、$I$:損失勾配とする。 若干の式変形を行なうと次のとおりとなる。 $$ \begin{align} \dfrac{\displaystyle \sum^n T_n}{\rho g } = I \sum^n A_n \end{align} $$