ラダーのみを用いて旋回する航空機の垂直尾翼・上反角設計

ラダーのみを用いて旋回する航空機について、垂直尾翼容積比 \(V_v\) と等価上反角 \(\Gamma_{\mathrm{eff}}\) を軸にした設計チャートの作り方を説明する。

GitHub - mtkbirdman/OpenVSP

Contribute to mtkbirdman/OpenVSP development by creating an account on GitHub.

github.com

はじめに

ラダーのみを用いて旋回する航空機では、ラダー入力によってヨーイング角速度と横滑りを作り、ヨーとロールの連成項や上反角効果を通じてロールを作る。

したがって、このような航空機の垂直尾翼・上反角の設計では、スパイラル安定性や Dutch roll 特性だけでなく、旋回性能、横転性能、横風突風応答性能、抵抗性能などをまとめて見つつ、バランスの良い垂直尾翼と上反角を設定する必要がある。

そこで、本記事では、ラダーのみを用いて旋回する航空機の垂直尾翼・上反角設計について、横軸に垂直尾翼容積比 \(V_v\)、縦軸に等価上反角 \(\Gamma_{\mathrm{eff}}\) を取った \(V_v-\Gamma_{\mathrm{eff}}\) チャートを用いて検討してみる。

↓対象とする機体はこちら

鳥コン滑空機を想定したSampleGliderのOpenVSPモデル

鳥コン滑空機を想定したSampleGliderのOpenVSPモデルおよび数値例を示す。

mtkbirdman.com

2026.07.16

↓最初から理論を追いたい人はこちらから

航空機のラダーのみ定常旋回

ラダーのみを用いた定常旋回について説明する。

mtkbirdman.com

2026.07.11

それではいってみよう。

\(V_v-\Gamma_{\mathrm{eff}}\) チャートの考え方

\(V_v-\Gamma_{\mathrm{eff}}\) チャートでは、垂直尾翼と上反角の組合せを

横軸に垂直尾翼容積比：

\begin{align} V_v=\frac{S_vl_v}{Sb} \end{align}

縦軸に等価上反角：

\begin{align}
\Gamma_{eff} =\frac{ \displaystyle\int_0^{b/2}y c(y)a(y)\Gamma(y)\,dy }{ \displaystyle\int_0^{b/2}y c(y)a(y)\,dy }
\end{align}

を取った 2 次元平面として考える。

ここで、等価上反角 \(\Gamma_{\mathrm{eff}}\) は、スパン方向に分布するたわみ角を、ローリングモーメントアームで重みづけした平均値として定義する。

\(V_v-\Gamma_{\mathrm{eff}}\) チャート上で、スパイラル安定性、旋回性能、横転性能、横風突風応答性能、抵抗性能などを重ねて見ることで、ラダーのみ旋回に必要な複数の条件を同時に満たす垂直尾翼と上反角の存在領域を探すことができる。

OpenVSP による設計空間の生成

\(V_v-\Gamma_{\mathrm{eff}}\) チャートを作成する流れは以下の通り。

• OpenVSP の基準形状を準備する。
• 垂直尾翼容積比 \(V_v\) を変化させた形状を作る。
• 翼端たわみ量 \(w_{\mathrm{tip}}\) から主翼のたわみ角分布を作る。
• たわみ角分布を OpenVSP の Dihedral へ反映する。
• VSPAERO で安定微係数を取得する。
• .stab ファイルからラダーのみ旋回の評価指標を計算する。
• \(V_v-\Gamma_{\mathrm{eff}}\) 平面に等値線を描く。

1 つずつ順番に説明していく。

基準形状に必要な要件

\(V_v-\Gamma_{\mathrm{eff}}\) チャートを作るには、まず OpenVSP の基準形状が必要である。

本記事では、次のような0Gの基準形状を想定する。

• 胴体、主翼、水平尾翼、垂直尾翼を持つ。
• 主翼、水平尾翼、垂直尾翼は VSPAERO の薄翼解析対象として扱う。
• ラダー、エレベータ、エルロンの Control Surface Group が定義されている。
• VSPAERO の Reference Geom が主翼に設定されている。
• 安定微係数解析で使う Set が整理されている。

↓安定微係数の計算と使用するモデルの詳細についてはこちら

OpenVSP の Python API で安定微係数を計算する

OpenVSPのPythonAPIで、VSPAEROの安定微係数を計算する方法について説明する

mtkbirdman.com

2026.07.15

基準形状で確認する量

基準形状から、次の量を取得する。

\begin{align} S &= S_{\mathrm{ref}}, \\ b &= b_{\mathrm{ref}}, \\ l_v &= x_{\mathrm{ac},v}-x_{\mathrm{cg}}, \\ S_v &= \text{垂直尾翼面積} \end{align}

ここで、\(S\) は基準主翼面積、\(b\) は基準翼幅、\(l_v\) は重心から垂直尾翼空力中心までの距離である。

垂直尾翼容積比は

\begin{align} V_v=\frac{S_vl_v}{Sb} \end{align}

で定義する。

この式を \(S_v\) について解くと、目標の垂直尾翼容積比 \(V_v\) に対して必要な垂直尾翼面積は

\begin{align} S_v=\frac{V_vSb}{l_v} \end{align}

である。

垂直尾翼容積比 \(V_v\) を変化させる方法

今回は、基準形状の主翼面積 \(S\)、主翼翼幅 \(b\)、重心位置 \(x_{\mathrm{cg}}\)、垂直尾翼位置から \(l_v\) を求め、

• アスペクト比固定
• スパン固定

の 2 通りの方法で、\(V_v\) に対応する垂直尾翼の形状を決定する。

目標面積は

\begin{align} S_{v,\mathrm{target}}=\frac{V_vSb}{l_v} \end{align}

である。

垂直尾翼の現状面積を \(S_{v,0}\) とすれば、面積倍率 \(\lambda_S\) は

\begin{align} \lambda_S=\frac{S_{v,\mathrm{target}}}{S_{v,0}} \end{align}

である。

アスペクト比固定の場合は、スパン方向倍率 \(\lambda_b\) とコード方向倍率 \(\lambda_c\) を同じにして拡大縮小する。

垂直尾翼面積は、代表的にスパンとコードの積に比例するので、

\begin{align} \lambda_S=\lambda_b\lambda_c \end{align}

である。

アスペクト比を保つために

\begin{align} \lambda_b=\lambda_c \end{align}

とおくと、

\begin{align} \lambda_S &= \lambda_b^2, \\ \lambda_b &= \sqrt{\lambda_S}, \\ \lambda_c &= \sqrt{\lambda_S} \end{align}

となる。

この場合、アスペクト比は

\begin{align}
AR_{v,\mathrm{new}}
=\frac{(\lambda_b b_v)^2}{\lambda_S S_v}
=\frac{\lambda_S b_v^2}{\lambda_S S_v}
=AR_{v,0}
\end{align}

