第 24 回数値流体力学シンポジウム
B8-1

再使用観測ロケット転回時の数値シミュレーション
Numerical flow simulation of a reusable sounding rocket during turn over
○ 葛生和人, JAXA/JEDI, 神奈川県相模原市中央区由野台 3-1-1, E-mail:kuzuu.kazuto@jaxa.jp
北村圭一, JAXA/JEDI, 神奈川県相模原市中央区由野台 3-1-1, E-mail:kitamura.keiichi@jaxa.jp
藤本圭一郎, JAXA/JEDI, 茨城県つくば市千現 2-1-1, E-mail:fujimoto.keiichiro@jaxa.jp
嶋英志, JAXA/JEDI, 神奈川県相模原市中央区由野台 3-1-1, E-mail:shima.eiji@jaxa.jp
Kazuto Kuzuu, Japan Aerospace Exploration Agency, 3-1-1 Yoshinodai, Sagamihara-shi, Kanagawa
Keiichi Kitamura, Japan Aerospace Exploration Agency, 3-1-1 Yoshinodai, Sagamihara-shi, Kanagawa
Keiichiro Fujimoto, Japan Aerospace Exploration Agency, 2-1-1 Sengen, Tsukuba-shi, Ibaraki
Eiji Shima, Japan Aerospace Exploration Agency, 3-1-1 Yoshinodai, Sagamihara-shi, Kanagawa
Flow around a reusable sounding rocket is simulated using unstructured compressible CFD solver, LS-FLOW. Unlike
landing of the space shuttle, this rocket must turn over from gliding position to vertical landing position. It is expectedthat this motion greatly influences the aerodynamic characteristics of the vehicle. In this study, we focus on theanalysis of dynamic characteristics during turn over. An employed numerical code is based on a cell-centered finite
volume compressible flow solver applied to a moving grid system. With respect to a turbulence model, we employed
DDES model in order to investigate massive separated flow. Furthermore, we also carried out the analysis of static
characteristics for comparison. Through this simulation, we discuss the effect of turn over motion on aerodynamiccharacteristics.

１．はじめに垂直離着陸型再使用観測ロケットは、大幅な輸送コスト低減を目指した将来型宇宙輸送システムの 1 つとして期待されている。
JAXA 宇宙科学研究所（ISAS）は、本ロケットの技術的課題、システム要求の抽出などを目的として、これまで、小型再使用ロケット実験機（RVT）の飛行実験を行ってきた。
このような垂直離着陸型ロケットシステムを構築するには、解決すべき多くの工学的課題が存在する。空気力学的課題もそのような検討項目の 1 であろう。たとえば、本ロケットは、地球帰還時、最終的な着陸態勢をとるために、飛行姿勢を滑空姿勢から垂直姿勢に変更する必要が生ずる。その時に発生する流れ場は当然大規模剥離を伴う非定常現象となる。空力設計の観点から、あらかじめそのような流れ場を理解し、その時の空力特性を把握しておくことは必要不可欠なことである。
本研究では、そのような再使用観測ロケットの飛行フェーズの1 つと目されるノーズエントリ帰還飛行時の機体転回運動に着目し、数値シミュレーションを実施することとした。解析は、機体に生ずる非定常流れ場の理解と動的空力特性の把握を目指したものである。
本シミュレーションを行うにあたっては、JAXA 情報計算工学センター（JEDI)の CFD 解析ツール、LS-GRID および LS-FLOWを用いることとした。これらは、いずれも当センターが航空機・宇宙輸送機設計用に開発してきた次期基盤 CFD 解析ツールである。ここで、LS-GRID は直交・物体適合ハイブリッド非構造格子自動格子生成ツール、LS-FLOW は任意多面体非構造圧縮性ソルバーである(1) ( 2)。
なお、解析対象となる流れ場は、その飛行環境から非定常大規模剥離流れを伴うと同時に、機体姿勢も時々刻々変化することとなる。LS-FLOW の開発に当たっては、これまで空力設計ツールとしての実用性を目指して、DES/DDES (3) (4) 乱流モデル、移動格子などを導入し、それらの検証を行ってきた。本シミュレーションではそれらの機能を用いて解析を実施する。

