— 67 —SENAC Vol. 47, No. 2（2014. 4）[民間企業利用サービス利用成果]

民間航空機開発における大規模空力弾性解析シミュレーションの適用
森野裕行
三菱航空機株式会社

三菱航空機株式会社では、平成 20 年度第 2 期の先端的大規模計算シミュレーションプログラム利用サービスに始まり、平成 23 年度からは民間企業利用サービスの枠組みで、東北大学サイバーサイエンスセンター所有のベクトル計算機 SX-9 を利用している。この大規模計算機を利用することで、数値流体力学（CFD）による空力弾性解析シミュレーションを MRJ 開発で本格的に活用することが可能となった。これにより、フラッタに代表される空力弾性不安定問題に対する設計リスクを最小限に抑え、飛行安全上のリスク低減を図ることができた。本稿では、この空力弾性解析シミュレーションの内容と、MRJ 開発における適用状況を紹介する。

1.

はじめに

三菱航空機株式会社では、平成 20 年度より YS-11 以来半世紀ぶりとなる国産旅客機、三菱リージョナルジェット（MRJ:図 1）の開発を進めている。MRJ の開発では、環境負荷低減のため、同クラスの現行ジェット旅客機の燃費に対して、機体の軽量化・低抵抗化と新エンジンの搭載を含めて 2 割以上の燃費削減を目標としており、これを支える要素技術の開発に東北大学や宇宙航空研究開発機構と連携した産官学の共同研究を積極的に活用してきた。その内容は、空気力学、空力弾性、材料／構造、装備、飛行制御等の各要素技術や、多分野統合最適化／多目的設計探査の研究など、広範囲にわたっている。本稿では、東北大学との共同研究で開発した数値流体力学（CFD）を用いた空力弾性解析コードについて、東北大学所有の SX-9 による解析例と MRJ 開発における適用状況を紹介する。

図 1 三菱リージョナルジェット（MRJ）

— 68 —

2.

SENAC Vol. 47, No. 2（2014. 4）

民間航空機開発における空力弾性設計

航空機の燃費削減を実現するためには、機体の低抵抗化と軽量化が必要となる。航空機の軽量化を進めると機体構造の剛性は低下することになるため、高速飛行時においてフラッタと呼ばれる自励振動現象の発生が懸念される。フラッタの発生は最悪の場合機体構造の破壊につながるため、民間航空機の安全性を担保する型式証明を取得するためには、解析、風洞試験、飛行試験の3 本柱で 15%マージンを含めた運用領域内でフラッタ等の空力弾性上の不安定現象が起こらないことを示す必要がある。
フラッタに関する風洞試験、飛行試験は開発フェーズの最終段階で実施され、その実施条件も限定されることから、開発フェーズの各段階で設計進捗に応じた空力形状・構造特性（マス、剛性分布など）を反映した解析を行い、軽量化（低剛性化）によるフラッタ特性変化を高い精度で確認しておくことが設計上のキーポイントとなる。
航空機設計におけるフラッタ解析については商用ソフト NASTRAN の使用実績が世界的にも多く、当社でも NASTRAN をフラッタ解析の主力ソフトとして利用している。ただし、NASTRAN のフラッタ解析では線形揚力面理論に基づく非定常空力計算を適用しているため、衝撃波や粘性に起因する空気力の非線形特性が支配的となる遷音速条件ではフラッタ解析精度が低下する。解析精度の低下は構造剛性の設計要求に安全マージンを加えることに繋がるため、機体構造に対して重量ペナルティを負う可能性がある。そのため当社では、Euler／Navier-Stokes 方程式を用いた CFD 解析をフラッタなどの空力弾性解析に適用する研究開発を東北大学との共同研究で行い、MRJ の空力弾性設計に活用している。しかしながら、このような CFD を利用した空力弾性解析は、非定常計算で 1 ケースあたり数千～数万ステップの時間進行計算が必要となることから、一民間企業の計算機環境では設計ツールとして積極利用するのが大変厳しいのが現状である。研究目的の解析であれば数ケースで十分であるが、航空機設計における空力弾性解析となると、運用領域を包含する飛行条件及び想定されるあらゆる故障条件を考慮して数百～数千規模の解析ケースが必要となる。これら全てを CFD 空力弾性解析で実行するのは現実的ではないものの、NASTRAN による解析では精度低下が懸念される特定のケースを抽出して CFD 空力弾性解析を適用できれば、フラッタに対する設計リスクを最小限に抑えることができる。当社では、平成 20 年度第 2 期の先端的大規模計算シミュレーションプログラム利用サービスに始まり、平成 23 年度からは民間企業利用サービスの枠組みで東北大学のベクトル計算機 SX-9 を利用している。この大規模計算機を活用することで、計算機環境の問題を解決し、MRJ 開発における CFD 空力弾性解析の本格適用を実現することができた。以降に CFD 空力弾性解析コードの内容とその解析例を紹介する。