となり、基準形状の垂直尾翼アスペクト比を保ったまま垂直尾翼面積を変更できる。

一方、スパン固定の場合は、垂直尾翼のスパンを変えずにコード長だけを変える。

したがって、

\begin{align} \lambda_b=1 \end{align}

とすると、面積倍率 \(\lambda_S\) は

\begin{align} \lambda_S=\lambda_b\lambda_c=\lambda_c \end{align}

となる。

このようにして、チャートの横軸である \(V_v\) を直接指定し、その値に対応する OpenVSP 形状を生成する。

等価上反角

長スパン・弾性翼では、荷重によって翼が上方にたわむ。

たわみ角を \(\theta_b(y)\) とすれば、たわみによる上反角は次の式で表される。

\begin{align} \Gamma_{\mathrm{elastic}}(y) \simeq \frac{dw}{dy} =\theta_b(y) \end{align}

元の幾何上反角を \(\Gamma_{\mathrm{rigid}}(y)\) とすれば、たわみ後の局所上反角は

\begin{align} \Gamma_{\mathrm{eff}}(y) =\Gamma_{\mathrm{rigid}}(y) +\theta_b(y) \end{align}

である。

長スパン・弾性翼では、主翼のたわみにより上反角がスパン方向に分布するため、チャートの縦軸として使える代表的な上反角を単一の値として表す必要がある。

横滑り角 \(\beta\) と上反角 \(\Gamma(y)\) による局所揚力変化を

\begin{align} d(\Delta L) \simeq q_\infty c(y)a(y)\beta\Gamma(y)\,dy \end{align}

とすると、この揚力差がモーメントアーム \(y\) を持つため、ローリングモーメント係数への寄与は

\begin{align}
C_{l\beta,\Gamma} \simeq -\frac{2}{Sb} \int_0^{b/2} y c(y)a(y)\Gamma(y)\,dy
\end{align}

である。

このとき、スパン方向に一様な上反角 \(\Gamma_0\) を持つ場合は

\begin{align}
C_{l\beta,\Gamma} \simeq -\frac{2}{Sb} \Gamma_0 \displaystyle\int_0^{b/2}y c(y)a(y)\,dy
\end{align}

となる。

よって、スパン方向に分布する上反角 \(\Gamma(y)\) と同じ \(C_{l\beta,\Gamma}\) を与える一様上反角は

\begin{align}
-\frac{2}{Sb} \Gamma_0 \displaystyle\int_0^{b/2}y c(y)a(y)\,dy
=-\frac{2}{Sb} \int_0^{b/2} y c(y)a(y)\Gamma(y)\,dy
\end{align}

\begin{align}
\rightarrow \quad
\Gamma_0 =\frac{ \displaystyle\int_0^{b/2}y c(y)a(y)\Gamma(y)\,dy }{ \displaystyle\int_0^{b/2}y c(y)a(y)\,dy }
\end{align}

であり、これを 等価上反角 \(\Gamma_{eff}\) と定義する。

今回は、楕円コード長分布・桁径一定を仮定して、翼端たわみ量 \(w_{\mathrm{tip}}\) からたわみ形状を作成し、そのたわみをOpenVSPの0Gモデルに足し合わせて、等価上反角を計算した。

詳細な計算式は、Appendix を参照してほしい。

安定微係数の計算

↓安定微係数の詳しい計算方法はこちら。

OpenVSP の Python API で安定微係数を計算する

OpenVSPのPythonAPIで、VSPAEROの安定微係数を計算する方法について説明する

mtkbirdman.com

2026.07.15

本記事では、OpenVSP / VSPAERO により .stab ファイルを生成し、その .stab ファイルから必要な微係数を読み取るものとして扱う。

プログラムの説明

\(V_v-\Gamma_{\mathrm{eff}}\) チャート作成プログラムは、主に VvGammaChart.py に実装されている。

処理は、次の段階に分けられる。

基準形状の読み取り

read_reference_wing_summary() と read_wing_section_table() により、基準主翼の面積、翼幅、断面情報を取得する。

主翼の面積と翼幅は、垂直尾翼容積比

\begin{align} V_v=\frac{S_vl_v}{Sb} \end{align}

の基準量として使う。

垂直尾翼容積比の変更

update_vtail_volume() により、指定した \(V_v\) から目標垂直尾翼面積を求め、垂直尾翼形状を更新する。

\begin{align} S_{v,\mathrm{target}} =\frac{V_vSb}{l_v} \end{align}

ここで、垂直尾翼形状の変更方法は update_vtail_volume() の vtail_area_scale_mode で指定でき、

vtail_area_scale_mode="fixed_aspect_ratio" では

\begin{align} \lambda_b &= \sqrt{\lambda_S}, \\ \lambda_c &= \sqrt{\lambda_S} \end{align}

vtail_area_scale_mode="fixed_span" では

\begin{align} \lambda_b &= 1, \\ \lambda_c &= \lambda_S \end{align}

となる。

翼端たわみ量から上反角を作る

elastic_wing_deflection_distribution() により、翼端たわみ量 \(w_{\mathrm{tip}}\) からたわみ角分布 \(\theta_b(x)\) を作る。

apply_wing_deflection_as_dihedral() により、得られたたわみ角分布を OpenVSP の Dihedral に反映する。

VSPAERO 安定微係数解析

run_vspaero_stability_case() により、生成した .vsp3 に対して VSPAERO の安定微係数解析を実行する。

基準形状の検証と安定微係数解析については、AnalysisVSPAERO.py の validate_vsp3_for_stability_derivatives() と vsp_stability_derivatives() の考え方に従う。

後処理

postprocess_vv_gamma_cases() により、各ケースの .stab から次の量を計算する。

• \(w_{\mathrm{tip}}\) から計算した \(\Gamma_{\mathrm{eff}}\)
• \(C_{l\beta}\)、\(C_{n\beta}\) などの安定微係数
• spiral_margin
• 旋回性能指標
• 横転性能指標
• 横風突風応答指標

可視化

plot_vv_gamma_contour_panel() により、\(V_v-\Gamma_{\mathrm{eff}}\) 平面上に各指標の等値線を描く。

使い方

test_vv_gamma_chart.ipynb.py と plot_vv_gamma_chart.ipynb.py では、notebook からチャート作成処理を呼び出す例が示されている。

基本的な使い方は次の流れである。

入力条件を設定する

まず、基準 OpenVSP ファイル、出力ディレクトリ、\(V_v\) の候補、翼端たわみ量の候補を設定する。

base_vsp3_path = Path('../../models/BRGlider/BRGlider.0G.vsp3')
output_dir = Path('')

vv_values = np.linspace(0.001, 0.006, 11)
tip_deflections = np.linspace(0.0, 2.0, 21)