２．解析手法
本解析で用いた解析手法について述べる。
解析に使用した非構造圧縮性ソルバーLS-FLOWは、直交・物体適合ハイブリッド非構造格子への対応を前提とした任意多面体
非構造圧縮性ソルバーであり、支配方程式は3次元非定常圧縮性
Navier- Stokes方程式∂∂t

F Q dSQ dV

F Q dS1に基づいている。空間離散化にはセル中心有限体積法を適用している。ここで、Q は保存量ベクトル、Fe(Q), Fv(Q) は非粘性流束ベクトル、粘性流束ベクトルである。
計算は、時間刻みの離散化に3点後退差分近似、時間発展陰解法にLU-SGS法を用いた。その他の数値計算スキームは、再構築用勾配計算法として非粘性項にLSQ法、粘性項にGreen-Gauss法、勾配制限関数にVenkatakrishnan(5)、数値流束スキームは、非粘性項にSLAU(6)、粘性項にWangの方法(7)を採用した。
また、流れ場は大規模剥離を伴う非定常流れとなることを想定し、DDESモデルを適用した。なお、筆者らはDES/DDESモデル検証のため、LS-FLOWを用いて円柱の流れ場解析を行っており、そこでストローハル数など流れ場の非定常性が定量的にも妥当な値を示していることを確かめている(8)。
次に、今回機体転回時の非定常解析に適用されるLS-FLOW移動格子機能について述べる。支配方程式は、従来のソルバーの支配方程式(1)に移動座標系を考慮した項を加えた3次元非定常圧縮性Navier -Stokes方程式∂∂tQ dVF QQ∙dS

F Q dS2である。ここで は面要素dSの移動速度である。ただし、本機能
1

Copyright © 2010 by JSFM

第 24 回数値流体力学シンポジウム
B8-1
結果を示す。設定した迎角は 5.6 度、42.2 度、84.5 度、126.7 度、169.0 度である。なお、これらの迎角は、転回時の流れの可視化スナップショットとの比較のため、機体姿勢が一致するように設定したものである。
Fig.2 に解析ステップ数に対する揚力係数、抵抗係数、モーメント係数の変動履歴を示す。迎角は 42.2 度である。解析結果
は乱流モデル DDES を用いたため非定常的な変動履歴をとる。
Fig.3(a)(b)は、機体ベース部より非定常的に後流側に放出される渦の様子を示したものである。なお、静的特性値は Fig.2 に示したような周期的特性が認められて以降の平均値として求める。

をLS-FLOWへ実装するにあたり、移動格子は並進または回転するものに限定し、格子自体の変形は伴わないものとしている。実装のための離散化および保存量の線形化は以下のとおりである。

D ∆Q

12FQF

∆QQλ ∆Q

∈2F

QFQV
∆t

12

λ dS

∆Q
Q dSdS3∈

D14∈

なお、移動格子を適用した時のLS-FLOW解析の妥当性に関しては、NACA0012振動翼の空力解析を通してその定量的な検証を行っているので文献(9)を参照されたい。

0.9

cl0.8

cd0.7cm0.6３．解析条件
解析対象となる機体形状および解析条件について述べる。
今回の解析は、機体転回という特殊な飛行環境における流れ場の基本構造を捉えること、およびその時の動的空力特性を求め静的特性との違いを把握することを主眼としている。そのため、機体形状としては付属突起物のない単純形状とした。また、計算格子は、LS-GRIDにて作成した直交・物体適合ハイブリッド非構造格子である。計算領域はx-z中心断面を対称面とする半裁モデルとし、格子点数は約160万点である。機体形状および計算格子を
Fig.1に示す。

