(57)【要約】【課題】  解析に用いる式の自由度を低減させる。【解決手段】  ｘｙｚ座標空間内にてｚ軸について軸対称とした円筒状のタンク１の解析モデルを設定する。流体運動の解析領域とするタンク１の周囲の空間については、ｘ＝ｘｍｉｎに設定した爆風源における地面上に位置する点と、タンク１の底面の中心とを結ぶ直線を軸とする円筒座標を設定する。これにより、流体運動の解析領域を、ｘ軸方向におけるタンク１が存在する区間においてタンク１の外面に対応する位置に内半径ｒｍｉｎが変位する非一様の空間として捉え、流体の運動を、この非一様の空間に爆風源で発生した爆風の圧力波が伝播する問題として捉えることで、流体運動の解析領域における流体の流速成分と動圧について、一般化座標に関する時間の常微分方程式を求め、この常微分方程式を解くことで、タンク１の外面に対応する位置に作用する流体の圧力を求めるようにする。【選択図】図１
2015022599.tif 000002
【特許請求の範囲】
【請求項１】
コンピュータにより、爆風によってタンクに作用する圧力を予測する爆風によりタンクに作用する圧力の予測方法であって、
ｘｙｚ座標空間内に、ｚ軸について軸対称な円筒状のタンクの解析モデルを設定し、
流体運動の解析領域とする前記タンクの周囲の空間について、ｘ軸方向の或る位置に設定する爆風源における地面上に位置する点と、前記タンクの底面の中心とを結ぶ直線を軸とし、該軸に直交する半径方向距離と、該軸を中心とする周方向角座標と、ｘ座標とによる円筒座標
【数１】
2015022599.tif 000044
を設定して、前記流体運動の解析領域は、前記ｘ軸方向における前記タンクが存在する区間で、前記タンクの外面に対応する位置に内半径ｒｍｉｎを有するものとして、前記円筒座標により
【数２】
2015022599.tif 000045
と表し、
該流体運動の解析領域における流体の流速成分と動圧を、下記のように表示し、
【数３】
2015022599.tif 000046
前記式（３）の座標に依存する関数を、下記のように与え、
【数４】
2015022599.tif 000047
次に、流体の支配方程式を、以下のように変分表示し、
【数５】
2015022599.tif 000048
連続条件は、爆風源のないときの次式によって与え、
【数６】
2015022599.tif 000049
次いで、静的密度の座標微分を省略することにより、次の式を得るようにし、
【数７】
2015022599.tif 000050
更に、断熱の関係式を用いて密度変動を動圧変動で表し、変分表示して、爆風源としてＶｘ．ｓｒｃ（ｔ）を導入して次式を得るようにし、
【数８】
2015022599.tif 000051
その後、前記式（３）を、式（７）、式（８）、式（９）、式（１２）に代入し、変分δａ１ｋｎｌ−δａ４ｋｎｌの係数を０と置くことにより、一般化座標に関する時間の常微分方程式を、次のマトリックス方程式の形に導き、
【数９】
2015022599.tif 000052
前記式（１３）を解くことにより、流体運動の解析領域における流体の運動を定めて、前記タンクの外面に対応する位置に作用する流体の圧力を求めるようにする
ことを特徴とする爆風によりタンクに作用する圧力の予測方法。
【請求項２】
請求項１の爆風によりタンクに作用する圧力の予測方法により導かれた一般化座標に関する時間の常微分方程式
【数１０】
2015022599.tif 000053
と、該常微分方程式を導出する過程で用いた流体運動の解析領域における流体の流速成分と動圧を表示した式、
【数１１】
2015022599.tif 000054
を導入し、
タンクのシェルの解析に関して、前記流体運動の解析領域における流体の流速成分と動圧を求める際に設定した円筒座標
【数１２】
2015022599.tif 000055
とは別の円筒座標を以下のように設定し、
【数１３】
2015022599.tif 000056
更に、タンクのシェルについて、以下の有限要素を設定し、
【数１４】
2015022599.tif 000057
次いで、前記タンクのシェルのラグランジュアン汎関数の変分を、前記式（１３）を基に、次のように表し、
【数１５】
2015022599.tif 000058
前記式（１４）における固有値問題
【数１６】
2015022599.tif 000059
を解くことによって、次の変数変換を行い、
【数１７】
2015022599.tif 000060
前記式（１６）及び前記式（３）を、前記式（１４）に代入して、Ｅｍｑｘ（ｔ）、Ｅｍｑｙ（ｔ）に関する常微分方程式を次の形に得て、
【数１８】
2015022599.tif 000061
前記常微分方程式（１７）と（１８）を解くことにより、爆風の動圧をタンクのシェルが受けるときに該タンクのシェルに生じる変位と応力を求めて、該タンクに生じる振動を予測するようにする
ことを特徴とする爆風によるタンク振動の予測方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、タンクの近くで爆発が生じた際に、爆風によりタンクに作用する圧力を予測する方法と、該タンクに作用する圧力の予測方法を基にして、前記タンクに生じる振動を予測するための爆風によるタンク振動の予測方法に関するものである。
【背景技術】
【０００２】
大型のタンクに関する振動予測は、工学における重要課題であり、これまでは、内部に貯留された液体によって生じる振動に関し、多くの研究が行われてきている。
【０００３】
ところで、近年では、大型のタンクの近くで万一爆発が起きた際に、該タンクに作用する爆風の圧力や、爆風により該タンクに生じる振動の問題についての関心が高まっている。特に、紛争が懸念される地域に大型のタンクを設置する際には、前記のような問題が非常に重要になることがある。
【０００４】
しかし、タンク周りの流体運動は、一般に、解析的に表すことが難しいという特徴があり、そのため、従来は、タンクに作用する爆風動圧の評価は、数値流体力学（ＣＦＤ）のような数値的解析方法によって行われてきている（たとえば、非特許文献１参照）。
【先行技術文献】
【非特許文献】
【０００５】
【非特許文献１】島村和夫、他４名，「円筒タンクの耐爆風設計への数値シミュレーションの適用事例」，２００５年６月，石川島播磨技報，第４５巻，第２号，ｐ．８０−８５
【発明の概要】
【発明が解決しようとする課題】
【０００６】
ところが、前記非特許文献１に示された手法は、タンクに作用する爆風動圧の評価として非常に有効であるが、自由度が高い数値的解析を行う必要がある。
【０００７】
そこで、本発明は、解析を行う際の自由度を従来に比して低減させることができる爆風によりタンクに作用する圧力の予測方法、及び、該圧力の予測方法を基に、前記タンクに生じる振動を予測するための爆風によるタンク振動の予測方法を提供しようとするものである。
【課題を解決するための手段】
【０００８】
本発明は、前記課題を解決するために、請求項１に対応して、コンピュータにより、爆風によってタンクに作用する圧力を予測する爆風によりタンクに作用する圧力の予測方法であって、ｘｙｚ座標空間内に、ｚ軸について軸対称な円筒状のタンクの解析モデルを設定し、流体運動の解析領域とする前記タンクの周囲の空間について、ｘ軸方向の或る位置に設定する爆風源における地面上に位置する点と、前記タンクの底面の中心とを結ぶ直線を軸とし、該軸に直交する半径方向距離と、該軸を中心とする周方向角座標と、ｘ座標とによる円筒座標
【数１】
2015022599.tif 000003
を設定して、前記流体運動の解析領域は、前記ｘ軸方向における前記タンクが存在する区間で、前記タンクの外面に対応する位置に内半径ｒｍｉｎを有するものとして、前記円筒座標により
【数２】
2015022599.tif 000004
と表し、該流体運動の解析領域における流体の流速成分と動圧を、下記のように表示し、
【数３】
2015022599.tif 000005
前記式（３）の座標に依存する関数を、下記のように与え、
【数４】
2015022599.tif 000006
次に、流体の支配方程式を、以下のように変分表示し、
【数５】
2015022599.tif 000007
連続条件は、爆風源のないときの次式によって与え、
【数６】
2015022599.tif 000008
次いで、静的密度の座標微分を省略することにより、次の式を得るようにし、
【数７】
2015022599.tif 000009
更に、断熱の関係式を用いて密度変動を動圧変動で表し、変分表示して、爆風源としてＶｘ．ｓｒｃ（ｔ）を導入して次式を得るようにし、
【数８】
2015022599.tif 000010
その後、前記式（３）を、式（７）、式（８）、式（９）、式（１２）に代入し、変分δａ１ｋｎｌ−δａ４ｋｎｌの係数を０と置くことにより、一般化座標に関する時間の常微分方程式を、次のマトリックス方程式の形に導き、
【数９】
2015022599.tif 000011
前記式（１３）を解くことにより、流体運動の解析領域における流体の運動を定めて、前記タンクの外面に対応する位置に作用する流体の圧力を求めるようにする、爆風によりタンクに作用する圧力の予測方法とする。
【０００９】
又、請求項２に対応して、請求項１の爆風によりタンクに作用する圧力の予測方法により導かれた一般化座標に関する時間の常微分方程式
【数１０】
2015022599.tif 000012
と、該常微分方程式を導出する過程で用いた流体運動の解析領域における流体の流速成分と動圧を表示した式、
【数１１】
2015022599.tif 000013
を導入し、タンクのシェルの解析に関して、前記流体運動の解析領域における流体の流速成分と動圧を求める際に設定した円筒座標
【数１２】
2015022599.tif 000014
とは別の円筒座標を以下のように設定し、
【数１３】
2015022599.tif 000015
更に、タンクのシェルについて、以下の有限要素を設定し、
【数１４】
2015022599.tif 000016
次いで、前記タンクのシェルのラグランジュアン汎関数の変分を、前記式（１３）を基に、次のように表し、
【数１５】
2015022599.tif 000017
前記式（１４）における固有値問題
【数１６】
2015022599.tif 000018
を解くことによって、次の変数変換を行い、
【数１７】
2015022599.tif 000019
前記式（１６）及び前記式（３）を、前記式（１４）に代入して、Ｅｍｑｘ（ｔ）、Ｅｍｑｙ（ｔ）に関する常微分方程式を次の形に得て、
【数１８】
2015022599.tif 000020
前記常微分方程式（１７）と（１８）を解くことにより、爆風の動圧をタンクのシェルが受けるときに該タンクのシェルに生じる変位と応力を求めて、該タンクに生じる振動を予測するようにする、爆風によるタンク振動の予測方法とする。
【発明の効果】
【００１０】
本発明によれば、以下のような優れた効果を発揮する。
（１）請求項１の爆風によりタンクに作用する圧力の予測方法によれば、流体運動の解析領域について、新たな円筒座標を設定することにより、流体運動の支配方程式を、該円筒座標で表される流体運動の解析領域について、内半径がｘ軸方向にて部分的に変化する問題に変換することができる。したがって、前記流体運動の支配方程式は、時間に関する常微分方程式に帰着させることができ、且つ該常微分方程式は、１００次元程度と、従来のＣＦＤのような数値解析的手法に比して大幅に低次元となるようにさせることができる。
（２）又、請求項２の爆風によるタンク振動の予測方法によれば、前記（１）の効果を得ることができる爆風によりタンクに作用する圧力の予測方法によって時間に関する常微分方程式に帰着させた流体運動の支配方程式を基にして、タンクのシェルに生じる振動を求めることができる。
【図面の簡単な説明】
【００１１】
【図１】本発明の爆風によりタンクに作用する圧力の予測方法及び爆風によるタンク振動の予測方法の実施の一形態を示すもので、（ａ）はｘｙｚ座標空間におけるタンクの解析モデルと流体運動の解析領域を示す図、（ｂ）は前記流体運動の解析領域に設定するｘ軸を中心とする円筒座標を示す図である。
【図２】タンクのシェルの運動解析に用いるための別の円筒座標を示す図である。
【図３】本発明の有効性について検証した数値例題で用いた爆風源での動圧変動の時刻歴を示す図である。
【図４】本発明についての数値計算の結果として、タンク屋根上の無次元圧力について、ｙ＝０の位置に沿う分布を示す図である。
【図５】非特許文献１に示されたＣＦＤによる図４に対応する解析結果を示すものである。
【図６】タンク屋根に作用する圧力の解析結果を基に求めたタンク屋根の変位と応力の結果を示すもので、（ａ）はタンク屋根上のｙ＝０の位置における法線方向変位、（ｂ）は子午線方向の膜応力、（ｃ）は周方向の膜応力、の分布をそれぞれ示す図である。
【発明を実施するための形態】
【００１２】
以下、本発明を実施するための形態を図面を参照して説明する。
【００１３】
図１（ａ）（ｂ）乃至図６（ａ）（ｂ）（ｃ）は本発明の爆風によりタンクに作用する圧力の予測方法、及び、爆風によるタンク振動の予測方法の実施の一形態を示すものである。
【００１４】
本発明の爆風によりタンクに作用する圧力の予測方法（以下、本発明の圧力予測方法と略称する）では、先ず、図１（ａ）に示すように、ｘｙｚ座標空間内に、解析対象とするタンク１をモデル化する。該タンク１は、球面屋根を有する円筒タンクであり、ｚ軸について軸対称となるようにしてある。
【００１５】
図１におけるａはタンク円筒部（胴部）の半径、Ｈ１はタンク円筒部２の高さ、Ｈは球面状のタンク屋根３の頂上のｚ座標である。
【００１６】
又、後述するように前記非特許文献１に示された数値的解析方法による結果との比較検証のために、ｘ＝ｘｍｉｎの位置に、ｘ軸に垂直な面状の爆風源（面爆源）があるものとする。図１（ａ）における白抜きの矢印は、前記爆風源で発生した爆風の圧力波の伝播方向を示している。
【００１７】
ｚｘ平面に関する問題の対称性により、流体運動の解析は次の領域で行う。
【数１９】
2015022599.tif 000021
【００１８】
ところで、前述したように、従来、タンク周りの流体運動は、一般に、解析的に表すことが難しいとされている。
【００１９】
そのために、本発明の圧力予測方法では、曲線座標の使い方を工夫し、固有値に相当するパラメータの座標依存性を考慮してガレルキン法（モード展開法）を用いることによって、数値的方法に比べ自由度を大きく低減した解析モデルを導くようにする。
【００２０】
そのために、本発明の圧力予測方法では、流体運動の解析領域として、前記面状の爆風源における地面上に位置する点（前記ｘｙｚ座標空間では、（ｘｍｉｎ，０，０））と、前記タンク１の底面の中心（前記ｘｙｚ座標空間では、（０，０，０））とを結ぶ直線（ｘ軸）を軸とする円筒座標を導入する。
【００２１】
これにより、図１（ｂ）に示すように、前記タンク周りの流体運動の解析領域は、ｘ軸からの距離（ｘ軸に直交する半径方向距離）、ｘ軸周りのｘｙ平面から測った周方向角座標、ｘ座標としての
【数２０】
2015022599.tif 000022
を用いて次のように表される。
【数２１】
2015022599.tif 000023
【００２２】
ここで、ｒｍａｘは、前記流体運動の解析領域の外半径であり、該解析領域の外縁の位置を示している。又、ｒｍｉｎは、前記流体運動の解析領域の内半径であり、該内半径ｒｍｉｎのｘ依存性は、タンク１の存在によって生じる。すなわち、前記ｘ軸方向において前記タンク１が存在する個所（ｘの絶対値がａ以下となる範囲）では、該内半径ｒｍｉｎは、前記タンク１の外面の位置を示している。したがって、前記ｒｍａｘとｒｍｉｎの関数形は幾何学的計算によって決めることができることは明らかである。なお、前記ｘの絶対値がａよりも大きい範囲では、前記内半径ｒｍｉｎは、タンク１のサイズに比べて非常に小さい値に設定するものとする。
【００２３】
ところで、前記タンク１の存在によってｘ軸方向に非一様形状となる領域に解析的手法を適用する際に遭遇する問題は、固有値に相当するパラメータが一定でなく、周方向及び軸方向座標に依存することである。これは、従来、一定固有値が解析的手法の適用可能条件と考えられていることによる。
【００２４】
そこで、本発明の圧力予測方法では、前記の適用可能条件の限界を解くために、固有値の周方向及び軸方向座標に関する微分を考慮したガレルキン法を実行するようにする。
【００２５】
前記の円筒座標を用いることにより、解（流体の流速成分と動圧）を下記のように表示する。
【数２２】
2015022599.tif 000024
【００２６】
前記式（３）の座標に依存する関数は、下記のように与える。
【数２３】
2015022599.tif 000025
この際、下記の２点に注意する。
【数２４】
2015022599.tif 000026
【００２７】
次に、流体の支配方程式（基礎式）を、以下のように変分表示（仮想仕事表示）する。
【数２５】
2015022599.tif 000027
【００２８】
連続条件は、爆風源のないときに次式によって与えられる。
【数２６】
2015022599.tif 000028
【００２９】
次いで、静的密度の座標微分を省略することにより、次のようになる。
【数２７】
2015022599.tif 000029
【００３０】
更に、断熱の関係式を用いて密度変動を動圧変動で表し、変分表示して爆風源（ソース）としてＶｘ．ｓｒｃ（ｔ）を導入すると、次式が得られる。
【数２８】
2015022599.tif 000030
【００３１】
前記式（３）を、式（７）、式（８）、式（９）、式（１２）に代入し、変分
δａ１ｋｎｌ−δａ４ｋｎｌ
の係数を０と置くと、一般化座標に関する時間の常微分方程式が、次のマトリックス方程式の形に導かれる。
【数２９】
2015022599.tif 000031
【００３２】
よって、前記式（１３）を解くことによって、前記流体運動の解析領域における流体運動が定められる。これにより、該流体運動を基に、前記タンク１の外面の位置に作用する圧力（予測値）が求まる。
【００３３】
このように、本発明の圧力予測方法によれば、流体運動の解析領域について、前記円筒座標
【数３０】
2015022599.tif 000032
を新たに設定するようにしてあるので、流体運動の支配方程式を、該円筒座標で表される流体運動の解析領域について、内半径ｒｍｉｎがｘ軸方向にて部分的に変化する問題に変換することができる。
【００３４】
したがって、本発明の圧力予測方法によれば、前記流体運動の支配方程式を、時間に関する常微分方程式に帰着させることができ、且つ該常微分方程式は、１００次元程度と、従来のＣＦＤのような数値解析的手法に比して大幅に低次元となるようにすることができる。
【００３５】
続いて、本発明の爆風によるタンク振動の予測方法について説明する。
【００３６】
本発明の爆風によるタンク振動の予測方法（以下、本発明のタンク振動予測方法と略称する）では、先ず、タンク１のシェルの解析に関して、前記流体運動の解析領域に設定した円筒座標とは別の以下の円筒座標を、図２に示すように設定する。
【数３１】
2015022599.tif 000033
【００３７】
更に、タンク１のシェルについて、図２に示すように、以下の有限要素を設定する。
【数３２】
2015022599.tif 000034
【００３８】
これにより、タンク１のシェルのラグランジュアン汎関数の変分（運動方程式に対応）は、前記式（１３）を基にして、次のように書ける。
【数３３】
2015022599.tif 000035
又、前記式（１４）の左辺第３項は、前記式（１３）による流体運動解析によって求めた圧力のなす仮想仕事を表している。
【００３９】
前記式（１４）における固有値問題
【数３４】
2015022599.tif 000036
を解くことによって、次の変数変換を行う。
【数３５】
2015022599.tif 000037
【００４０】
前記式（１６）及び式（３）を、式（１４）に代入して、Ｅｍｑｘ（ｔ）、Ｅｍｑｙ（ｔ）に関する常微分方程式を次の形に得る。
【数３６】
2015022599.tif 000038
【００４１】
したがって、本発明のタンク振動予測方法によれば、前記常微分方程式（１７）と（１８）を解くことにより、前記本発明の圧力予測方法によって求めた爆風の動圧をタンク１のシェルが受けるときに、該タンク１のシェルに生じる変位と応力を求めることができる。このため、前記変位と応力から、前記タンク１に生じる振動を予測することができる。
【００４２】
次いで、本発明の圧力予測方法、及び、本発明のタンク振動予測方法の有効性について数値例題により検証した結果を示す。
【００４３】
数値計算諸元は以下の表のとおりとした。
【表１】
2015022599.tif 000039
【００４４】
又、タンク１のシェルにおけるＺ１≦Ｚ≦Ｚ２でのパラメータは、以下の表のとおりとした。
【表２】
2015022599.tif 000040
【００４５】
爆風源での動圧変動の時刻歴が次式及び図３で示される場合について計算を行った。
【数３７】
2015022599.tif 000041
ここで、ｐｍａｘ＝１６．６ｋＰａ、ｔｍａｘ＝０．１５８５ｓである。規定モード関数は、ｋ＝２、ｎ＝２，ｌ＝１３，ｍ＝３，ｑ＝２０までを考慮した。なお、ここに示すものは計算の一例であり、特定機種を示すものではない。
【００４６】
図４は、数値計算の結果として、タンク屋根３上の無次元圧力について、ｙ＝０の位置に沿う分布を示すものである。各曲線の数値は、タンク円筒部（側壁）２における爆風源に最も近い個所に爆風の圧力波が到達する時間、
ｔ＝（ａ−ｘｍｉｎ）／（ｐ０γ／ρ０）１／２＝０．１１９ｓ
からの経過時間を示すものである。
【００４７】
図５は、比較として、非特許文献１に示されたＣＦＤによる同様の条件の下での解析結果を示すものである。
【００４８】
前記本発明に基づく図４の結果は、図５の結果とよく一致しており、次のことが分かる。
【００４９】
（Ａ）図４、図５共に、最大圧力の生じる位置が時間の経過と共に下流方向へ移動していく様子がシミュレートされている。
【００５０】
（Ｂ）最大圧力の生じる位置が下流方向に移動する速度は、図４と図５で概ね一致している。
【００５１】
（Ｃ）図４、図５共に、上流側で無次元圧力が１より大きくなり得る傾向、すなわち、爆風源よりも大きい動圧がタンク１に作用し得る傾向が示されている。
【００５２】
（Ｄ）図４と図５では、前記（Ｃ）の傾向の強さも概ね一致し、又、下流に移るに従って圧力が減少する度合いも、図４と図５では近似的に等しくなっている。
【００５３】
又、前記タンク屋根３に作用する圧力の解析結果を基に、該タンク屋根３の変位と応力について調べた。これは、一般に、タンク１では、タンク円筒部（側壁）２に比して、タンク屋根３の方が柔軟であるため、該タンク屋根３が爆風による圧力変動の影響を受け易いためである。
【００５４】
図６（ａ）（ｂ）（ｃ）は、タンク屋根３上のｙ＝０の位置について、法線方向変位、子午線方向（ｙ＝０における鉛直面とタンク屋根３との交線に沿う方向）の膜応力、周方向の膜応力の分布をそれぞれ示すものである。なお、膜応力は、表側と裏側の応力を加えて２で割ることによって計算するものとする。
【００５５】
図６（ａ）（ｂ）（ｃ）の結果から、下記の点が確認される。
【００５６】
（Ｅ）最大圧力の生じる位置が前記図４に示したように時間の経過と共に下流方向へ移動することに伴い、変位及び応力が最大となる位置も、時間の経過と共に下流方向に移動する。
【００５７】
（Ｆ）図６（ｃ）に示した周方向の膜応力は、図６（ｂ）に示す子午線方向の膜応力よりも大きい。
【００５８】
更に、曲げ応力は、シェル（曲面殻）の場合、端部のみ集中して大きくなり、図６（ｂ）（ｃ）のように示すと、両端ｘ／ａ＝±１の位置でのみ大きくなる。そこで、以下の表に、タンク屋根３上のｙ＝０の位置における子午線方向の曲げ応力を示す。
【表３】
2015022599.tif 000042
【００５９】
同様に、以下の表は、タンク屋根３上のｙ＝０の位置における周方向の曲げ応力を示すものである。
【表４】
2015022599.tif 000043
【００６０】
前記各表に示された子午線方向と周方向の各曲げ応力の結果からは、次のことが分る。
【００６１】
（Ｇ）初期には上流側の端ｘ／ａ＝−１における応力が大きく、時間と共に下流側の端ｘ／ａ＝１における応力が大きくなる。これは、最大圧力の生じる位置が前記図４に示したように、時間の経過と共に下流方向へ移動することを良く反映している。  
【００６２】
（Ｈ）子午線方向の曲げ応力の方が、周方向の曲げ応力よりも大きくなっている。
【００６３】
前記数値計算では、流体に関するモード次数は、ｋ＝２，ｎ＝２，ｌ＝１３までを考慮している。長手方向（ｘ軸方向）のモード次数ｌは、図４及び図６（ａ）（ｂ）（ｃ）に示すように、爆風の圧力波が到達していない初期の下流領域での圧力が非常に小さいことをシミュレートする上で有効である。
【００６４】
又、タンクの振動に関するモード次数は、ｍ＝３，ｑ＝２０まで考慮している。時間の１階常微分方程式系で表した場合の状態変数の個数は、流体に関し２０８、全体で３０８であり、自由度に応じた計算時間をＣＦＤに比して大きく低減するモード展開法の適用が可能となった。なお、前記数値例題についての計算時間は、現在の標準的レベルのパーソナルコンピュータで２分程度にすることができた。
【００６５】
なお、本発明は前記実施の形態のみに限定されるものではなく、本発明を適用するタンク１は、タンク円筒部２と球面状のタンク屋根３を備えていれば、前記各表に示した数値計算諸元及びシェルのパラメータ以外のいかなる諸元を有するタンク１であってもよい。
【００６６】
又、本発明に従って爆風によりタンク１に作用する圧力の予測を行うため、又は、爆風によるタンク振動の予測を行う数値計算を実施する際、タンク１に作用させる爆風のモデルとしては、爆風源での動圧変動の時刻歴が図３及び式（１９）に示した以外のものであってもよい。
【００６７】
本発明に従って爆風によりタンク１に作用する圧力の予測を行うため、又は、爆風によるタンク振動の予測を行う数値計算を実施する際の前記各モード次数ｋ，ｎ，ｌ，ｍ，ｑは、任意に設定してよい。
【００６８】
その他本発明の要旨を逸脱しない範囲内で種々変更を加え得ることは勿論である。
【符号の説明】
【００６９】
１  タンク、２  タンク円筒部、３  タンク屋根
【図１】
2015022599.tif 000062
【図２】
2015022599.tif 000063
【図３】
2015022599.tif 000064
【図４】
2015022599.tif 000065
【図５】
2015022599.tif 000066
【図６】
2015022599.tif 000067