次に、VSPAERO 解析条件を設定する。

flight_condition = {
    'alpha_deg': 0.0,
    'mach': 0.0,
    'reynolds': 1.0e6,
}

さらに、形状変更に必要な設定を与える。

geometry_config = {
    'lv': 4.5,
    'wing_area': 18,
    'wing_span': 27,
    'xcg': 1.4,
    'wing_name': 'WingGeom',
    'vtail_name': 'VTailGeom',
    'n_span': 101,
    'vtail_area_scale_mode': 'fixed_aspect_ratio',
}

VSPAERO 安定微係数スイープを実行する

sweep = run_vv_wtip_stability_sweep(
    base_vsp3_path,
    vv_values,
    tip_deflections,
    flight_condition,
    geometry_config,
    output_dir,
    validate_base_model=True,
    fixed_wake_flag=False,
    wake_num_iter=10,
    ncpu=8,
    verbose=1,
)

この処理により、\(V_v\) と翼端たわみ量の組合せごとに .vsp3 と .stab が作成される。

1 形態 4～5 分で、200 形態ほど計算すると 12 時間くらいかかる。

後処理を行う

既存の .stab 群から、等価上反角、安定微係数、旋回トリム、6DOF ロール応答指標を計算する。

results = postprocess_vv_gamma_cases(
    sweep,
    mass=mass,
    inertia=inertia,
    delta_r=math.radians(10.0),
    target_delta_phi=math.radians(2.0),
    turn_trim_mode='gliding',
    turn_trim_phi=math.radians(2.0),
    output_csv_path=output_dir / 'vv_gamma_metrics.csv',
    history_output_dir=output_dir / '6dof_history',
    write_6dof_history=False,
    plot_6dof_history=False,
    verbose=1,
)

だいたい 1 形態 10～15 秒で、200 形態計算するのに 30 分ほどかかる。

チャートを描く

後処理結果 vv_gamma_metrics.csv を読み込み、必要な列を等値線として描く。

fig, ax = plt.subplots()
plot_vv_gamma_contour(results, 'spiral_margin', ax=ax)
plt.show()

複数のパネルへ 1 指標ずつ描く場合は、plot_vv_gamma_contour_panel() を使う。

複数指標の等値線を 1 つの \(V_v-\Gamma_{\mathrm{eff}}\) 平面へ重ねる場合は、plot_vv_gamma_contour_lines() を使う。

考察

ここからは、実際に \(V_v-\Gamma_{\mathrm{eff}}\) チャートを使って、ラダーのみで旋回する航空機の垂直尾翼・上反角の設計について考えてみる。

まずは、スパイラル安定性、旋回性能、横転性能、横風突風応答、抵抗性能などの複合指標を見る前に、各安定微係数とラダー効きを \(V_v-\Gamma_{\mathrm{eff}}\) チャート上で確認する。

安定微係数・ラダー効き

横・方向の線形空力係数は、入力を \(\beta\)、\(\hat p\)、\(\hat r\)、\(\delta_r\)、出力を \(C_Y\)、\(C_l\)、\(C_n\) として整理する。

入力	横力	ローリングモーメント	ヨーイングモーメント
横滑り角 \(\beta\)	\(C_{Y\beta}\)	\(C_{l\beta}\)	\(C_{n\beta}\)
ロールレート \(\hat p\)	\(C_{Y\hat p}\)	\(C_{l\hat p}\)	\(C_{n\hat p}\)
ヨーレート \(\hat r\)	\(C_{Y\hat r}\)	\(C_{l\hat r}\)	\(C_{n\hat r}\)
ラダー舵角 \(\delta_r\)	\(C_{Y\delta_r}\)	\(C_{l\delta_r}\)	\(C_{n\delta_r}\)

↓それぞれの安定微係数の詳細についてはこちら。

鳥コン滑空機における横・方向の微小擾乱理論

鳥コン滑空機のような大きな上反角を持つ機体に対する横・方向の微小擾乱理論について説明する

mtkbirdman.com

2021.06.01

ここで、各記号の意味は以下の通り。

記号	意味	記号	意味
\(a_1\)	主翼の局所揚力傾斜	\(S_f\)	垂直尾翼面積
\(c(y)\)	コード長	\(l_f\)	重心から垂直尾翼空力中心までの前後方向距離
\(\Gamma(y)\)	上反角	\(z_f\)	重心から垂直尾翼空力中心までの上下方向距離
\(z(y)\)	上下位置	\(a_f\)	垂直尾翼揚力傾斜

それでは順番に説明していく。

横滑り角による作用

今回のチャートでは、

\begin{align} -0.213923 &\le C_{Y\beta}\le -0.015450, \\ -0.222588 &\le C_{l\beta}\le 0.001263, \\ -0.023719 &\le C_{n\beta}\le 0.015016 \end{align}

である。

参考記事より、サイドウォッシュ角を無視すると、横滑り角に関する3係数は、

\begin{align} C_{Y\beta} &= -\frac{2}{S}\int_0^{b/2}a_1\sin^2\Gamma\,c\,dy -\frac{S_f}{S}a_f, \\ C_{l\beta} &= -\frac{2}{Sb}\int_0^{b/2} a_1\sin\Gamma \left(y\cos\Gamma+z\sin\Gamma\right)c\,dy -\frac{S_f}{S}a_f\frac{z_f}{b}, \\ C_{n\beta} &= -\frac{2}{Sb}\int_0^{b/2} \left(a_1\sin\alpha+C_L\cos\alpha+C_D\sin\alpha\right) \sin\Gamma\,yc\,dy +\frac{S_f}{S}a_f\frac{l_f}{b} \end{align}

である。

ここで、\(\sin\Gamma\simeq\Gamma_{\mathrm{eff}}\)、\(\cos\Gamma\simeq1\)、\(|z_f/b|\ll1\)、\((z\sin\Gamma/y\cos\Gamma)\ll1\) を用いると

\begin{align} C_{Y\beta} &\simeq -\frac{2}{S} \left( \int_0^{b/2}a_1c\,dy \right) \Gamma_{\mathrm{eff}}^2 -a_f\frac{b}{l_f}V_v, \\ C_{l\beta} &\simeq -\frac{2}{Sb} \left( \int_0^{b/2}a_1yc\,dy \right) \Gamma_{\mathrm{eff}}, \\ C_{n\beta} &\simeq -\frac{2}{Sb} \left[ \int_0^{b/2} \left(a_1\sin\alpha+C_L\cos\alpha+C_D\sin\alpha\right) yc\,dy \right] \Gamma_{\mathrm{eff}} +a_fV_v \end{align}

なので、概念的には