0.50.4
0.3
0.2
0.1
0
-0.1
050000100000step

150000200000

Fig.2 History of aerodynamic coefficient (AOA=42.2degs)(a) 18000step
Fig. 1 Configuration of body and mesh
流速条件は、転回運動直前の飛行条件として標準大気中80m/sの飛行速度を想定し、一様流マッハ数M∞=0.235とした。また、レイノルズ数はRe=1.0×107としている。転回特性としては、実機ベースで300度/sと150度/sの回転角速度を想定した。これは、一様流速と機体全長に基づく無次元角速度で0.445、0.2225に相当する。
この条件は機体の終端速度でM∞=0.4～0.5、局所流速ではそれ以上の速度となる可能性があるため、勾配制限関数Venkatakrishnanを採用した。なお、転回の範囲は迎角0度から180度までとし、回転中心は機体先端より60%の位置とした。

(b) 19000step
Fig.3 Pressure and eddy viscosity distribution (AOA=42.2degs)(Color:pressure, B/W:eddy viscosity)
Fig.4(a)(b)に各迎角に対する静的空力特性解析値を示す。

４．結果および考察
４．１ 静的空力特性
最初に、迎角を固定したときの静的空力特性と流れ場の解析
2

Copyright © 2010 by JSFM

第 24 回数値流体力学シンポジウム
B8-1
次に静的解析と実験の比較を示す。
空力解析結果はいずれも非定常な CN（横力係数）、CA（軸力係数）の平均値、実験データは ISAS において実施された同形状機体模型の風洞試験結果である。図に示したように静的な空力特性は迎角の小さい範囲では実験結果と比較しても妥当な値を得ている。しかしながら、大迎角（迎角 90 度を中心とする±約 30 度範囲）では機体の横力が大きく異なってくる。これは、90 度付近の迎角での流れの条件は円柱周りの流れの条件に近くなり、今回のような半裁モデルでは大迎角時における剥離位置が的確に捉えられてないものと思われる。Fig.5 に示したように 84.5 度の高迎角解析結果で後流側に渦が止まっている様子から見ても、半裁モデルの対称条件が本来の非定常渦の生成を正しく再現していない可能性が大きい。横力を正しく捉えるためには機体全周の流れの計算が必要となるであろう。

0.50.40.3
0.20.10-0.1-0.2-0.3-0.4

X=20% from nose
X=40% from nose
Fig.5 Velocity vector field and stream lines on the body surface
AOA=84.5degrees

Cm static cal.
Cm exp.

0306090 120degrees

150４．２ 動的空力特性
機体転回時の動的空力特性の解析結果を示す。
Fig.6 は迎角 0 度の転回直後より 180 度まで反転した時の CN（横力係数）、CA（軸力係数）の履歴を表したものである。回転速度 300 度/s、150 度/s の結果に加えて静的解析結果も示している。なお、転回計算は、迎角 0 度の状態で一様流による初期流れを生成した後、指定した角速度で計算格子を回転させはじめている。
静的結果と比較して転回時の特性で特徴的なのは、CN が転回開始直後より急激に増加していることである。これは、一様流に加えて機体の回転速度成分の空気力が加味されるためであり、回転速度が大きいほどその影響は顕著である。同様に、CAに関しても値は回転速度とともに増加する傾向がある。

180(a) Moment coefficients

5CA static cal.
CA exp.
CN static cal.
CN exp434CA fast rot.
CN fast rot.
CA slow rot.
CN slow rot.
CA static
CN static

32120-10306090 120degrees

1501180(b) Axial and Normal force coeffients
Fig.4 Static aerodynamic characteristics

0-10306090120degrees

150180Fig.6 History of axial and normal force coeffients

3Copyright © 2010 by JSFM