3.

空力弾性解析コード

開発した空力弾性解析コードは、東北大学の非構造格子ソルバーTAS（Tohoku University
Aerodynamic Simulation）コード[1-5]をベースとしており、流体側の計算では Euler／
Navier-Stokes 方程式をセル節点有限体積法で解く。計算格子の移動変形に伴う 3 次元圧縮性
Euler 方程式は以下のようになる。

QdV   F Q   U g Q   ndS  0 t t  t (1)Q   , u , v, w, e

(2)T民間航空機開発における大規模空力弾性解析シミュレーションの適用— 69 — u   v  w  u 2  p   uv  uw  F Q    uv i   v 2  p  j   vw  k  2 uw   vw  w  p e  p w e  p u   e  p v (3)ri n 1  ri n

t(4)

Ungi

 1
p    1e   u 2  v 2  w 2  2


(5)

ここで Q は保存量ベクトル、F Q  は非粘性流束ベクトル、U g i は各格子点の移動速度である。
セル境界での流束計算には HLLEW 近似リーマン解法[6]で高次精度化を図り、空間精度を悪化させず、かつ収束性の高い Venkatakrishnan の流束制限関数を用いている。時間積分には LU-SGS 法と3 点後退差分による Newton 反復を組み合わせた 2 次精度の陰解法を用いた。また計算格子の変形には Murayama らの方法[7]を採用した。
構造側の弾性振動方程式は次式で与えられる。
nM d K d   F 

(6)

ここで M は質量マトリクス、 K は剛性マトリクス、 d は変位、 F は空気力である。式(6)にレーリー・リッツ法を適用して運動方程式を構築し、振動変位 d を次式のように幾つかの固有振動モードの重ね合わせで表現するモーダル法により解析を行う。

d   q(7)

ここで  は固有振動モードを表す変換マトリクス、q は一般化変位である。式(6)及び式(7)より、最終的な支配方程式は次のようになる。

0S   2i0 1
S T
0
 i F / M i (8)S  q q (9)ここで  i は i 番目振動モードの固有角振動数である。
流体と構造の連成には Fully Implicit Coupling 法[8]を採用し、時間積分の際に内部反復を行なうことで時間遅れのない完全な時間 2 次精度を確保した。
開発した解析コードはベクトル計算機向けにチューニングされており、効率の良い計算領域分割による MPI 並列化を用いて大規模並列計算を実現している。

4.

空力弾性解析例

当社では、MRJ の開発試験として、主翼と尾翼を対象としたフラッタ風洞試験を実施している。
本章では、風洞試験（WTT）で使用した弾性模型のフラッタ解析を開発コードで実施した解析例を紹介する。また、風洞試験データを利用した開発コードの精度検証についても紹介する。

— 70 —

SENAC Vol. 47, No. 2（2014. 4）