\begin{align} C_{Y\beta} &= C_0+C_{\Gamma^2}\Gamma_{\mathrm{eff}}^2+C_vV_v, \\ C_{l\beta} &= C_0+C_\Gamma\Gamma_{\mathrm{eff}}, \\ C_{n\beta} &= C_0+C_\Gamma\Gamma_{\mathrm{eff}}+C_vV_v \end{align}

となる。

よって、これらの安定微係数について、\(V_v-\Gamma_{\mathrm{eff}}\) の設計空間における感度は、以下のようになる。

• \(C_{Y\beta}\)：\(V_v\) または \(\Gamma_{\mathrm{eff}}\) を増やすと、全体として負の方向へ変化する。
• \(C_{l\beta}\)：主として \(\Gamma_{\mathrm{eff}}\) によって変化し、\(\Gamma_{\mathrm{eff}}\) を増やすと負の方向へ変化する。
• \(C_{n\beta}\)：\(V_v\) を増やすと正の方向へ、\(\Gamma_{\mathrm{eff}}\) を増やすと負の方向へ変化する。\(C_{n\beta}\) に対する寄与のうち、垂直尾翼の寄与が方向安定、主翼の寄与が方向不安定として作用するので、上反角と垂直尾翼のバランスによって、全機として方向安定にも方向不安定にもなりうる。

ロールレートによる作用

今回のチャートでは、

\begin{align} -0.451055 &\le C_{Y\hat p}\le -0.019076, \\ -0.776207 &\le C_{l\hat p}\le -0.771120, \\ -0.0878647 &\le C_{n\hat p}\le -0.0870005 \end{align}

である。

参考記事より、ロールレートに関する3係数は、

\begin{align} C_{Y\hat p} &= -\frac{4}{Sb}\int_0^{b/2} a_1y\sin\Gamma\cos\Gamma\,c\,dy -\frac{S_f}{S}a_f\frac{2z_f}{b}, \\ C_{l\hat p} &= -\frac{4}{Sb^2}\int_0^{b/2} a_1y\cos\Gamma \left(y\cos\Gamma+z\sin\Gamma\right)c\,dy -\frac{S_f}{S}a_f\frac{2z_f^2}{b^2}, \\ C_{n\hat p} &= -\frac{4}{Sb^2}\int_0^{b/2} \left(a_1\sin\alpha+C_L\cos\alpha+C_D\sin\alpha\right) \cos\Gamma\,y^2c\,dy +\frac{S_f}{S}a_f\frac{2z_fl_f}{b^2} \end{align}

である。

ここで、\(\sin\Gamma\simeq\Gamma_{\mathrm{eff}}\)、\(\cos\Gamma\simeq1\)、\(|z_f/b|\ll1\)、\((z\sin\Gamma/y\cos\Gamma)\ll1\) を用いると

\begin{align} C_{Y\hat p} &\simeq -\frac{4}{Sb} \left( \int_0^{b/2}a_1yc\,dy \right) \Gamma_{\mathrm{eff}}, \\ C_{l\hat p} &\simeq -\frac{4}{Sb^2} \int_0^{b/2}a_1y^2c\,dy, \\ C_{n\hat p} &\simeq -\frac{4}{Sb^2} \int_0^{b/2} \left(a_1\sin\alpha+C_L\cos\alpha+C_D\sin\alpha\right) y^2c\,dy \end{align}

なので、概念的には

\begin{align} C_{Y\hat p} &= C_0+C_\Gamma\Gamma_{\mathrm{eff}}, \\ C_{l\hat p} &= C_0, \\ C_{n\hat p} &= C_0 \end{align}

となる。

これを踏まえると、これらの安定微係数について、\(V_v-\Gamma_{\mathrm{eff}}\) の設計空間における感度は、以下のようになる。

• \(C_{Y\hat p}\)：主として \(\Gamma_{\mathrm{eff}}\) によって変化し、\(\Gamma_{\mathrm{eff}}\) を増やすと負の方向へ変化する。
• \(C_{l\hat p}\) および \(C_{n\hat p}\)：今回の設計範囲では変化が小さい。

ヨーレートによる作用

今回のチャートでは、

\begin{align} 0.008091 &\le C_{Y\hat r}\le 0.158519, \\ 0.220813 &\le C_{l\hat r}\le 0.224332, \\ -0.011632 &\le C_{n\hat r}\le -0.007140 \end{align}

である。

参考記事より、ヨーレートに関する3係数は、

\begin{align} C_{Y\hat r} &= \frac{8}{Sb}\int_0^{b/2} C_Ly\sin\Gamma\,c\,dy +\frac{S_f}{S}a_f\frac{2l_f}{b}, \\ C_{l\hat r} &= \frac{8}{Sb^2}\int_0^{b/2} C_Ly \left(y\cos\Gamma+z\sin\Gamma\right)c\,dy +\frac{S_f}{S}a_f\frac{2l_fz_f}{b^2}, \\ C_{n\hat r} &= \frac{8}{Sb^2}\int_0^{b/2} \left(C_L\sin\alpha-C_D\cos\alpha\right)y^2c\,dy -\frac{S_f}{S}a_f\frac{2l_f^2}{b^2} \end{align}

である。

ここで、\(\sin\Gamma\simeq\Gamma_{\mathrm{eff}}\)、\(\cos\Gamma\simeq1\)、\(|z_f/b|\ll1\)、\((z\sin\Gamma/y\cos\Gamma)\ll1\) を用いると

\begin{align} C_{Y\hat r} &\simeq \frac{8}{Sb} \left( \int_0^{b/2}C_Lyc\,dy \right) \Gamma_{\mathrm{eff}} +2a_f\frac{b}{l_f}V_v, \\ C_{l\hat r} &\simeq \frac{8}{Sb^2} \int_0^{b/2}C_Ly^2c\,dy, \\ C_{n\hat r} &\simeq \frac{8}{Sb^2} \int_0^{b/2} \left(C_L\sin\alpha-C_D\cos\alpha\right)y^2c\,dy -2a_f\frac{l_f}{b}V_v \end{align}

なので、概念的には

\begin{align} C_{Y\hat r} &= C_0+C_\Gamma\Gamma_{\mathrm{eff}}+C_vV_v, \\ C_{l\hat r} &= C_0, \\ C_{n\hat r} &= C_0+C_vV_v \end{align}

となる。

これを踏まえると、これらの安定微係数について、\(V_v-\Gamma_{\mathrm{eff}}\) の設計空間における感度は、以下のようになる。

• \(C_{Y\hat r}\)：\(V_v\) または \(\Gamma_{\mathrm{eff}}\) を増やすと、全体として正の方向へ変化する。
• \(C_{l\hat r}\)：今回の設計範囲では変化が小さい。
• \(C_{n\hat r}\)：主として \(V_v\) によって変化し、\(V_v\) を増やすと負の方向へ変化する。

ラダー舵角による作用

今回のチャートでは、

