大阪湾の縦式単節静振に就いて、日高孝次
大阪湾を長軸の一端から他端に向かって次第に深くなっている楕円形の海として取り扱う。大阪湾の長径を 61.0km、平均の深さを27m、ｇ＝9.779m/s2とすれば縦式静振の周期として単節静振では 118.3min., 双節静振では 69.6min. を得た。尚縦式単節静振では大阪湾北東部の振幅は、その南西部の振幅の約3倍になるという結論を得た。これらは実測値と良く一致している。
基本方程式
∂2ζ／∂ｔ2＝（ｇ／ｗ）∂（ｗｈ∂ζ／∂ｘ）／∂ｘ (1)
ここでζ：水位、ｈ：水深、ｗ：水路幅、ｘ：楕円の長径方向の座標。
水路を長軸ａ、短軸ｂの楕円と考え、
ｗ＝２ｂ（１−（ｘ／ａ）2）1/2 (2)
水深は直線的に変化すると考える。
ｈ＝ｈ0（１−ｘ／ａ） (3)
ｈ0は平均水深。
振動を周期Ｔの単弦運動とする。
ζ＝η（ｘ）ｃｏｓ（２πｔ／Ｔ） (4)
ここで
ξ＝ｘ／ａ (5)
と置けば方程式(1)は
（１−ξ2）（ｄ2η／ｄξ2）−（１＋２ξ）（ｄη／ｄξ）＋（ｋａ）2（１＋ξ）η＝０ (6)
ただし
ｋ2＝４π2／（ｇｈ0Ｔ2） (7)
方程式(6)は（ｋａ）2のある特定な値に対してのみ収斂する解を与える。（固有値問題）
境界条件
（ｄη／ｄξ）ξ=+1≠∞ (8)
（ｄη／ｄξ）ξ=-1＝０ (9)
方程式(6)を解くため、次のような冪級数展開を用いる。
ｄη／ｄξ＝Ａ0＋Ａ1ξ＋Ａ2ξ2＋・・・ (10)
と置く。これを(6)に代入すれば
（ｋａ）2Ａ-1−Ａ0＋Ａ1＝０ (11)
（ｋａ）2Ａ-1−（（ｋａ）2−２）Ａ0−Ａ1＋２Ａ2＝０ (12)
およびｎ≧１に対し
（ｋａ）2Ａn-1／ｎ−（ｎ＋２−（ｋａ）2／（ｎ＋１））Ａn−Ａn+1＋（ｎ＋２）Ａn+2＝０ (13)
これを解くのに日高は以下のようにしている
Ａn+1／Ａn＝Ｎn+1 (14)
と置き、(11)と(12)からＡ1を消去したものから
（１−（ｋａ）2）／（２Ｎ1）＋１＋Ｎ2＝０ (15)
および(13)から
（ｋａ）2）／（ｎＮn）−（ｎ＋２−（ｋａ）2）／（ｎ＋１））−Ｎn+1＋（ｎ＋２）Ｎn+1Ｎn+2＝０ for ｎ≧１ (16)
また(11)から
Ａ-1＝（１−Ｎ1）Ａ0／（ｋａ）2） (17)
及び(14)によって
Ａ1＝Ｎ1Ａ0
Ａ2＝Ｎ1Ｎ2Ａ0
・・・・
Ａn＝Ｎ1Ｎ2・・ＮnＡ0
であるから、
ｄη／ｄξ＝Ａ0（１＋Ｎ1ξ＋Ｎ1Ｎ2ξ2＋・・・） (19)
等が得られる。
今、条件(8)によって(19)はξ＝１において収斂しなければならないから
Limn→∞Ｎn+1
なる極限値は１より小さい一定値に収斂しなければならない。今(16)からこれは０か１であることは容易に判る【原文ママ】。したがって
Limn→∞Ｎn＝０ (21)
なる極限値を採用しなければならない。
このようにしてＮnを求める。そして条件(9)にしたがって
１−Ｎ1＋Ｎ1Ｎ2−Ｎ1Ｎ2Ｎ3＋・・・＝０ (24)
を満足するように（ｋａ）2を選ぶ。
ここで、ｋａの大体の値は次の様にして目安を付けた。
深さが一様な矩形の海では
ｋａ＝π／２，π，・・・
本問題の様に（水深が）変わる矩形の海では
Ｊ0'（２1/2ｋａ）＝０
の解で
ｋａ＝1.355, 2.480, ・・・
となるから、矩形の場合は深さ（が変わること）の影響は各 16%および 26%である。楕円形の場合も等深のときのｋａは判っているからそれらより 16%および 26%小さい値を用いれば大体の値として考えることができる。
以上が日高が用いた計算アルゴリズムであるが、この計算を実際に行ってみよう。これは、Ｎnが（式(16)と同じ事であるが）
Ｎn＝（（ｋａ）2／ｎ）／（ｎ＋２−（ｋａ）2／（ｎ＋１）＋Ｎn+1−（ｎ＋２）Ｎn+1Ｎn+2） (22)
および式(24)を満たすような（ｋａ）2の固有値および対応する固有ヴェクトル｛Ｎn，ｎ＝１，２，・・・｝を求める問題である。現在ではこれはたとえば Microsoft Excelなどの「ゴールシーク」機能を用いれば比較的簡単に追算できる。
Oosakawan.xls はこれを実行するためのワークシートである。
[B1]には（ｋａ）2が与えられ、"Lambda"という名前が付されている。[B2]にはｋａが（[B1]から）求められる。
[D ]の列には順にＮ1，Ｎ2，Ｎ3，・・・が式(22)により求められる。ただし、最後の２つはＮn+1，Ｎn+2を０としている。[E ]の列には積Ｎ1Ｎ2・・Ｎnが、[F ]列には式(24)の右辺が求められる。
計算を行うためには、まず[B1]に適当な値（日高に従えば等深楕円の場合より16%および26%小さい値）を与える。次に、「ツール／ゴールシーク」を起動し、「数式入力セル」として[F ]列の最後のものを選び、「目標値」を０とし、「変化させるセル」として[B1]を選んで実行すれば正しい解が得られる。
以下に結果を示す。「Excelの結果(1)」は、日高と同様、単節静振は5次まで、双節静振は9次までを求めた場合、「Excelの結果(2)」はいずれも20次までの場合である。
単節静振のＮ1だけが日高の結果と少し差がある。それ以外はほとんど一致している。単節静振のＮ1は単なる計算間違いまたは誤植であろうか。
単節静振
| 日高の結果 | Excelの結果(1) | Excelの結果(2)
|
|---|
| （ｋａ）2 | 2.755 | 2.754 | 2.754
|
| ｋａ | 1.660 | 1.659 | 1.659
|
| Ｎ1 | 1.5040 | 1.5435 | 1.5436
|
| Ｎ2 | 0.43246 | 0.43218 | 0.43204
|
| Ｎ3 | 0.20945 | 0.20933 | 0.20935
|
| Ｎ4 | 0.12640 | 0.12633 | 0.12510
|
| Ｎ5 | 0.08424 | 0.08419 | 0.08367
|
単節静振のη／Ｈ
| ξ | 日高の結果 | Excelの結果(2)
|
|---|
| 1.0 | 1.000 | 1.000
|
| 0.8 | 0.668 | 0.668
|
| 0.6 | 0.399 | 0.400
|
| 0.4 | 0.186 | 0.186
|
| 0.2 | 0.019 | 0.019
|
| 0.0 | -0.108 | -0.108
|
| -0.2 | -0.203 | -0.201
|
| -0.4 | -0.267 | -0.266
|
| -0.6 | -0.308 | -0.307
|
| -0.8 | -0.330 | -0.329
|
| -1.0 | -0.337 | -0.336
|
双節静振
| 日高の結果 | Excelの結果(1) | Excelの結果(2)
|
|---|
| （ｋａ）2 | 7.966 | 7.967 | 7.967
|
| ｋａ | 2.822 | 2.823 | 2.823
|
| Ｎ1 | -0.579309 | -0.57869 | -0.57867
|
| Ｎ2 | 7.012410 | 7.01917 | 7.01945
|
| Ｎ3 | 0.940249 | 0.94040 | 0.94040
|
| Ｎ4 | 0.456507 | 0.45656 | 0.45656
|
| Ｎ5 | 0.280798 | 0.28083 | 0.28083
|
| Ｎ6 | 0.192979 | 0.19300 | 0.19299
|
| Ｎ7 | 0.141750 | 0.14176 | 0.14176
|
| Ｎ8 | 0.109244 | 0.10925 | 0.10895
|
| Ｎ9 | 0.086747 | 0.08675 | 0.08654
|