4.1 MRJ 主翼フラッタ風洞試験
風洞試験には、弾性翼、パイロンナセル部、胴体フェアリングで構成されるセミスパン模型を使用し、試験は三菱重工業（株）名古屋航空宇宙システム製作所（名航）の高速風洞（ブローダウン式、計測部 0.6m×0.6m）で実施した。模型の洞内設置状態を図 2 に、模型の 3 面図を図 3 に示す。弾性翼の主構造はアルミ製の梁要素とリブ・ウェブ要素で構成される。主翼の弾性特性は梁要素で模擬し、空力形状は樹脂で整形した。パイロンナセル部はアルミ製パイロンと光造形で製作したフロースルーナセルで構成される。この弾性翼模型は風洞壁に片持ち支持するのではな
く、板要素を有するソフト支持装置に固定する方式とし、機体が空中状態にある支持条件での振動特性を再現可能とした。模型のマス分布と剛性分布は、想定される主翼フラッタモードのうち、最も重要となるエンジン系モードを再現するよう設計した。このエンジン系フラッタモードは主翼の曲げ 1 次モード（WB1）とナセルのピッチングモード（NCP）の連成により発生する。模型の構造解析モデルを NASTRAN で作成し、振動解析を行った。振動解析結果を表 1、図 4 に示す。代
表的な振動モードについて、解析結果は振動試験（GVT）結果と良好に対応しており、作成した構造解析モデルの妥当性が確認された。
風洞試験では、マッハ数を固定して風洞総圧（動圧）を段階的に上昇させ、フラッタが発生するポイントを慎重に調査し、複数のマッハ数条件でフラッタポイントを取得した。図 5 に模型のフラッタ境界について、風洞試験結果と解析結果を比較した例を示す。図中の縦軸・横軸の数値については、公開できない情報のため削除しているが、試験のマッハ数条件は遷音速領域であることに留意されたい。解析は開発コードの Euler 計算（非粘性流）で実施した。解析に使用した計算格子と定常圧力分布を図 6 に示す。解析で予測したフラッタポイント（安定から不安定に移行するポイント）は試験結果と良好に対応しており、開発コードは主翼の遷音速フラッタを高精度に予測できることが確認された。

図2

MRJ 主翼フラッタ模型 洞内設置状態図3

MRJ 主翼フラッタ模型 3 面図

表 1 MRJ 主翼フラッタ模型 代表モードの固有振動数

ModestWing 1 Bending (WB1)
stFuselage 1 Bending (FB1)
Nacelle Pitching (NCP)

GVT

NASTRAN

[Hz]

[Hz]

82.0
152.0
180.0

81.9
150.0
179.6

GVT/
NASTRAN
[-]
1.00
1.01
1.00

民間航空機開発における大規模空力弾性解析シミュレーションの適用

GVT

図4

NASTRANMRJ 主翼フラッタ模型 振動モード形の比較（NCP）

図5

MRJ 主翼フラッタ模型 フラッタ境界— 71 —

— 72 —

SENAC Vol. 47, No. 2（2014. 4）

図6

MRJ 主翼フラッタ模型

Euler 計算格子と定常圧力分布

4.2 MRJ 尾翼フラッタ風洞試験
風洞試験では、エレベータ付きの水平尾翼と胴体フェアリングで構成されるセミスパン模型を使用し、試験は主翼フラッタ風洞試験と同様、三菱重工業（株）名航の高速風洞で実施した。模型の洞内設置状態を図 7 に、模型の 3 面図を図 8 に示す。水平尾翼とエレベータの主構造はアルミ製の桁リブ要素と CFRP 外板で構成され、空力形状はバルサ材と樹脂で整形した。エレベータ部
は 3 ヵ所のヒンジ金具で水平尾翼と結合し、翼根部のトルクチューブを介して風洞壁に固定する方式とした。模型のマス分布、剛性分布は想定される尾翼フラッタモードのうち、最も重要となる舵面フラッタを再現するよう設計した。この舵面フラッタはエレベータの回転モード（ELV-R）と水平尾翼の捩り 1 次モード（HT1）の連成により発生する。模型の構造解析モデルを NASTRANで作成し、振動解析を実施した。振動解析結果を表 2、図 9 に示す。代表的な振動モードについて、解析結果は振動試験結果と良好に対応しており、作成した構造解析モデルの妥当性が確認された。
風洞試験では、マッハ数を固定して風洞総圧（動圧）を段階的に上昇させ、フラッタが発生するポイントを慎重に調査し、複数のマッハ数条件でフラッタポイントを取得した。図 10 に模型のフラッタ境界について、風洞試験結果と解析結果を比較した例を示す。図中の縦軸・横軸の数値については、公開できない情報のため削除しているが、試験のマッハ数条件は遷音速領域であることに留意されたい。解析は開発コードの Navier-Stokes 計算（粘性流）で実施した。解析に使用した計算格子と定常圧力分布を示す。解析で予測したフラッタポイントは試験結果と良好に対応しており、開発コードは尾翼の遷音速舵面フラッタを高精度に予測できることが確認された。