\begin{align} 0.010920 &\le C_{Y\delta_r}\le 0.060707, \\ -0.000612 &\le C_{l\delta_r}\le 0.002702, \\ -0.010628 &\le C_{n\delta_r}\le -0.001633 \end{align}

である。

参考記事より、ラダー舵角に関する3係数は、

\begin{align} C_{Y\delta_r} &= \frac{S_f}{S}a_f\tau, \\ C_{l\delta_r} &= \frac{S_f}{S}a_f\tau\frac{z_f}{b}, \\ C_{n\delta_r} &= -\frac{S_f}{S}a_f\tau\frac{l_f}{b} \end{align}

である。

ここで、\(\sin\Gamma\simeq\Gamma_{\mathrm{eff}}\)、\(\cos\Gamma\simeq1\)、\(|z_f/b|\ll1\)、\((z\sin\Gamma/y\cos\Gamma)\ll1\) を用いると

\begin{align} C_{Y\delta_r} &\simeq a_f\tau\frac{b}{l_f}V_v, \\ C_{l\delta_r} &\simeq 0, \\ C_{n\delta_r} &\simeq -a_f\tau V_v \end{align}

なので、概念的には

\begin{align} C_{Y\delta_r} &= C_0+C_vV_v, \\ C_{l\delta_r} &= 0, \\ C_{n\delta_r} &= C_0+C_vV_v \end{align}

となる。

これを踏まえると、これらの安定微係数について、\(V_v-\Gamma_{\mathrm{eff}}\) の設計空間における感度は、以下のようになる。

• \(C_{Y\delta_r}\)：主として \(V_v\) によって変化し、\(V_v\) を増やすと正の方向へ変化する。
• \(C_{n\delta_r}\)：主として \(V_v\) によって変化し、\(V_v\) を増やすと負の方向へ変化する。
• \(C_{l\delta_r}\)：\(C_{Y\delta_r}\) および \(C_{n\delta_r}\) と比べて小さい。VSPAERO の計算結果では小さな非ゼロ値が現れる。

各種設計指標

各安定微係数の \(V_v-\Gamma_{\mathrm{eff}}\)チャート上での感度が確認できたので、いよいよラダーのみで旋回する航空機の垂直尾翼・上反角設計の設計指標について、\(V_v-\Gamma_{\mathrm{eff}}\)チャート上でどのような値をとっているのかを確認する。

この章で確認する指標は以下の通り。

見る項目	主指標	簡易指標
安定性	\(\delta_{r,\mathrm{trim}}\)	\(C_{l\beta}C_{n\hat r}-C_{n\beta}C_{l\hat r}\)
旋回性能	\(\left|\phi\right|_{\max,\delta_r}\)	\(\displaystyle \left|K_\phi\right|\simeq\left|-\frac{C_{n\delta_r}}{C_{l\beta}C_{n\hat r}-C_{n\beta}C_{l\hat r}}\left(\frac{4m}{\rho Sb}C_{l\beta}+C_{Y\beta}C_{l\hat r}\right)\right|\)
横転性能	\(\displaystyle \frac{\Delta\phi/\Delta t}{\delta_r}\)	\(K_{\mathrm{rudder,roll}}=(C_{l\hat r}-C_{l\beta})C_{n\delta_r}\)
横風突風応答	\(\Delta\phi(\tau_e)\)	\(\displaystyle K_{\mathrm{gust,roll}}=\frac{8I_z}{\rho Sb^3}C_{l\beta}+C_{l\hat r}C_{n\beta}\)

それでは順番に説明していく。

安定性

安定性では、スパイラルモードの安定性と、小バンクのラダーのみ定常旋回に必要なラダー舵角を確認する。

スパイラルモードの判別式は、

\begin{align} D_{\mathrm{sp}} =C_{l\beta}C_{n\hat r} -C_{n\beta}C_{l\hat r} \end{align}

である。

今回のチャートでは、

\begin{align} -0.003321 \le D_{\mathrm{sp}} \le 0.006916 \end{align}

となった。

\(C_{l\beta}\) は主として \(\Gamma_{\mathrm{eff}}\)、\(C_{n\beta}\) と \(C_{n\hat r}\) は主として \(V_v\) で変化するため、中立線 \(D_{\mathrm{sp}}=0\) は両設計変数の組合せで決まる。

小バンクのラダーのみ定常旋回に必要なラダー舵角は、6 自由度定常トリムから得た \(\delta_{r,\mathrm{trim}}\) で確認する。

今回、トリムが成立した 157 ケースでは、

\begin{align} -26.35\ \mathrm{deg} \le \delta_{r,\mathrm{trim}} \le 21.08\ \mathrm{deg} \end{align}

であった。

ちなみに、ラダーのみ定常旋回に必要なラダー舵角の簡易式は、

\begin{align} \delta_r \simeq -\frac{\Omega_0b}{2V_0} \frac{ C_{l\beta}C_{n\hat r} -C_{n\beta}C_{l\hat r} }{ C_{l\beta}C_{n\delta_r} } \end{align}

である。

\(\delta_r\) の値を見てみると、\(D_{sp}=0\) の付近で \(\delta_r=0\) になっている。

そこから\(V_v\) を増やすと、主に \(C_{n\beta}\)、\(C_{n\hat r}\)、\(C_{n\delta_r}\) が変化し、特に\(|C_{n\delta_r}|\) の増加は必要ラダー舵角を小さくする方向に働くが、方向安定とヨー減衰も同時に変わるため、必要舵角はラダー効きだけでは決まらない。

一方で、\(\Gamma_{\mathrm{eff}}\) を増やすと、主に \(C_{l\beta}\) の絶対値が増えるが、\(C_{l\beta}\) は簡易式の分子と分母の両方に入るため、\(\delta_{r,\mathrm{trim}}\) の \(\Gamma_{\mathrm{eff}}\) 感度は単純な反比例にはならない。

したがって、必要ラダー舵角は、スパイラルモードの判別式の値をベースに、\(V_v\) による方向安定・ヨー減衰・ラダー効きの変化と、\(\Gamma_{\mathrm{eff}}\) による上反角効果の変化を合わせて判断する必要がある。

↓詳細な導出はこちら。

航空機のラダーのみ定常旋回

ラダーのみを用いた定常旋回について説明する。

mtkbirdman.com

2026.07.11

旋回性能

旋回性能では、最大ラダー舵角の範囲内で定常的につり合える最大バンク角を確認する。

数値計算における旋回性能の指標は、ラダー舵角を正負の限界値に固定して定常トリムを解いたときの最大バンク角である。

\begin{align}
|\phi|_{\max,\delta_r} =\max \left( |\phi(-\delta_{r,\max})|, |\phi(+\delta_{r,\max})| \right)
\end{align}

片側でトリムが成立しない場合は、成立した側だけを候補とする。

