船舶・設備

【1Dシミュレーション②】Pythonで作る揚水発電所とロケットエンジンの1Dモデル

Pythonで作る1Dモデル 揚水発電所とロケットエンジンを解く

ロケットエンジンの条件を変えてみる

前のページのモデルを使って、次の7つを順に調べます。以下の数値は、すべてこのモデルでの計算結果です。計算は、どれも1ケース数秒です(最後の通し計算だけ1〜2分)。

  1. 酸燃比がずれたとき(オリフィスや噴射器の出来)
  2. GGの弁を開く順番
  3. 主弁を開く順番と速さ
  4. 起動トルクの強さ
  5. 推力を絞る
  6. エンジンを止める手順
  7. 起動から燃焼終了までを通して解いてみる

酸燃比がずれたとき

Merlinでは、酸燃比を推進剤の配管の太さで決め、バタフライ弁で微調整します(Wikipedia)。つまり、配管やオリフィスの出来が、そのまま酸燃比に効きます。そこで、LOX側のオリフィスや噴射器の流れやすさ(流量係数×面積)が、設計から最大±20%ずれた場合を計算しました。数%の製作公差だけでなく、詰まりや部品の付け間違いまで含めた幅です。

LOX側のオリフィス・噴射器の流れやすさのずれと推力の関係
LOX側の流れやすさが±20%ずれたときの推力

GGのLOXオリフィスが20%広いだけで、推力は845 kNから1,028 kNへ22%増えます。20%狭ければ722 kNへ15%減ります。GGに流れるのは全体の7.1%ですが、そのO/Fがタービンの出力を通じてエンジン全体を左右します。主燃焼室のLOX噴射器を同じだけずらしても、推力は828〜860 kNで、変化は3%以内です。

なお、このモデルでは、タービン出口の圧力を0.2 MPaで一定としています。実機では排気の出口の大きさが決まっているので、タービンを通る流量が増えれば出口の圧力も上がり、タービンの出力の増え方は抑えられます。GGのずれに対する推力の変化は、実際にはこれより小さめになります。

GGのLOXオリフィスのずれとGG出口温度の関係
GGのLOXオリフィスが±20%ずれたときのGG出口温度

GGのLOXオリフィスを±20%ずらしたとき、GGの出口温度(=タービン入口温度)は約800〜1,000 K(O/F 0.26〜0.39)の範囲で動きます。O/F 0.30より下は教科書に値がないので、表の傾向を延ばした点を足して補っています。Huzel & Huangの教科書は、タービン入口温度の実用上の限界を約1,200 K(1,700 °F)としています。±20%のずれなら、まだ余裕があります。

主燃焼室側のずれは、GGにも影響します。主燃焼室のLOX噴射器をずらしただけで、GGのO/Fは0.31〜0.35に動きました。ポンプの吐出圧を通じて全体がつながっているためで、こうした間接的な影響は、システム全体を解かないと見えません。

全体の7%ほどしか流れないGGが、エンジン全体を握ってるんだね。

GGの弁を開く順番

次は起動の順番です。GGの燃料弁は0.15秒のままにして、LOX弁を開くタイミングだけを前後に数ms〜数十msずらします。

GGの燃料弁とLOX弁を開く順番によるGG出口温度の時間変化の違い
GGの弁を開く順番と、GG出口温度の時間変化
GGのLOX弁のタイミングのずれとGG出口温度のピークの関係
GGのLOX弁のタイミングと、GG出口温度のピーク

LOX弁が5 ms早いだけで、GG出口温度のピークは約1,300 Kまで上がり、限界の目安の1,200 Kを超えます。10 msなら約2,500 K、20 msなら約3,300 Kです。燃料が先なら、温度は設計値付近(910 K)のままです。

たった5 msのずれで、400 Kも上がるんだ。。。

LOXが先に入ると、燃料が届いた直後のGGは極端な酸素過剰です。そこから燃料が増えていく途中で、燃焼温度が最高になるO/F 3付近を必ず通過します。温度だけを見れば、燃料過剰で運転するGGは、燃料を先に入れる順番が安全です。