第 24 回数値流体力学シンポジウム
B8-1
くの風下側と機体後方の風上側で圧力上昇が認められる。また、機体側面の圧力分布も大きく異なっており、この結果は、横力特性への影響の可能性も示唆している。今回の解析は半裁モデルであるため、両者の横力特性の差を評価するには至っていないが、今後全機モデルの検討も必要となるであろう。

0.5
Cm fast rot.

0.4

Cm slow rot

0.3

Cm static

0.2
0.1
0-0.1
-0.2

(a) α=5.6 degrees

-0.3
-0.4
-0.5
0306090120150180degrees(b) α=42.2 degrees

Fig.7 History of moment coefficients
一方、回転速度 150 度/s と 300 度/s の結果を比べると、回転速度が遅いほど特性履歴の変動が顕著になっている。４．１で示したように、固定迎角の解析では後流の乱れを反映した空力特性の非定常性が確認されているが（Fig.2,3）、この時の振動特性は無次元周期で 4.5～5.0sec（ストローハル数で約 0.7）となる。
これは、低回転速度 150 度/s の迎角 20 度～50 度間に見られる振動周期と一致する。すなわち、固定迎角で既に認められた非定常特性が、転回時の空力特性変動にもそのまま反映しているものと考えられる。
Fig.7 はモーメント特性の変化を表したものである。
回転速度 300 度/s でのモーメント特性は、揚力や抵抗と同様に大きな回転速度の影響を受け、なおかつモーメントは終始回転方向とは逆の負の側に生ずる結果となっている。しかし、回転速度 150 度/s ではその影響が弱まり、静的特性に近づく傾向が見られる。また、CN 履歴で見られた変動履歴中の振動はモーメント履歴にも認められる。
一方、迎角 0 度から 180 度までのモーメント特性履歴の外形に注目すると、回転速度の増加に従って全体として負の側にシフトしているのがわかる。このことは機体転回時のモーメントダンピング特性と捉えることができる。モーメントダンピング
特性は、動特性として把握すべき重要な特性量であり、今後、定量的な評価が必要となるであろう。なお、回転に伴う位相遅れ現象に関しては顕著な傾向は見られなかった。

(c)α=84.5 degrees(d)α=126.7 degrees(e)α=169.0 degrees
Fig. 8 Comparison of pressure distributions(left:Turn Over cal., right:Fixed AOA cal.)

４．３ 流れ可視化
４．２で示した空力特性が実際にどのような流れ場によってもたらされているのかを明らかにするため、機体周りの流れ場を可視化する。
Fig.8 は中心断面および機体表面上の圧力分布である。図は各迎角において転回時の動的解析結果（300 度/s）と固定迎角の静的解析結果を比較したものである。いずれの迎角においても、転回の影響により機体周りおよび表面の圧力分布が大きく異なっている。特に転回初期において、モーメント特性にダンピング効果を与えるような圧力分布の差異が見られる。すなわち、機体先端近

次に、転回時と静止時の各迎角の機体付近の流線と機体表面上の流線により流れの 3 次元構造を調べる。
Fig.9 は各迎角時の速度ベクトル場に基づいて計算された瞬間流線である。ただし、実際の流れ場は非定常であるため、ここで示されている流線はその時の流れの 3 次元構造を示す指標と考える。なお、転回時の流線で一部空間流線から外れたような機体表面上の流線のずれが見られるが、これはノンスリップ条件で回転する機体表面の流れの状況を反映したものである。
4

Copyright © 2010 by JSFM

図から明らかなように大迎角姿勢時に機体後方に剥離渦が生じているが、機体静止時に比べて、転回時の方が後流側の乱れは小さくなっている。すなわち、転回時の流れ場では機体背面に生成される縦渦が比較的安定した状態で保持されているのに対して、固定迎角ではそのような縦渦が崩壊から大きな非定常流れ場へ発達している。
（Fig.9(b)-(d)）また、この傾向は機体周りの圧力分布にも反映されている。
（Fig.8(c)-(d)）元来、後流側で大きな乱れが生ずるのは大迎角流れの特徴であるが、転回運動がそのような非定常流れを抑制する方向に働いているものとしたら、そのことは本機の空力安定性能を見積もる上で看過することのできない重要な問題である。回転速度との関係も含めて、今後はより詳細な解析と定量的な評価が求められるところである。