一方、ラダー舵角に対する定常バンク角の簡易指標は、

\begin{align}
\left|K_\phi\right|
\simeq
\left|-\frac{ C_{n\delta_r} }{ C_{l\beta}C_{n\hat r} -C_{n\beta}C_{l\hat r} } \left( \frac{4m}{\rho Sb}C_{l\beta} +C_{Y\beta}C_{l\hat r} \right)\right|
\end{align}

である。

今回、ラダー最大舵角における定常旋回が成立した 182 ケースでは、\(\left|K_\phi\right|\) が発散したケースを除いて

\begin{align}
0.0241\ \mathrm{deg} \le |\phi|_{\max,\delta_r} \le 7.495\ \mathrm{deg}
\end{align}

\begin{align}
0.017700 \le \left|K_\phi\right| \le 0.6266
\end{align}

であった。

今回、数値計算指標と簡易指標を同時に得られた 182 ケースにおける 相関係数は、

\begin{align}
r\left(|\phi|_{\max,\delta_r}, \left|K_\phi\right|\right)=0.989
\end{align}

であり、設計平面上の傾向はよく対応した。

\(|\phi|_{\max,\delta_r}\) は、\(D_{sp}=0\) のライン上で大きくなっている。

\(V_v\) を増やすと \(C_{n\delta_r}\) と \(C_{n\beta}\) の絶対値が増え、ラダー効きと方向安定が同時に変化し、\(\Gamma_{\mathrm{eff}}\) を増やすと \(C_{l\beta}\) の絶対値が増え、横滑りからバンク角を作る経路が強くなる。

よって、旋回性能を高めるには、スパイラル安定性が中立に近づけ、\(V_v\)を大きく、\(\Gamma_{eff}\)を大きくすればよい。

ただし、\(K_\phi\) はスパイラル判別式を分母に含み、スパイラル中立付近では値が発散しているため、この領域では定常トリム計算の結果を優先する必要がある。

↓詳細はこちら。

航空機のラダーのみ旋回における旋回性能

ラダーのみを用いて旋回する航空機の旋回性能について説明する。

mtkbirdman.com

2026.07.11

横転性能

横転性能では、ラダー入力後の有限時間内にどれだけバンク角を作れるかを確認する。

6 自由度有限時間応答から、平均バンク角速度をラダー舵角で正規化した指標を、

\begin{align}
\frac{\Delta\phi/\Delta t}{\delta_r}
\end{align}

とする。

簡易ラダー横転指数は、主要な二つのロール生成経路をまとめて、

\begin{align}
K_{\mathrm{rudder,roll}} =\left( C_{l\hat r} -C_{l\beta} \right) C_{n\delta_r}
\end{align}

である。

今回、応答計算が成立した 174 ケースでは、

\begin{align}
-0.07277 \le \frac{\Delta\phi/\Delta t}{\delta_r} \le -0.009356
\end{align}

\begin{align}
-0.004677 \le K_{\mathrm{rudder,roll}} \le -0.000371
\end{align}

であった。

両指標を同時に得られた 174 ケースの Pearson 相関係数は

\begin{align}
r\left(
\frac{\Delta\phi/\Delta t}{\delta_r},
K_{\mathrm{rudder,roll}}
\right) =0.973
\end{align}

であり、設計平面上の傾向はよく対応した。

\(V_v\) を増やすと \(C_{n\delta_r}\) の絶対値が増え、ラダー入力からヨーレートを作る経路が強くなり、\(\Gamma_{\mathrm{eff}}\) を増やすと \(C_{l\beta}\) の絶対値が増え、横滑りからロールを生成する経路が強くなる。

したがって、横転性能を高めるには、\(V_v\)を大きく、\(\Gamma_{eff}\)を大きくすればよい。

↓詳細な導出はこちら。

航空機のラダーのみで旋回における横転性能

ラダーのみで旋回する航空機の横転性能について説明する。

mtkbirdman.com

2026.07.12

横風突風応答

横風突風応答では、1-cosine 横風突風を受けたときのバンク応答と、安定微係数だけから作る簡易指標を確認する。

数値計算における横風突風応答の指標には、突風通過終了時のバンク角変化

\begin{align}
\Delta\phi(\tau_e)_{max,gust}\ \mathrm{[deg]}
\end{align}

を用いる。

簡易横風突風ロール指数は、横滑り角から直接ロールする経路と、方向安定からヨーレートを介してロールする経路をまとめて、

\begin{align}
K_{\mathrm{gust,roll}} =\frac{8I_z}{\rho Sb^3}C_{l\beta} +C_{l\hat r}C_{n\beta}
\end{align}

である。

今回の範囲は、

\begin{align}
-3.047 \mathrm{[deg]} \le \Delta\phi(\tau_e)_{max,gust} \le 0.974 \mathrm{[deg]}
\end{align}

\begin{align}
-0.009390 \le K_{\mathrm{gust,roll}} \le 0.003321
\end{align}

である。

両指標の相関係数は、

\begin{align}
r\left(
K_{\phi,g}, K_{\mathrm{gust,roll}}
\right)=0.995
\end{align}

であり、今回の設計平面では非常によく対応した。

\(\Gamma_{\mathrm{eff}}\) を増やすと \(C_{l\beta}\) の絶対値が増え、横風から直接ロールする経路が強くなり、\(V_v\) を増やすと \(C_{n\beta}\) と \(C_{n\hat r}\) が変化し、ヨー経由のロール応答とヨー減衰が同時に変わる。

よって、横風突風に対する煽られやすさを小さくするには、\(V_v\)を大きく、\(\Gamma_{eff}\) を小さくすればよい。

↓詳細な導出はこちら。

航空機の横風突風応答

航空機の横風突風応答について説明する。

mtkbirdman.com

2026.07.12

抵抗性能

最後に、シンプルな抗力係数と、ラダーのみ定常旋回における横滑り角をを確認する。

抗力係数の範囲は、

\begin{align}
0.021355 \le C_{D,\mathrm{base}} \le 0.021748
\end{align}

であり、今回の \(V_v-\Gamma_{\mathrm{eff}}\) の範囲では約2%の変動幅がある（たかが2%とみるか、400mの飛距離に対する8m分とみるか）。

なお、この抗力値は、Re数による有害抗力の補正を行う前の値である。

OpenVSPのVSPAEROにおける有害抗力CDoの計算方法

OpenVSPのVSPAEROにおける有害抗力CDoの計算方法について説明する

mtkbirdman.com

2024.12.31

また、垂直尾翼を大きくした設計では、表面摩擦の増加、支持部管掌、隙間、構造重量増加 等の影響も考慮する必要がある。

---

定常旋回中の横滑り角の大きさも抵抗性能の指標になる。

バンク角2°における定常旋回の横滑り角の範囲は、