ただし、ピークの値そのものは目安です。起動前からGGに入っているガス(今回は1気圧の窒素と仮定)が、最初にできた燃焼ガスを薄めるので、その量でピークが大きく変わります。10 msのずれの場合、窒素を半分にするとピークは約3,000 K、2倍にすると約1,900 Kです(図の細い線)。5 msのずれでは、窒素が2倍あるとピークは約1,100 Kで、1,200 Kを超えません。

薄まったガスの温度は、ガスの比熱の仮定にも左右されます。また、ピークのときのGGの圧力は約0.4〜0.9 MPaと低く、実際の燃焼温度は、表の値(12.9 MPaでの値)より低くなります(「エンジンを止める手順」で数字を挙げます)。高温が続く時間も短く、基準の窒素の量でも、10 msのずれで1,200 Kを超えているのは約20 msです。数msのずれでも温度が跳ね上がる、という傾向として見てください。

また、このモデルには、着火の遅れや、液滴が蒸発して燃えるまでの時間は入っていません。実機では、着火が遅れて溜まった推進剤が一気に燃える「ハードスタート」もあります。燃料を早く入れすぎるとその原因になるので、最後は燃焼試験で確かめる領域です。

主弁を開く順番と速さ

主燃料弁と主LOX弁でも、同じ計算をしました。

主燃料弁と主LOX弁を開く順番による主燃焼室の酸燃比の時間変化の違い
主弁を開く順番と、主燃焼室の酸燃比の時間変化

主LOX弁が50 ms早いと、主燃焼室のO/Fは酸素過剰の側から始まり、燃焼温度が最高になる範囲を通って設計値に近づきます。酸素過剰の高温ガスは、燃焼室の壁や噴射器を傷める原因になります。LOX弁が遅い場合は、燃料過剰の側から近づきます。

主LOX弁のタイミングのずれと回転数のピークの関係
主LOX弁のタイミングと、回転数のピーク

ところが回転数は逆の傾向です。主LOX弁が遅れるほど、LOXポンプに負荷がかかり始めるのも遅くなります。回転数のピークは、50 msの遅れで設計値の6%増し、100 msで14%増し、150 msでは24%増し(約44,800 rpm)まで上がります。負荷の軽いLOXポンプの吐出圧が上がり、GGに流れ込むLOXが増えて、GGの温度が上がる(150 msの遅れで約1,070 K)ことも、タービンの出力を押し上げています。LOX弁が早ければ主燃焼室が酸素過剰になり、遅ければ回転数が行き過ぎるので、どちらにずれても別の問題が出ます。実機でどちらを先にするかは、このモデルにない着火の遅れや推進剤の溜まりにも左右されるので、温度と回転数だけでは決まりません。

主弁が全開になる時刻と回転数のピークの関係
主弁を開き切る時刻と、回転数のピーク

主弁を開き切る時刻を基準の0.5秒から遅らせると、回転数のピークは0.9秒で設計値の18%増し、1.6秒では33%増し(約47,800 rpm)まで上がります。

主弁をゆっくり開いたときのタービン出力とポンプ動力の時間変化
主弁を0.9秒で開き切る場合の、タービンの出力とポンプの動力(ポンプは燃料とLOXの2台の合計)

原因はタービンとポンプの動力の差です。GGが先に立ち上がってタービンの出力が増える一方、主弁はまだ開き切っていません。ポンプを流れる量がまだ少なく、ポンプが吸収する動力も小さいので、タービンの出力に追いつきません。その差(水色の部分)が軸を加速します。起動トルクは回転数が上がると弱まり、設計回転数で0になるので、設計回転数を超えたあとの加速は、すべてこの差によるものです。しかも回転数が上がるとポンプの吐出圧が上がり、GGに流れ込む推進剤が増えて、タービンの出力がさらに増えます。このとき吐出圧は、設計値の約1.6倍まで上がります。

このモデルでは、主弁を速く開くほどターボポンプには楽になります。ただし速く開けば、配管の中の液体が急に動かされて起きる圧力波(水撃)や、着火時の圧力の急上昇といった、このモデルに入っていない現象が厳しくなります。

速く開けたほうがターボポンプには楽、っていうのは意外だね。ゆっくり開けたほうが優しそうなのに。

