― 16 ―SENAC Vol. 57, No. （1 2024. 1）［共同研究成果］

レベルセ ット関数 を用いた 振動する 膜翼の流 体解析
佐々木大輔、船田光星
大阪公立大学高橋俊
東海大学（現所属：宇宙航空研究開発機構）
川本裕樹
東海大学
中矢敦士
金沢工業大学
火星探査航空機を想定した膜翼周りの空力係数を推算するために、等間隔直交格子法において非圧縮性非定常流体解析を実施した。翼形状をレベルセット関数により表現することで、変動する翼を効率的に解析することが可能である。本研究では、強制振動により変形する円弧形状の解析を行い、変動する翼における空力係数を推算すると共に、翼の変動が流れ場に与える影響について調べた。

1．

緒言

火星は他の太陽系惑星に比べて地球に近いが、未だに謎に包まれている部分が多い。そのため、現在でも多くの探査機が送り込まれ調査が進められている。火星は人類が移住可能な惑星として注目されており、どのような資源が存在しているか、生物存在の痕跡がないか、など、火星の環境について調査することは非常に重要である。現在活動している火星探査機は軌道衛星上から地形を観測するオービター、着陸地点から周囲を探索可能なランダーとローバーである。これらの探査機が抱える問題点として、オービターは地表の詳細な情報が確認できないこと、ランダーとローバーは探索範囲が狭いことが挙げられる。これらの課題を克服する新たな探索手法として航空機利用の検討が進められている。航空機による探査では、渓谷や谷など高低差のある場所でも制約なしに移動でき、広く探索が可能となる。
しかし、火星へ航空機を送り込むにあたって解決すべき大きな課題がある。火星の希薄な大気条件において揚力を得るためには面積の大きな翼を必要とし、火星探査用の航空機は大型化せざるを得ない。しかし、火星に輸送する際のペイロードには制約があり、折りたためる、丸められる等、小型化して輸送する必要がある。その一案として、膜翼を主翼に用いることの有効性が検討されている。山本らは、柔軟膜翼の風洞実験を通して、最大揚力係数、最大揚抗比が向上することを示した[1]。その要因は、翼上面の非定常流れ場に応じて膜が振動するためであると述べているが、その詳細なメカニズムについては明らかになっていない。
そこで、本研究では変形した膜翼を模して、強制振動により変形する円弧形状の翼モデルの流体解析を行う。流体解析には等間隔直交格子を採用し、翼形状はレベルセット関数を使用して表現することで、格子の再生成を伴うことなく変動する翼の解析が可能である。解析結果により得られた翼周りの流れ場の変化を可視化し、形状の変動が流れ場へ与える影響を明らかにする。

レベルセット関数を用いた振動する膜翼の流体解析

2．

― 17 ―

解析手法

本解析における支配方程式は式(1)の連続の式と式(2)の 2 次元非圧縮性 Navier-Stokes 方程式とする。
(1)𝛻𝛻𝛻𝛻 𝛻 𝛻𝛻𝛻𝛻 = 0𝜕𝜕𝜕𝜕𝜕𝜕𝜕𝜕
𝜕𝜕𝜕𝜕𝜕𝜕𝜕𝜕

1𝜇𝜇𝜇𝜇+ 𝑈𝑈𝑈𝑈𝑈 𝛻𝛻𝛻𝛻𝑈𝑈𝑈𝑈 = − 𝛻𝛻𝛻𝛻𝛻𝛻𝛻𝛻 + 𝛻𝛻𝛻𝛻 2 𝑈𝑈𝑈𝑈𝜌𝜌𝜌𝜌𝜌𝜌𝜌𝜌

(2)

ここで、U は流体の速度、p は圧力、μ は粘性係数、ρ は密度である。式(1)と式(2)の時間発展には fractional step 法を用いており、離散化は完全等間隔直交格子上で行う。 空間の離散化手法
として式(2)の対流項に 5 次精度 WENO 法、粘性項に 2 次精度中心差分近似を用い、時間の離散化に対流項では 3 次精度 TVD Runge-Kutta 法、粘性項では 2 次精度 Crank-Nicolson 法を適用した。また、圧力の Poisson 方程式には SOR 法を用いている。
物体境界はレベルセット法を用いた埋め込み境界法で定義する[2、3]。セルの中心座標から物体境界までの垂直距離であるレベルセット関数で計算領域内のセルを流体セル、物体セル、物体と流体の中間に位置するゴーストセルに分類し、ゴーストセルの値に任意の境界条件を与える。本
研究ではゴーストセルの値の決定にイメージポイントを用いる。イメージポイントは物体表面から法線方向に伸びるプローブの先端にある点であり，その長さは格子幅の 1.45 倍とした。ゴース
トセルの速度 UGC はイメージポイントの値 UIP を用いて式(3)により決定する。
𝑑𝑑𝑑𝑑 +𝑑𝑑𝑑𝑑𝐺𝐺𝐺𝐺𝐺𝐺𝐺𝐺(𝑈𝑈𝑈𝑈𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 − 𝑈𝑈𝑈𝑈𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 )𝑑𝑑𝑑𝑑𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝑈𝑈𝑈𝑈𝐺𝐺𝐺𝐺𝐺𝐺𝐺𝐺 = 𝑈𝑈𝑈𝑈𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 − 𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼(3)

ここで、イメージポイントの値 UIP はその点を囲むセルから線形内挿により決定する。dGC はゴーストセルから物体表面までの法線方向距離、dIP はプローブの長さ、UIB は物体の移動速度である。ゴーストセルの圧力はイメージポイントに内挿された値をゴーストセルに用いることで物体近傍での勾配が 0 となるように決定した。