\begin{align}
1.827\ \mathrm{deg} \le \beta_{\delta_r,\lim} \le 37.745\ \mathrm{deg}
\end{align}

であった。

大きな横滑り角を必要とする設計は、抗力と飛行姿勢の観点から採用しにくい。

垂直尾翼と上反角設計のまとめ

ここまでの結果をまとめると、ラダーのみで旋回する航空機の垂直尾翼と上反角設計は、次の 5 つの性能のバランスで決まる。

• スパイラル安定性
• 旋回性能
• 横転性能
• 横風突風応答
• 抵抗性能

設計変数ごとの効き方をまとめると、次のようになる。

上反角の大きさは、横転性能と横風突風応答のトレードオフで決まる。

横転性能を高めるには大きな上反角が有利だが、横風突風に対する耐性を高めるには小さな上反角が有利である。

したがって、上反角は、垂直尾翼とのバランスでスパイラルモードを中立からやや安定側に置きつつ、横転性能が必要量を満たす範囲に設定する。

垂直尾翼の大きさは、操縦性と抵抗特性のトレードオフで決まる。

ラダーのみ旋回を成立させるには、垂直尾翼をある程度大きくして、方向安定、ラダー効き、ヨー減衰を確保する必要がある。

一方で、過大な垂直尾翼は抗力と重量を増やす。

したがって、フライトプランをよく練り、旋回性能、横転性能、横風突風応答、抵抗性能のどれを優先するかを明確にした上で、\(V_v\) を決める。

このように、\(V_v-\Gamma_{\mathrm{eff}}\) チャートは、ラダーのみで旋回する航空機に必要な操縦性、安定性、突風耐性を同じ平面上で比較する設計ツールとして活用できる。

おわりに

この記事では、ラダーのみを用いて旋回する航空機について、垂直尾翼容積比 \(V_v\) と等価上反角 \(\Gamma_{\mathrm{eff}}\) を軸にした設計チャートの作り方を整理した。

\(V_v-\Gamma_{\mathrm{eff}}\) チャートを使うと、旋回性能、横転性能、突風応答性能、スパイラル安定性を同じ平面上で比較できる。

↓関連記事

鳥コン設計者向け

鳥人間コンテストの機体設計に携わる人に読んでほしい記事

mtkbirdman.com

自家用操縦士（上級滑空機）向け

「自家用操縦士（上級滑空機）向け」の記事一覧です。

mtkbirdman.com

Appendix

ここでは、翼端たわみ量 \(w_{\mathrm{tip}}\) から代表的なたわみ形状を作るため、楕円揚力分布を受ける桁径一定の片持ち翼を考える。

楕円揚力分布を受ける片持ち翼

半スパンを \(s\)、無次元スパン座標を

\begin{align} x=\frac{y}{s} \end{align}

とする。

単位スパンあたり揚力を次でおく。

\begin{align} \ell(y)=\ell_0\sqrt{1-\left(\frac{y}{s}\right)^2} \end{align}

全翼揚力を \(L_w\) とすると、

\begin{align} L_w &=2\int_0^s\ell(y)\,dy \\ &=2\ell_0s\int_0^1\sqrt{1-x^2}\,dx \\ &=2\ell_0s\frac{\pi}{4} \\ &=\frac{\pi}{2}\ell_0s \end{align}

したがって、

\begin{align} \ell_0=\frac{2L_w}{\pi s} \end{align}

である。

片持ち翼の断面 \(y\) における曲げモーメントは、翼端側の分布荷重によって

\begin{align} M(y)=\int_y^s\ell(\eta)(\eta-y)\,d\eta \end{align}

である。

無次元座標で書くと、

\begin{align} M(x)=\ell_0s^2\int_x^1\sqrt{1-\xi^2}(\xi-x)\,d\xi \end{align}

である。

積分を実行すると、

\begin{align} M(x)=\ell_0s^2\left[ \frac{x^2\sqrt{1-x^2}}{6} +\frac{x\sin^{-1}x}{2} -\frac{\pi x}{4} +\frac{\sqrt{1-x^2}}{3} \right] \end{align}

となる。

桁径一定の場合のたわみ角分布

桁径一定、すなわち断面二次モーメント \(I\) がスパン方向に一定の場合を考える。

Euler-Bernoulli はりとして、

\begin{align} EI\frac{d^2w}{dy^2}=M(y) \end{align}

を用いる。

たわみ角は

\begin{align} \theta_b(y)=\frac{dw}{dy}=\int_0^y\frac{M(\eta)}{EI}\,d\eta \end{align}

である。

無次元座標で書くと、

\begin{align} \theta_b(x) =\frac{\ell_0s^3}{EI} \int_0^x \left[ \frac{\xi^2\sqrt{1-\xi^2}}{6} +\frac{\xi\sin^{-1}\xi}{2} -\frac{\pi\xi}{4} +\frac{\sqrt{1-\xi^2}}{3} \right]d\xi \end{align}

である。

この積分を実行すると、

\begin{align} \theta_b(x) =\frac{\ell_0s^3}{EI} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right] \end{align}

である。

\(\ell_0=2L_w/(\pi s)\) を代入すると、

\begin{align} \theta_b(x) =\frac{2L_ws^2}{\pi EI} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right] \end{align}

である。

翼端では \(x=1\) であるから、

\begin{align} \theta_{b,\mathrm{tip}} &=\frac{\ell_0s^3}{EI}\left(\frac{5\pi}{32}-\frac{4\pi}{32}\right) \\ &=\frac{\ell_0s^3}{EI}\frac{\pi}{32} \end{align}

である。

\(\ell_0=2L_w/(\pi s)\) を代入すると、

\begin{align} \theta_{b,\mathrm{tip}}=\frac{L_ws^2}{16EI} \end{align}

である。

さらに \(s=b/2\) を代入すれば、

\begin{align} \theta_{b,\mathrm{tip}}=\frac{L_wb^2}{64EI} \end{align}

となる。

桁径一定の場合のたわみ分布

たわみは、

\begin{align} w(y)=\int_0^y(y-\eta)\frac{M(\eta)}{EI}\,d\eta \end{align}

である。

無次元座標で書くと、

\begin{align} w(x) =\frac{\ell_0s^4}{EI} \int_0^x (x-\xi) \left[ \frac{\xi^2\sqrt{1-\xi^2}}{6} +\frac{\xi\sin^{-1}\xi}{2} -\frac{\pi\xi}{4} +\frac{\sqrt{1-\xi^2}}{3} \right]d\xi \end{align}

である。

この積分を実行すると、

\begin{align} w(x) =\frac{\ell_0s^4}{EI} \left[ \frac{x(4x^2+3)}{48}\sin^{-1}x -\frac{\pi x^3}{24} +\frac{\sqrt{1-x^2}(6x^4+83x^2+16)}{720} -\frac{1}{45} \right] \end{align}

である。