ポンプに早く負荷がかかるほど、回りすぎる前にブレーキが効くんです。

起動トルクの強さ

起動の条件の最後は、スピンスタートの強さです。止まっているときの起動トルクを200〜800 N·mで変えました。

起動トルクの強さと回転数のピークの関係
起動トルクの強さと、回転数のピーク

基準の300 N·m以下なら、回転数の行き過ぎは出ません。その代わり、弱いほど立ち上がりは遅くなり、推力が90%に達する時刻は、基準の0.46秒に対して、200 N·mでは0.53秒です。逆に強くすると、400 N·mで設計値の約2%、600 N·mで6%、800 N·mでは9%の行き過ぎが出ます。800 N·mのとき、ポンプの吐出圧は設計値の28%増しです。基準の300 N·mは、調べた範囲では、行き過ぎを出さずに、いちばん速く立ち上がる強さです。

行き過ぎの分を押しているのは、起動トルクそのものではありません。起動トルクが強いと、主弁が開き切る前に回転数が上がり、GGも早く立ち上がるので、ここでもタービンの出力がポンプの動力を上回ります。

起動の手順は、温度・回転数・圧力・立ち上がりの速さのせめぎ合いで決まります。1回数秒で計算できるので、何百通りでも試せるのが1Dシミュレーションの強みです。

推力を絞る

ガスジェネレータサイクルのエンジンは、GGに流す量でタービンの出力が決まり、それがエンジン全体の流量を決めます。GGのオリフィスを、燃料側もLOX側も同じ割合で絞ったときの定常値を計算しました。

GGオリフィスの開き具合と推力の関係
GGへの流量を絞ったときの推力

GGのオリフィスを40%まで絞ると、推力は設計値の46%になります。30%まで絞ると、設計値の33%(276 kN)まで下がりました。絞るほど燃焼圧が下がり、海面ではノズルの膨張が過剰になるので、比推力は282 sから256 sへ落ちます。タービン出口の圧力を一定としているぶん、絞ったときの推力の落ち方は、「酸燃比がずれたとき」で書いたとおり大きめに出ています。

30%まで絞ったときの燃焼圧は約4.8 MPaです。流れが出口まで壁に沿うとすると、出口でのガスの圧力は大気圧の約1/3まで下がります。ノズルの中の圧力が大気圧よりここまで低くなると、流れは出口の手前で壁からはがれます。このモデルでは、壁の圧力が大気圧の0.4倍まで下がった位置で流れがはがれるとし、そこから先の壁には大気圧がかかるとして、推力を計算しています。0.4倍という値は、Wikipediaにある、出口の圧力についての目安です(約30〜45%の幅があります)。このモデルでは、同じ値をノズルの途中の壁の圧力に当てはめました。この扱いに切り替わるのは、燃焼圧が約5.7 MPaより下がったときです。

実機のMerlin 1Dは、推力を最大の40%まで絞れます(Wikipedia)。このモデルでは、オリフィスを30〜40%まで絞ったあたりに当たります。

この絞りが生きるのが、Falcon 9の1段目の着陸です。最後は中央の1基だけを、推力を絞って噴射します(3基で減速を始めて、途中で外側の2基を止めることもあります)。推進剤をほとんど使い切った機体は軽く、いちばん絞っても推力が重さを上回るので、ホバリングはできません。そこで、着地の瞬間に速度がちょうどゼロになるように、噴射を始めるタイミングと推力を調整しています。

エンジンを止める手順

停止(シャットダウン)には、起動以上に気を使う面があります。起動では低い圧力と回転数から立ち上げていきますが、停止では、ターボポンプが全速で回り、燃焼室も高圧という状態から、短い時間でエンジンを止めなければなりません。

Huzel & Huangの教科書は、先にGGを止めてから主燃焼室を止めること、燃料過剰のまま止まるように弁を閉じることを、停止の手順の基本としています。解析で調べる項目としては、弁を閉じる順番とタイミング、配管の圧力の急上昇、主燃焼室やGGの温度の急上昇、停止の指令のあとに残る推力を挙げています。これに沿って、手順を次のように仮定しました(実機の停止手順も公表されていません)。

  • 0秒:停止の指令。GGのLOX弁と主LOX弁を閉じ始める
  • 0.02秒:GGの燃料弁と主燃料弁を閉じ始める(LOX弁より20 ms遅らせて、燃料過剰のまま止める)
  • GGの弁は0.05秒、主弁は0.15秒かけて閉じる(GGの弁と主弁は同時に閉じ始めるので、GGが先に止まる)

