強風飛来物衝突確率評価コードTOMAXIの開発とその適用例
公開日:
カテゴリ: 第15回
強風飛来物衝突確率評価コード TOMAXI の開発とその適用例
Development of Wind-borne Missile Strike Probability Evaluation Code TOMAXI and Its Application
電力中央研究所
江口
譲
Yuzuru EGUCHI
Member
電力中央研究所
平口
博丸
Hiromaru HIRAKUCHI
Member
電力中央研究所
村上
貴裕
Takahiro MURAKAMI
Member
電力中央研究所
杉本
聡一郎
Soichiro SUGIMOTO
Member
電力中央研究所
服部
康男
Yasuo HATTORI
Member
Injection of objects into high wind and flight in the air are numerically simulated using a random orientation model in order to mimic the tumbling of the objects. The computed random trajectories are used to evaluate annual strike probability on a target by considering the annual exceedance probability of the wind speed at the points where the objects are initially placed. Assuming that the direction of tornado translational path or the direction of straight wind is stochastically random, such numerical scheme is compiled in a computer code, named as TOMAXI. The TOMAXI code enables us to evaluate annual strike probability for a vertical plate or a horizontal plate of unit area placed in an arbitrary position. As an example of the application, the code has been applied to an automobile initially staying on the ground to compute the spatial distribution of annual strike probability for a vertical and horizontal plate of unit area.
Keywords: High Wind, Missile Strike, Probability, Tumbling, Tornado, Typhoon
１．序論
新規制基準では、竜巻が設計・運転で考慮すべき自然現象として新たに設定され、「原子力発電所の竜巻影響評価ガイド」（以下「竜巻評価ガイド」）が原子力規制委員会によって制定された[1]。この竜巻評価ガイドの設計飛来物速度の例示値は、決定論的な評価方法に基づき計算されたものである。一方、飛来物による損傷確率P は、P
＝（竜巻来襲確率）×（飛来物存在確率）×（浮上飛散確率）×（衝突確率）×（衝突による損傷発生確率）の形式で表現することができる。竜巻来襲確率は 10-5/年以下のオーダーであるので、竜巻飛来物による損傷発生は確率的にさらに小さいものとなる。米国規制ガイドNUREG1.76 [2]の設計飛来物速度はこのような確率的な側面を部分的に考慮して設定されていることを Antaki[3] が指摘している。NUREG1.76 [2]の設計飛来物速度は、日本の「竜巻評価ガイド」[1] の設計飛来物速度よりも小さいが[4]、その理由の一つとしてこのような確率的な考慮の有無が影響しているものと考えられる。
て、飛来物が安全関連施設へ衝突する確率が十分に小さ いならば、その影響を無視することができるものと考え られる。飛来物の衝突確率の評価法については、Twisdale ら[5]、Goodman and Koch [6]、Hope ら[7]が解析モデルや解析コードを開発している。著者らは、竜巻や台風など の強風の方向依存性が統計的に一様である場合を想定し、比較的簡易な評価方法を提案した[8]-[10]。その方法は、 強風飛来物の年間衝突確率評価コード TOMAXI[9]として取り纏められている。本論文ではその概要と必要とな る入力データを説明し、具体的な適用例を示す。
２．計算モデル
2.1 衝突確率の全体計算フロー
被衝突体（ターゲット）が飛来物の初期位置（原点） から見て位置r にあるとし、飛来物の初期位置での上空風速V（地上10m での3 秒平均風速相当）の年間超過確率がH(V)で表されるものとする。この場合、上空風速V の確率密度関数p(V)は以下で表される。
竜巻飛来物による安全関連施設の損傷確率は、飛来物が安全関連施設へ衝突する確率がその上限となる。従っ
p(V )
?H (V )
?? ?V
(1)
連絡先: 江口 譲、〒270-1194 千葉県我孫子市我孫子1646、一般財団法人 電力中央研究所、E-mail: eguchi@criepi.denken.or.jp
確率密度関数p(V)は種々の強風特性（竜巻強度、竜巻半径、移動速度や台風の特性など）をすべて包含した値であるため、確率密度関数p(V)は強風特性i の強風ごとに分解することができる。つまり、確率密度関数p(V)は強風特性i に関する確率密度関数pi(V)の総和として、以下の関係で表される。
N
2.2 竜巻ハザード計算方法の概要[12], [13],[15]
以下では、竜巻に限定して竜巻風速のハザード曲線の計算法の概要を説明する。ある地点において、1 年間の期間において竜巻による風速がV を超える事象をTVと表すと、事象TVが発生する確率H(TV)= H(V)は以下で表される。
?
p(V ) ? ? pi (V )
i?1
(2)
H (TV ) ? ? P(TV | N )Pa (N )
N ?0
(5)
ここで、N は強風特性i の分類総数を表し、i=1, 2,…, N である。いま、強風特性i の強風によって飛来物の初期位置での上空風速がV となる場合に、飛来物が浮上・飛散し位置r にある単位面積の平板状のターゲットに衝突する条件付確率をq i(V, r)で表すことにする。このとき、位置r にある単位面積の平板状のターゲットにこの飛来物が衝突する年間確率Ph(r)は以下で表すことができる。
ここで、N はこの地点と同様の竜巻発生特性を有する地域で1 年間に発生する竜巻の数を表し、P(TV |N)はこの地域で1 年間にN ヶの竜巻が発生した場合に事象TVが発生する条件付確率を表す。また、Pa(N)はこの地域で1 年間にN ヶの竜巻が発生する確率を表し、以下のポリヤ分布でよく表現されることが知られている[13]。
? N? N ? 1 N ?1
Vm N
? ?
Ph (r) ? 0pi (V )qi (V ,r)dV
i?1
(3)
Pa (N ) ?
(1 ? ??)
N !
? ?(1 ? ?i)
i?1
(6)
ここで、Vm は上空風速V の最大値である。いま、条件付衝突確率q i(V, r)が強風の特性（竜巻強度、竜巻半径、移動速度など）に顕著に依存しないと仮定すると、すべ
(
ここで、総積?の添え字の上限がN?1= ?1, 0 の場合は総積の値は1 をとるものとする。?は竜巻発生数N の平均値、?は竜巻発生数N の分散? 2 と以下の関係にある[14]。
ての強風特性i についてq i(V, r)= q(V, r)が近似的に成立す
? ? 1 ? 2
?1)
(7)
ることになる。この関係を用いると、年間衝突確率 Ph(r)
の評価式(3)は以下のように表される。
? ?
なお、? 2= ?の場合は? = ?となり、式(6)はポアソン分布
Ph (r) ?
Vm
0 p(V )q(V ,r)dV
?
(4)
に収斂する。
一方、条件付確率P(TV
|N)は以下のように表される。
従って、年間衝突確率Ph(r)を評価するには、上空風速V の確率密度関数p(V)と条件付衝突確率q (V, r)が必要と
P(TV
| N ) ? 1 ? (1 ? R)N
(8)
なる。年間衝突確率評価コード TOMAXI では、Fig.1 に示すように、条件付衝突確率q (V, r)の計算には竜巻飛来物浮上・飛散解析コード TONBOS[11]などを用い、確率
ここで、R(V)は１つの竜巻が発生した場合に、当該地点での上空風速がV を超える平均的な確率であり、以下で表される。
密度関数p(V)の計算にはTOWLA[12]などを用いる
R(V ) ? A?V
S
(9)
ここで、S はこの地域の面積を表す。また、A>V はこの地域で発生する１つの平均的な竜巻の中で風速がV を超える領域の面積を表す。
式(6), (8)を式(5)に代入すると、以下が得られる。
?
H (T ) ? 1 ? ?(1 ? R)N P (N )
Va
N ?0
(10)
Fig.1 Computing flow of annual strike probability in TOMAXI
?
ここで、? Pa (N ) ? 1 を考慮した。
N ?0
ところで、関数f(x)=(1+ ? x)?1/?をテーラー展開すると、以下が得られる。
? (?x)N N ?1
k：上向き鉛直方向単位ベクトル
f (x) ? ?
N ?0
?(1 ? ?i)
i?1
N !
(11)
q：単位揚力方向ベクトル
この関係式を用いると、式(10)右辺第２項は以下のように表される。
?? 1
s：単位横力方向ベクトル
?：空気密度
CD：飛来物の抗力係数
?(1 ? R)N P (N ) ? (1 ? ??R) ?
a
N ?0
(12)
CL：飛来物の揚力係数CS：飛来物の横力係数
従って、ある地点において、事象 TV が発生する確率
H(TV) = H(V)は以下で表される。
A：飛来物の代表見附面積（=Ld）
? 1
H (T ) ? H (V ) ? 1 ? (1 ? ??R(V )) ?
V
(13)
ここで、ベクトルvre, q, s は互いに直交し，物体の姿勢と相対風速ベクトルに依って変化するものである。
なお、?R が1 より十分小さい場合はH(V)=?R(V)で近似的に表される。?、?およびR を算出するには、竜巻等の記録や竜巻特性等の統計モデルに基づき数値解析する必要がある[12], [15]。
抗力・揚力・横力の更新間隔については、Twisdale ら [5] は経験的に求めた一定値を使用しているが、著者ら[17] は相対風速の大小が自励回転の周期に影響を及ぼすことを考慮し、以下の回転の基礎式に基づく物体姿勢更新方法を新たに提案した。
浮上・飛来計算モデル[10], [17]
?2?
?C A? ? v ? v 2
? k DW
(16)
強風により空中に浮上・離脱した物体は、物体が浮上
または支持部から離脱する際の初期姿勢の違いや空中での自励回転(autorotation)等の影響を受けて、一般的に不規則な回転を伴いながら飛散する。そこで、著者らはTwisdale ら [5]が提案したランダム回転モデルを新たに
?t2mL
ここで、k は定数を、?は回転角を、L は飛来物の最長寸法を表す。また、[CDA]は丸山ら[18]の評価モデルで用いられている飛行定数[CDA]/m と同様に、L, d, b を飛来物の各辺の寸法として以下の式で定義されるものである。
TONBOS[11]に導入し、TONBOS[(ver.4)とした[10], [17]。[C A] ? CLd Ld ? CLb Lb ? Cdbdb
このモデルでは、乱数を用いてある時間間隔ごとに物体D3
(17)
の空中姿勢をランダムに更新し、各時刻の物体姿勢に対応する抗力・揚力・横力を作用させて、並進運動のみを解析する。以上のモデル化に基づく並進運動に関する方程式は以下のように記述される。
なお、CLd 等は各方向から風を受けた場合の抗力係数を表す。式(16)を?=0 を初期条件として時間積分し、?= ?となった時点で物体姿勢を更新し、更新後は?=0 にリセットして同様の更新を繰り返す。また、風速に関しても乱
m dVM
dt
? f ? mgk(14)
数を用いて同じ更新間隔で変動させた。
その他のモデル化等については、Table1 に示すものを
ここで、空気力ベクトルf は以下のように表される[5]。
2
用いた。モデルの詳細については、文献[10], [17]に記載し
f ? 1
2
? VW ? VM
A?CDvre ? CLq ? CSs?
(15)
ている。
なお、各変数は以下で定義されるものである。
モデル化項目
モデル化法
飛来物初期姿勢
乱数によるランダムな初期姿勢
空中での飛散挙動
ランダム回転モデルによる並進運
動方程式を時間積分、式(14)
空力係数
クロスフロー理論
空中姿勢と風速変動の更新間隔
回転の方程式,式(16)を時間積分
し、回転角?が?となるごとに姿勢と風速変動を更新する
m：飛来物の質量t：時間 g：重力加速度f：空気力ベクトルVW：風速ベクトル
VM：飛来物速度ベクトルVW ?VM：相対風速ベクトル
vre：単位相対風速（抗力方向）ベクトル
Table 1 Physical model and numerical methods
TOMAXI の確率計算モデル[9]
Fig.2 のように原点に１つの飛来物があるものとし、進路方向が?の１つの竜巻Ａが通過する場合を考える。この竜巻Ａによって原点での上空風速がV であるとする。この竜巻Ａと同じ竜巻特性と進路方向?を有し、発生点が異なる竜巻Ｂなどによっても原点での上空風速がV であるとする。このような竜巻の発生点の集合（領域）を?T(?) で表すことにする。この場合、領域?T(?)を発生点とする竜巻によって、原点の飛来物は種々の方向へ異なる軌跡で飛散する。
いま、Fig.3 のように原点を発生点としてX 方向に進む竜巻によって、上空風速がV である地点の集合（領域）
?M にある飛来物が飛散する現象を考える。この飛散現象はFig.2 の飛散現象と相対的に同一となる。つまり、両者の軌跡の座標は平行移動と回転移動によって変換可能となる。また、この変換は任意の進路方向?（Fig.2）に対して可能であるので、すべての進路方向?の竜巻による飛散現象は、Fig.3 の飛散現象に集約できる。
安全上重要な機器・構造物などのターゲットが飛来物からr（水平距離L、地上高H）の位置にあるものとする。また、竜巻の進路方向は統計的に一様と仮定すると、このターゲットへの単位面積当たりの衝突確率q(V, r)は、Fig.4 に示すようなリング状のターゲットの単位面積当たりの衝突確率と同一となる。つまり、Fig.5 のように各飛来物に対応した仮想的なリング状のターゲットへの衝突確率を計算し、その面積で除すことにより、衝突確率q(V, r)を計算することができる。例えば、ターゲットが鉛直面である場合、Fig.6 のように各飛来物からの軌跡を(r,z)面に射影して、この軌跡がリング状のターゲットを通過する割合を計算し、これを面積で除すことにより、以下のように条件付の衝突確率q(V, r)を計算することができる。
Fig.2 Spatial relation between a missile and tornadoes
Fig.3 Transformed relation between missiles and a tornado
q(V ,r) ?
mstrike
2? L?hNall
(17)
ここで、?h はターゲットの高さを、Nallは領域?M内にある飛来物の総数を、mstrike は軌跡がリング状のターゲットを通過する個数を表す（mstrike≦Nall）。式(17)の関係を用いると、年間衝突確率Ph(r)は式(4)で求められる。ターゲットの位置r（水平距離L、地上高H）を空間的に変化させることにより、年間衝突確率 Ph(r)の空間分布を計算することができる。
Fig.4 Plan view of a target location relative to a missile
Fig.5 Plan view of virtual target locations relative to missiles
Fig.6 Trajectories of missiles from a region ?M
３．適用例
地上に飛来物として自動車（寸法：4.4mx1.7mx1.5m, 質量：1,140kg）が１台あるものとし、最大竜巻風速 89m/s
（移動速度 18m/s，旋回速度 72m/s，竜巻半径 45.7m）の竜巻がランダムに来襲するものとした。まず、一辺が竜巻半径の４倍の正方形領域に約 50,000 台の自動車を等間
隔で配置して飛散軌跡を 2.3 節に記載した方法で計算した。この軌跡データをTOMAXI の入力データとして、2.4 節の方法で条件付の衝突確率 q(V, r)を計算した。一方、竜巻風速V の年超過確率H(V)はH(V)=1.52 exp(-0.10466V) で表されるものと仮定した。
この年超過確率H(V)と条件付の衝突確率q(V, r)を用いて、年間衝突確率Ph(r)を式(4)で求めた。その結果、単位面積の水平面に対する年間衝突確率の空間分布はFig.7(a) に、鉛直面はFig.7(b)に示すものとなった。
Fig.7(a) Annual strike probability for horizontal plate
Fig.7(a) Annual strike probability for vertical plate
４．結論
本論文では、竜巻や台風などの強風によって物体が浮上・飛散し、安全上重要な機器や構造物に衝突する年間確率を計算するための評価方法（理論的背景・解析フロー・解析コードの概要）を説明した。
特に、浮上・飛散する物体の姿勢や風速のランダム性を考慮して確率的に軌跡を評価する方法（TONBOS(ver.4) で用いられている方法）や竜巻ハザード曲線を評価するための方法について、概要を説明した。確率的な飛来物の軌跡を風速ハザード曲線と組み合わせ、さらに竜巻進行方向の一様性を仮定することにより、物体が平板に衝突する年間確率の空間分布を解析コード TOMAXI を用いて計算することができることを本論文で示した。
参考文献
原子力規制委員会，原子力発電所の竜巻影響評価ガイド（実用発電用原子炉の規制基準に関連する内規，平成25 年6 月19 日 原規技発第13061911 号にて制定， 平成26 年9 月17 日 原規技発第1409172 号にて改正）, 2013.
U.S. Nuclear Regulatory Commission, “Regulatory Guide 1.76, Design-Basis Tornado and Tornado Missiles for Nuclear Power Plants”, Revision 1, March 2007.
Antaki, G., “Protection of Nuclear Power Plants from Tornado Missiles”, Becht Engineering Blog, 2015. (https://becht.com/blog/protection-of-nuclear-power-plants-fr om-tornado-missiles)
江口 譲, 竜巻飛来物バザード評価におけるシミュレ
ーション利用の現状, 日本原子力学会 2016 秋の大会
（計算科学技術部会企画セッション）, PL2L03, 2016.
Twisdale, L.A., Dunn, W.L. and Davis, T.L., “Tornado Missile Transport Analysis”, Nuclear Engineering and Design, 51, 295-308, 1979.
Goodman, J. and Koch, J. E., “The Probability of a Tornado Missile Hitting a Target”, Nuclear Engineering and Design, 74, pp.125-155, 1983.
Hope, K.D., Povroznyk, N. and Schneider, R., “Tornado Missile Strike Calculator: An Excel-based Stochastic Model of Tornado-Driven Missile Behaviour for Use in High Winds PRA”, International Topical Meeting on Probabilistic Safety Assessment and Analysis (PSA 2015), Sun Valley, Idaho, USA, April 26-30, 12086, 2015.
Eguchi, Y., Murakami, T., Hirakuchi, H., Sugimoto, S. and Hattori, Y., “An Evaluation Method for Tornado Missile Strike Probability with Stochastic Correlation”, Nuclear Engineering Technology, 49, Issue 2, pp.395-403, 2017.
Eguchi, Y., Sugimoto, S., Hattori, Y., Murakami, T., and Hirakuchi, H., “An Estimation Method for Tornado Missile Strike Probability under Assumption of Statistically Isotropic Tornado Path Directions”, pp.241-246, Proc. of 2017 International Topical Meeting on Probabilistic Safety Assessment and Analysis (PSA 2017), Pittsburgh, USA,
2017.
Eguchi, Y., Sugimoto, S., Hattori, Y., Murakami, T. and Hirakuchi, H., “Computation of Annual Strike Probability of a Wind-borne Tumbling Missile using TOMAXI”, Proc. the 14th International Probabilistic Safety Assessment & Management (PSAM14), 2018 (to appear).
江口 譲，杉本 聡一郎，服部 康男，平口 博丸, 竜巻
による物体の浮上・飛来解析コードTONBOS の開発， 電力中央研究所 研究報告N14002, 2014.
平口 博丸, 野原 大輔,杉本 聡一郎，江口 譲，服部康男、沿岸立地原子力発電所の竜巻風速ハザードモデル TOWLA の開発、電力中央研究所 研究報告書O15005, 2016.
Wen, Y.K. and Chu, S.L., “Tornado risks and design wind speed”, Proceedings of American Society of Civil Engineering, Journal of Structural Division, Vol.99, No.ST12, pp.2409-2421, 1973.
印東 太郎、確率および統計、p.113、コロナ社、1957.
Garson, R.C., Morla-Catalan J. and Cornell C.A., “Tornado risk evaluation using wind speed profiles”, Journal of Structural Division, Proceedings of American Society of Civil Engineering, Vol.101, No.ST5, pp.1167-1171, 1975.
Twisdale, L.A. and Vickery, P.J., Comparison of Debris Trajectory Models for Explosive Safety Hazard Analysis, 25th DoD Explosive Safety Seminar, Anaheim, California, 1992.
江口 譲，村上 貴裕，杉本 聡一郎，服部 康男，平
口 博丸,ランダム回転モデルによる強風下での角柱飛散解析, 日本流体力学会2017, 2017.
丸山 敬，河井 宏允，奥田 泰雄，中村 修, 竜巻中の飛散物の特性に関する数値計算，京都大学防災研究所年報 57(B), pp.248-259, 2014.