3．

解析対象および解析条件

本研究では、図 1 に示す平板翼、上面が変形する翼（上面変形翼）、上下面の両面が変形する翼（両面変形翼）を計算対象として空力係数を算出する。上面変形翼では、平板翼の上面が徐々に膨らみ円弧翼となり、その後また平板翼へと戻る形状変形を繰り返す。両面変形翼は一定の厚みを維持しながら上下共に円弧形状となり、平板翼へと戻る変形を繰り返す。平板翼では翼弦長 cを 40×10−3 [m]、翼厚は 0.05c とする。上面変形翼では、平板翼と同じく翼弦長 40[mm]とし、変形後の円弧の厚みの最大値を 0.07c とする。両面変形翼では厚さ 0.05c を維持し、最大変形時は翼端を結ぶ弦から頂点までの距離を 0.07c とする。
翼形状の変形による影響を調べるために平板翼、変動する円弧翼ともにレイノルズ数
Re=1.0×104 で計算を実施する。変動する円弧翼では、平板翼の形状から翼弦長中心を頂点として楕円状に変形し、最大径を 0.07c とする形状まで変位を与える。その後、元の平板翼に形状に戻るといった形状変化を繰り返す。この際の振動数は 30[Hz]とする。今回の解析で使用した格子は、格子幅 0.125×10−3 [m]、総格子数 512 万(=3200×1600)である。図 2 に両面変形翼まわりの格子と変形形状を示す。

― 18 ―

SENAC Vol. 57, No. （1 2024. 1）

0.07c
0.07c

0.05c
(a)平板翼(b)上面変形翼図 1 解析対象の翼モデル(a)最小変形時

4．

図 2 両面変形翼の変形と格子(c)両面変形翼(b)最大変形時

解析結果

平板翼・上面変形翼・両面変形翼の解析をして得られた平均揚力係数と平均抗力係数を図 3 に示す。平板と上面変形翼を比較すると揚力係数の差は少なく、両面変形翼に比べて低い揚力係数を示す。一方、抗力係数は上面変形翼が低いのに対し、両面変形翼と平板翼は同様の値である。
図 4 に、迎角 3 度における上面変形翼と両面変形翼における主流方向流速分布を示す。両面変形翼では、下面も変動することで、下面の流れ場だけではなく、上面にも低速部が現れている。
また、同様に迎角 6 度における主流方向流速分布を図 5 に示す。図 5 より、最小変形時・最大変形時共に、両面変形翼の上面では低速域が上面変形翼に対して広くなっている。このことが、両面変形翼は揚力が大きく、かつ抗力も大きくなっている原因と考えられる。

(a) 揚力係数
図 3 時間平均化した空力係数（続く）

レベルセット関数を用いた振動する膜翼の流体解析― 19 ―(b) 抗力係数
図 3 時間平均化した空力係数(a)上面変形翼(b)両面変形翼図 4 迎角 3 度における瞬間の主流方向流速分布（左：最小変形時、右：最大変形時）

― 20 ―

SENAC Vol. 57, No. （1 2024. 1）(a)上面変形翼(b)両面変形翼図 5 迎角 6 度における瞬間の主流方向流速分布（左：最小変形時、右：最大変形時）

5.

結言

本研究では、レベルセット法を用いた埋め込み境界法を適用した等間隔直交格子法により、Re=1.0×104 における平板翼と変形する翼の流体解析を行った。平板翼と上面変形翼の円弧形状の際の流れ場を比較すると、上面変形翼は剥離領域が小さいことが確認できた。このことが上面変形翼が両面変形翼に比べて揚力・抗力共に低い要因である。今回、強制的に変形させた解析を実施したが、膜翼特有の現象を詳しく調査するためには、流体構造連成解析の実施に加えて、3 次元解析など、実際の膜翼に近い条件での解析が必要になる。レベルセット法を用いた埋め込み境界法を適用した等間隔直交格子法では、様々な形状変形に容易に対応可能であるため、本手法を拡張してより現実的な解析を実施する予定である。

謝辞
本研究は、令和 4 年度共同研究「直交格子法による移動境界問題の解法に関する研究」として実施した。研究にあたっては同センター関係各位に有益なご指導とご協力をいただいた。また、
本解析は、東北大学サイバーサイエンスセンターのスーパーコンピュータを利用することで実現することができた。
加えて、本解析の可視化には、FieldView CFD 社の FieldView を使用した。

レベルセット関数を用いた振動する膜翼の流体解析― 21 ―

参考文献[1] 山本健太郎、藤田昂志、伊神翼、高橋幸一、永井大樹，“低レイノルズ数における柔軟膜翼の非定常流れ場の PIV 解析”、第 53 回流体力学講演会 / 第 39 回航空宇宙数値シミュレーション技術シンポジウム講演論文集、東京、2021．
[2] Takahashi, S., Nonomura, T., Fukuda, K., “A numerical scheme based on an immersedboundary method for compressible turbulent flows with shocks: application to two-dimensionalflows around cylinders”, Journal of Applied Mathematics, Vol. 20, Article ID 252478, 2014.
[3] Kawamoto, Y., Takahashi, S., Ochiai, M., Azetsu, A., Yamamoto, K, “Prediction of oil
behavior in piston ring groove based on gas-liquid two-phase flow analysis”, Journal ofAdvanced Mechanical Design, Systems, and Manufacturing, Vol. 14, No. 6, 2020.