停止過渡での推力の時間変化(主弁を閉じる時間を変えた場合)
停止過渡:主弁を閉じる速さと、推力の下がり方

基準の手順では、推力は指令から0.14秒で10%まで下がり、0.15秒でほぼゼロになります。指令のあとに出る推力を時間で積分した力積は約67 kN·sで、停止前の推力(845 kN)の約0.08秒分です(どちらも海面での計算です)。この力積がばらつくと、そのぶん機体の速度もばらつきます。そのため力積は、小さく、しかも毎回同じであることが求められます。なお、燃焼圧が約5.7 MPaより下がってからの推力は、「推力を絞る」で説明した、ノズルの流れが壁からはがれる場合の計算によるものです。

主弁を0.05秒で閉じれば、力積は約28 kN·sまで減ります。ただし、全速で回っているポンプの出口を急にふさぐことになり、吐出圧は設計値の約1.1倍に上がります。実機ではこれに水撃が加わります。逆に0.5秒かけて閉じると、力積は約124 kN·sに増えます。

停止過渡でのターボポンプ回転数の時間変化(GGの弁を閉じるのが遅れた場合)
停止過渡:GGの弁を閉じるのが遅れたときの回転数

次に、GGの弁を閉じ始めるのを遅らせました。主弁が先に閉じて負荷を失ったポンプを、タービンが回し続けるので、回転数が跳ね上がります。ピークは、100 msの遅れで設計値の12%増し、150 msの遅れで29%増しです。150 ms遅れた場合、ポンプの吐出圧は設計値の約2.1倍になります。基準の手順では、回転数は指令のあと下がる一方ですが、0.2秒後でも、まだ設計回転数の半分で回っています。

GGのLOX弁と燃料弁を閉じる順番とGG出口温度のピークの関係
GGの弁を閉じる順番と、GG出口温度のピーク

GGの弁を閉じる順番も、温度に大きく効きます。LOX弁が先に閉じれば、GGの温度は下がるだけです。燃料弁が先に閉じると、GGの温度は5 msのずれで約1,420 Kに達し、目安の1,200 Kを超えます。10 msのずれでは約3,790 K、20 ms以上では燃焼温度の表の最高値(約3,800 K)まで上がります。起動のときは、同じ10 msのずれで約2,500 Kでした。燃料弁が閉じたあとも流れ込むLOXが、GGに残っている燃料過剰のガスと燃えて、中身のO/Fが、燃焼温度が最高になる3付近まで上がるためです。

起動のときと違って、停止のときのピークは、起動前のガスの仮定には左右されません。薄める窒素はもう残っておらず、GGはすでに約900 Kあります。しかもポンプの吐出圧が高いままなので、弁が閉じ切る直前まで、推進剤が勢いよく流れ込みます。

ただし、この温度は、流れ込んだLOXがすぐに混ざって燃え、壁へ熱が逃げないとした場合の値です。燃焼温度の表も、燃焼圧12.9 MPaでの値です。温度がピークになるとき、GGの圧力はそれより低く、10 msのずれでは約0.8 MPa、40 msのずれでも約9.0 MPaです。化学平衡計算では、O/F 3の燃焼温度は圧力が低いほど下がり、4 MPaで約3,620 K、1 MPaで約3,410 Kです(表の値より約180〜390 K低い値)。10 msのずれでは圧力が1 MPaより低いので、実際のピークは、ここに書いた値よりさらに低くなります。

主弁でも順番が効きます。基準ではLOX弁を先に閉じるので、主燃焼室は燃料過剰のまま推力が下がっていきます。燃料弁が20 ms先に閉じると、止まる直前の主燃焼室は酸素過剰になります。

止めるだけだと思ってたけど、弁の順番は起動より厳しいんだね。

しかも、ここまではタンクの圧力が一定という、いちばん簡単な条件での話です。実機では、ここにさらに条件が加わります。