\(\ell_0=2L_w/(\pi s)\) を代入すると、

\begin{align} w(x) =\frac{2L_ws^3}{\pi EI} \left[ \frac{x(4x^2+3)}{48}\sin^{-1}x -\frac{\pi x^3}{24} +\frac{\sqrt{1-x^2}(6x^4+83x^2+16)}{720} -\frac{1}{45} \right] \end{align}

である。

翼端たわみは \(x=1\) を代入して、

\begin{align} w_{\mathrm{tip}} =\frac{\ell_0s^4}{EI} \left(\frac{\pi}{32}-\frac{1}{45}\right) \end{align}

である。

\(\ell_0=2L_w/(\pi s)\) を代入すると、

\begin{align} w_{\mathrm{tip}} =\frac{2L_ws^3}{\pi EI} \left(\frac{\pi}{32}-\frac{1}{45}\right) \end{align}

である。

さらに \(s=b/2\) を代入すれば、

\begin{align} w_{\mathrm{tip}} =\frac{L_wb^3}{4\pi EI} \left(\frac{\pi}{32}-\frac{1}{45}\right) \end{align}

である。

翼端たわみ量からたわみ角分布を作る

チャート作成では、構造荷重や剛性の詳細指定よりも、翼端たわみ量 \(w_{\mathrm{tip}}\) を直接指定する方が扱いやすいため、前節までに含まれていた \(L_w/EI\) を、指定した翼端たわみ量 \(w_{\mathrm{tip}}\) で消去する。

翼端たわみの完全書き下し式より、

\begin{align} w_{\mathrm{tip}} =\frac{2L_ws^3}{\pi EI} \left(\frac{\pi}{32}-\frac{1}{45}\right) \end{align}

である。

これを \(L_w/EI\) について解くと、

\begin{align} \frac{L_w}{EI} =\frac{\pi w_{\mathrm{tip}}} {2s^3\left(\frac{\pi}{32}-\frac{1}{45}\right)} \end{align}

である。

この関係をたわみ角分布の \(L_w/EI\) に代入すると、

\begin{align} \theta_b(x) =\frac{w_{\mathrm{tip}}} {s\left(\frac{\pi}{32}-\frac{1}{45}\right)} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right] \end{align}

である。

さらに \(s=b/2\) を代入すれば、

\begin{align} \theta_b(x) =\frac{2w_{\mathrm{tip}}} {b\left(\frac{\pi}{32}-\frac{1}{45}\right)} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right] \end{align}

となる。

これにより、翼端たわみ量 \(w_{\mathrm{tip}}\) を直接入力として、弾性上反角分布を作ることができる。

VvGammaChart.py では、この考え方に対応して、elastic_wing_deflection_distribution() が翼端たわみ量からたわみ角分布を作り、apply_wing_deflection_as_dihedral() がそのたわみ角を OpenVSP の WingGeom の Dihedral へ近似的に加える。

楕円コード長分布・桁径一定を仮定したときの等価上反角

ここでは、前節で定義した重み付き平均に、楕円コード長分布と桁径一定のたわみ角分布を代入する。

半スパンを \(s\)、無次元スパン座標を

\begin{align} x=\frac{y}{s} \end{align}

とする。

楕円コード長分布を

\begin{align} c(x)=c_0\sqrt{1-x^2} \end{align}

とおく。

等価上反角の重みは、ロールモーメントアーム \(y\) と局所空力寄与 \(c(y)a(y)\) の積である。

局所揚力傾斜 \(a(y)\) を一定とすると、重みは

\begin{align} yc(y)a(y)\,dy \propto (sx)c_0\sqrt{1-x^2}\,s\,dx \propto x\sqrt{1-x^2}\,dx \end{align}

となる。

したがって、楕円コード長分布・局所揚力傾斜一定のもとで、等価上反角は

\begin{align} \Gamma_{\mathrm{eff}} =\frac{ \displaystyle\int_0^1x\sqrt{1-x^2}\,\theta_b(x)\,dx }{ \displaystyle\int_0^1x\sqrt{1-x^2}\,dx } \end{align}

である。

分母は

\begin{align} \int_0^1x\sqrt{1-x^2}\,dx=\frac{1}{3} \end{align}

である。

一方、桁径一定翼のたわみ角分布は、前節で導出したように

\begin{align} \theta_b(x) =\frac{2w_{\mathrm{tip}}} {b\left(\frac{\pi}{32}-\frac{1}{45}\right)} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right] \end{align}

である。

これを等価上反角の定義に代入すると、

\begin{align} \Gamma_{\mathrm{eff}} =\frac{ \displaystyle \int_0^1 x\sqrt{1-x^2} \frac{2w_{\mathrm{tip}}} {b\left(\frac{\pi}{32}-\frac{1}{45}\right)} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right]dx }{ \displaystyle\int_0^1x\sqrt{1-x^2}\,dx } \end{align}

である。

さらに分母の積分値を代入すると、

\begin{align} \Gamma_{\mathrm{eff}} =3\int_0^1 x\sqrt{1-x^2} \frac{2w_{\mathrm{tip}}} {b\left(\frac{\pi}{32}-\frac{1}{45}\right)} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right]dx \end{align}

となる。

この積分を実行すると、

\begin{align} \int_0^1 x\sqrt{1-x^2} \left[ \frac{4x^2+1}{16}\sin^{-1}x -\frac{\pi x^2}{8} +\frac{x\sqrt{1-x^2}(13+2x^2)}{48} \right]dx =\frac{128}{1575}-\frac{\pi}{60} \end{align}

である。

したがって、完全に整理すると、

\begin{align} \Gamma_{\mathrm{eff}} &= \frac{ 6\left(\frac{128}{1575}-\frac{\pi}{60}\right) }{ \left(\frac{\pi}{32}-\frac{1}{45}\right) } \frac{w_{\mathrm{tip}}}{b} \\ &= \frac{48(512-105\pi)}{35(45\pi-32)} \frac{w_{\mathrm{tip}}}{b} \end{align}

である。

数値的には

\begin{align} \Gamma_{\mathrm{eff}}\simeq2.284\frac{w_{\mathrm{tip}}}{b} \end{align}

である。

VvGammaChart.py の calculate_gamma_eff_from_tip_deflection() は、この定義に対応している。

まず elastic_wing_deflection_distribution() により、翼端たわみ量 \(w_{\mathrm{tip}}\) と半スパンから、楕円揚力分布・桁径一定モデルのたわみ角分布 \(\theta_b(x)\) を作る。

次に、

\begin{align} \Gamma_{\mathrm{eff}} =\frac{ \displaystyle\int_0^1x\sqrt{1-x^2}\,\theta_b(x)\,dx }{ \displaystyle\int_0^1x\sqrt{1-x^2}\,dx } \end{align}

を数値積分で評価する。