第 24 回数値流体力学シンポジウム
B8-1
転回による空力特性への影響を、動的空力特性の変動履歴として捉える事が出来た。
静動解析結果の流れの可視化を通して、流れ場の構造や非定常性、機体表面圧力の違いなどを確認することができた。
実験結果との比較では迎角の小さい範囲では定量的な一致は見たものの、大迎角時の特性としては定性的な評価にとどまった。

以上のことを踏まえて、今後は、大迎角時の定量的評価、横力特性の評価のため全機モデルのシミュレーションが必要となるであろう。
謝辞
実験データを提供して頂いた JAXA 宇宙科学研究所（ISAS）助教野中聡氏に感謝する。
参考文献(1) Fujimoto, K., and Fujii, K.,“Study on the Automated CFD Analysis
tools for Conceptual Design of Space Transportation Vehicles,”FEDSM 2007-37128, Proc. 5th ASME/JSME Joint Fluid
EngineeringConference, San Diego, California, July 2007.
(2) Fujimoto, K., and Fujii, K.,“Optimization-based Robust Feature
Preserving Technique for Bodyfitted Cartesian Grid Method,” 10thInternational Conference on Numerical Grid Generation in
Computational Field Simulations (CD-ROM), Crete, Greece,
September, 2007.
(3) Spalart, P. R., Jou, W. H., Strelets, M., and Allmaras, S. R.,“Comments of the Feasibility of LES for Wings, and on a Hybrid
RANS/LES Approach,” Advances in DNS/LES, 1st AFOSR Int.
Conf. on DNS/LES, (1997)(4) Spalart, P. R., Deek, S., Shur, M. L., Squires, K. D., Strelets, M., andTravin, A., “A new version of detached-eddy simulation, resistant toambiguous grid densities,” Theor. Comput. Fluid Dynamics, 20(2006) pp. 181-195.
(5) Venkatakrishnan, "Convergence to Steady State Solutions of theEuler Equations on Unstructured Grids with Limiters," Journal ofComp. Phys., Vol. 118, 1995, pp.120-130.
(6) 嶋, " AUSM 族全速度スキームと前処理付き陰解法," 第 22
回数値流体力学シンポジウム, C7-4 (2008).
(7) Wang, "A Quadtree-based adaptive Cartesian/ Quad grid flow
solver for Navier-Stokes equations," Computers and Fluids, Vol.27,(1998), pp.529-549.
(8) 葛生,北村,藤本,嶋,"直交・物体適合ハイブリッド非構造格子
ソルバーにおける DES,DDES の適用とその評価について,"第 23 回数値流体力学ｼﾝﾎﾟｼﾞｳﾑ, A2-4 (2009).
(9) 葛生,北村,藤本,嶋," 移動格子を適用した任意多面体非構造
圧縮性ソルバーによる振動翼解析," 第 42 回流体力学講演会/ANSS2010, 2E8 (2010).

(a) α=42.2 degrees(b)α=84.5 degrees(c)α=126.7 degrees(d)α=169.0 degrees
Fig. 9 Comparison of instantaneous streamlines(left:Turn Over cal., right:Fixed AOA cal.)
４．結論
任意多面体非構造圧縮性ソルバーLS- FLOW を用いて、再使用観測ロケットの機体転回運動時の数値シミュレーションを実施し、その動的な空力特性、非定常流れ場の解析を行った。
静的特性の解析結果や実験結果との比較を通して、以下のような結論を得た。


機体転回運動によるピッチングモーメントのダンピング特性が数値シミュレーションにより再現され、回転速度によって大きく変動することを確認した。
5Copyright © 2010 by JSFM