実機でさらに難しくなる点

このモデルに入っていないものもあります。弁を閉じたあとも、弁から噴射器までの配管や冷却通路に残った推進剤は流れ続けるので、実際の推力はもう少し尾を引きます。水撃も、配管の慣性を入れないと出てきません。

ポンプ入口の圧力も、実機では一定ではありません。実際のポンプ入口の圧力は、タンクの圧力(今回は0.3 MPaで一定と仮定)に、液面からポンプまでの液柱の圧力を足し、配管の損失と、流れの速さのぶん(動圧)を引いたものです。液柱の圧力は、液柱の高さと、機体にかかる加速度(地上で静止していれば1G)に比例します。たとえばLOXの液柱が10 mあれば、1Gで約0.11 MPa、3Gなら約0.34 MPaで、仮定したタンク圧と同じくらいの大きさです。

燃焼が進むと、タンクの残量が減って液面が下がる一方、機体は軽くなって加速度が上がります。機体が向きを変えれば、タンクの中の推進剤が揺れます(スロッシング)。エンジンを止めれば機体を押す力がなくなり、推進剤はタンクの底に押し付けられなくなります。このようにポンプ入口の圧力は時々刻々と変わり、下がりすぎると、ポンプの入口で液体が沸騰して泡ができ(キャビテーション)、揚程が落ちます。起動も停止も、この変動があっても確実に行えなければならないので、実機の手順づくりは、今回の計算よりさらに難しくなります。

起動から燃焼終了までを通して解いてみる

前編で、ロケットエンジンには「起動の1秒」と「燃焼の150秒」が同居していると書きました。同じモデルで、起動から推進剤が尽きるまでを通して計算します。タンクには150秒分の推進剤(LOX 30.7 t、燃料15.1 t)を積み、主燃焼室のLOX噴射器が設計より2%流れにくい、という小さなずれを入れました。

起動から燃焼終了までの推力の時間変化(対数時間軸)
起動から燃焼終了までの推力

時間刻み0.5 msで約30万ステップ、計算時間は1〜2分です。起動の1秒も、149秒後の燃焼終了も、同じ式を同じ時間刻みで解いた結果です。

燃焼中のタンクの推進剤残量の時間変化
タンクの推進剤の残量

噴射器の2%のずれで、エンジン全体のO/F(主燃焼室の2.36とGGの0.32を合わせた値)は2.03から2.01に下がります。その結果、149秒で燃料が先に尽き、LOXが約284 kg残りました。使い残した推進剤は、燃やせない重りとして運んだことになり、そのぶん打ち上げ能力が減ります。噴射器の2%のずれでも、こうして打ち上げ能力に響いてきます。

このモデルでは、どちらかのタンクが空になった時点で、すべての弁を閉じています。タンクの残量を見て、空になる前にエンジンを止めるといった手順は、モデルに入れていません。

ロケットの重さと打ち上げ能力の関係は、こちらの記事で紹介しています。

30万ステップでも1〜2分なんだね。ところで、時間刻みの0.5 msはどうやって決めたの?

いい質問です。実は、速い現象と遅い現象が同居するモデルには落とし穴があって、時間刻みを間違えると計算が発散します。次のページで、その落とし穴にわざと落ちてみます(汗)

※このページの値は、仮定を置いた簡易モデルの計算結果です。前のページの前提の表のうち非公表の値(燃焼圧は推定値)と、起動・停止の手順は仮定で、実機のMerlin 1Dの値や手順を表すものではありません。

次ページへ

次のページでは、計算が発散する例、燃焼振動、公開しているコードの使い方を紹介します。

1 2 3 4
ABOUT ME
Trude(まきブロ管理人)
Trude (まきブロ管理人)
多くの人に機械工学の面白さやエンジニアリングの楽しさを知ってもらうべく、解説や紹介記事を発信しています。

PR(広告を含みます)

オンラインゲームの通信が不安定なときの対策として、通信経路を最適化するソフト「ExitLag」があります。気になる方は公式サイト(紹介リンク)をご覧ください。

RELATED POST

COMMENT

メールアドレスが公開されることはありません。 ※ が付いている欄は必須項目です