— 73 —

民間航空機開発における大規模空力弾性解析シミュレーションの適用

図7

MRJ 尾翼フラッタ模型 洞内設置状態

表2

図8

MRJ 尾翼フラッタ模型 3 面図

MRJ 尾翼フラッタ模型 代表モードの固有振動数

ModestHT 1 Bending (HB1)
Elevator Rotation (ELV-R)
st
HT 1 Torsion (HT1)GVTNASTRAN

[Hz]

[Hz]

90.6
283.0
527.0

GVT

92.5
303.9
549.5

GVT/
NASTRAN[-]
0.98
0.930.96

NASTRAN

ELV-R

HT1

図9MRJ 尾翼フラッタ模型 振動モード形の比較（ELV-R、HT1）

— 74 —

SENAC Vol. 47, No. 2（2014. 4）

図 10 MRJ 尾翼フラッタ模型 フラッタ境界

図 11 MRJ 尾翼フラッタ模型 Navier-Stokes 計算格子と定常圧力分布

民間航空機開発における大規模空力弾性解析シミュレーションの適用

5.

— 75 —

おわりに民間企業利用サービスにより、世界でもトップレベルの計算機環境を利用することで、当社の計算機環境では困難であった CFD 空力弾性解析の本格活用を MRJ 開発で実現することができた。
解析コードの精度検証をフラッタ風洞試験で行い、主翼、尾翼で重要となる遷音速フラッタを高精度で予測できることが確認された。この解析コードを MRJ の設計段階から積極活用することにより、フラッタに代表される空力弾性不安定問題の設計リスクを最小限に抑え、飛行安全上のリスク低減を図ることができた。今後は飛行試験に向けた準備の中で解析コードを活用するとともに、飛行試験データを用いて解析コードのさらなる精度向上と高度化を図る予定である。

謝辞
本研究開発は、東北大学サイバーサイエンスセンターのスーパーコンピュータを利用することで実現することができた。解析コードの開発では、東北大学の大林研究室、旧中橋研究室にご協力いただいた。また、計算機利用と解析コードのチューニングにあたっては、同センター関係各位に有益なご指導とご協力をいただいた。ここに感謝の意を表します。
参考文献[1] Nakahashi, K., Togashi, F., Fujita, T. and Ito, Y., “Numerical Simulations on Separation ofScaled Supersonic Experimental Airplane from Rocket Booster at Supersonic Speed,”
AIAA Paper 2002-2843, June 2002.
[2] Murayama, M. and Yamamoto, K., “Comparison Study of Drag Prediction for the 3rd CFD
Drag Prediction Workshop by Structured and Unstructured Mesh Method,” AIAA Paper
2007-0258, June 2002.
[3] Ito, Y. and Nakahashi, K., “Surface Triangulation for Polygonal Models Based on CAD
Data,” International Journal for Numerical Methods in Fluids, Vol. 39, Issue 1, 2002.
[4] Shrov, D. and Nakahashi, K., “A Boundary Recovery Algorithm for Delaunay Tetrahedral
Meshing,” Proceedings of 5th International Conference on Numerical Grid Generation in
Computational Fluid Simulations, Mississippi State, Mississippi, 1996, pp.229-238.
[5] Ito, Y. and Nakahashi, K., “Improvements in the Reliability and Quality of Unstructured
Hybrid Mesh Generation,” International Journal for Numerical Methods in Fluids, Vol.
45, Issue 1, May 2004, pp.79-108.
[6] Obayashi, S., et al., “Convergence Acceleration of a Navier-Stokes Solver for Efficient
Static Aeroelastic Computations,” AIAA Journal, Vol.33, No.6, 1995, pp.1134-1141.
[7] Murayama, M., et al, “Unstructured Dynamic Mesh for Large Movement andDeformation,” AIAA paper 2002-0122.
[8] Melville, R. B., et al., “Implementation of a Fully-Implicit, Aeroelastic Navier-StokesSolver,” AIAA paper 97-2039.