風を読む！シリーズの解析結果の見方 #2 パラグライダー飛行と乱流

乱流

　以前の記事では流線が渦を巻いている箇所を紹介して、わかり易く「乱流」と呼んでいましたが、大きな渦だけが乱流ではありません。
　乱流とは、風が一見同じ方向に吹いているように見えても、煙や湯気がゆらゆらと揺れながら上がっていくように、細かな渦や突風があちこちで生まれ、ざわついている状態です。この「ざわつき」は規則性がなく、大きさや強さ、そして向きも絶えず変化します。

　パラグライダー飛行は基本、「山などの障害物の風下で飛ばないこと」が鉄則であるため、実際によく遭遇する乱流は2や3のケースです。

　とはいうものの、今のところ私の知識と手持ちのパソコンで扱えるのは上の1だけなので、当「風を読む！シリーズ」ではしばらくの間、山の風下に発生する乱流のみ扱います。
(将来的には、2や3も扱いたいと考えていますが、時期は未定です。）

パラグライダーにとって危険な乱流とその強さ

　パラグライダーが乱流に出会うと、以下のような2つの危険があり、それぞれどの程度の乱流の強さで、その危険が発生するかを考えてみます。　乱流は「渦や突風があちこちで生まれ、ざわついている状態」と説明しましたが、問題を簡単にするため、以下では乱流の強さ=突風の風速として議論を進めていきます。

　ここで紹介する危険となる「乱流の強さ」は、一般的に認知されている値ではなく、当ブログの「風を読む!シリーズ」で山の周りの風の流れの計算結果を確認するにあたっての目安です。
　まだ手探り状態ですので、今後変更する可能性があることをご承知おきください。

危険その1: パラグライダーが煽(あお)られる

山の風下の不規則で強い突風やローターに巻き込まれると、機体が大きく煽られて姿勢を乱し、操縦が難しい状況に陥ります。　
　しかし、強い突風ではなくても、例えばパラグライダーの沈下率と同等の突風が吹いた場合、通常の滑空飛行とは異なり機体が押し上げられたり押し下げられたりすることで翼の迎え角が大きくなったり小さくなったり、突風の翼への当たり方が左右で異なる場合は横揺れが起きたりと機体のバランスを乱される状況が発生し、機体を再度安定させる操作が必要になります。
　そこで、「煽られ」の強さの指標として、パラグライダーの沈下率と同等以上の垂直方向の突風の風速が一つの目安になるのではないかと考えています。
　パラグライダーの沈下率は、前回の記事で「各クラス共通で、沈下率は1[m/s]と仮定」することにしたので、当「風を読む！シリーズ」では、1[m/s]以上の垂直方向の突風が吹くと煽られる危険があると想定して計算結果を見ていきたいと思います。

危険その2: 翼が潰される

　強いサーマル(上昇気流)の中をくるくると旋回しているとき、翼端がサーマルの外に出るとその部分が潰され、翼すべてがサーマルの中に再度戻ると「バンッ」ととても大きな音とともに翼形状が復元し、そのまま飛行を継続できることもあります。
　しかし、左右のどちらか一方のみに突風が当たり片翼を大幅に潰されて意図しないスパイラル(らせん降下)に陥ってしまうこともあれば、翼の前縁を潰されて揚力がなくなり急激に高度を失うような場合もあり、とにかく数秒かつ一部でも飛ぶためにある翼がなくなるので、「翼が潰れる」という状況は全力で避けたいとても危険な状態と言えます。

　どんな突風が翼に当たると「翼が潰れる」状態になるかを理解するために、通常の滑空飛行中にどのように翼が膨らんでいるかを理解する必要があります。
　パラグライダーの翼は、下の図2-1の点線部分のように前縁に空気の取り入れ口(インテーク)があり、インテーク部分は滑空時に空気の淀み点(流速がゼロの点)となって翼内部の圧力$P_0$を保ちます。

　翼内部の圧力は淀み点の圧力$P_0$とほぼ等しく、数式で表すと下記のようになります。
ここで$P$は大気圧、$\rho$は空気密度、$v$は翼の滑空速度です。

$\textstyle P_0 = P + \dfrac{1}{2}\rho v^2$ (1)

　上の式(1)で、大気圧$P$は「静圧」と言い、$\dfrac{1}{2}\rho v^2$は「動圧」(またはラム圧)と呼ばれます。
動圧は空気の速度$v$の2乗に比例する圧力で、台風の時に外に出ると風に吹き飛ばされそうになったり、傘が裏返されたりしますが、あの力もとが動圧です。

　また、図2-1と上の式(1)から、滑空中の翼内圧は大気圧よりラム圧分だけ大きい圧力によって膨らんでいるということがわかります。
　ここで、大気圧$P$を基準値として0(ゼロ)に設定し、大気圧に対する圧力の差分をゲージ圧と言います。
上で見た翼内の動圧$\dfrac{1}{2}\rho v^2$も大気圧$P$に対するゲージ圧で、さらに、翼の周りの空気の流れの説明で、翼の上側は「負圧」で翼を吸いあげ、下側は「正圧」になっていて翼を押し上げているという話を聞いたことがあると思いますが、この「負圧」や「正圧」も大気圧に対するゲージ圧の話です。
　ゲージ圧には、大気圧との差分だと分かるように$“\Delta”$(デルタ)を先頭につけることにします。例えば翼内のラム圧は、$\Delta P_0 = \dfrac{1}{2}\rho v^2$のように記載します。

空気の密度1.166kg/m³@20℃とすると、滑空速度10m/sでの翼内のラム圧は以下のようになります。
$\Delta P_0 = \dfrac{1}{2} \times 1.166 \times 10^2 =58.3$Pa $\simeq 5.945$kgf/m²

　次に、翼断面で見た場合、翼に対してどの方向の突風が吹くと翼が潰れやすいかを考えてみます。

　翼の下方向からの突風に対しては、パラグライダーの翼は下方向にラインで引っ張られているので潰されることは無く、滑空方向からの突風はラム圧を上昇させるので翼を潰れにくくします。そのため、パラグライダーの翼が潰れるような突風の方向は、大きく分けて図2-2に示すMode1,2と3の様な翼の上方向からの風(下降気流)と、Mode 4,5と6の様な滑空方向と逆方向からの風(追い風)だと考えています。

　翼の上方向から下降してくる突風は、翼の上面に負圧を発生させる滑らかな流れを蹴散らして翼面に衝突すると、そこで動圧を発生します。この動圧が翼内のラム圧を上回ると、翼はその形状が維持できなくなり、「潰れ」が発生します。
　上方向からの突風による動圧についてもう少し詳しく見ると、翼のような曲面に風が当たる場合、動圧はその曲面に対して垂直方向の風速成分が動圧に寄与します。例えばMode 2のように真上から風速$v_G$の突風が吹く場合、図2-3に示すように、突風方向から見た頂点aでは風速$v_G$がそのまま動圧に寄与しますが、翼面に対して直角方向と突風の角度($\theta$とします)がある点b, cでは垂直成分が$v_G cos (\theta)$ですので、動圧は$\Delta P_G = \dfrac{1}{2}\rho\left( v_G cos(\theta)\right)^2$となります。

つまり、翼上方向からの突風により、翼断面が全体的に潰れるための条件は下の式のようになります。

$\Delta P_G > \Delta P_0$　⇒　$\dfrac{1}{2}\rho \left(v_G cos(\theta) \right)^2 > \dfrac{1}{2}\rho {v}^2$　⇒　$v_G cos(\theta)> v$　　　(2)

　ただ、$\theta = 0$の頂点aでも突風の風速$v_G$で翼形状が不安定になり、ここが起点になって翼形状の潰れが起きる可能性もあると考えられるので、安全を考慮し、翼の上方向から吹く突風(下降気流)の風速が翼の滑空速度と同じ($v_G = v$)になる場合を危険の目安とします。

$cos (0) = 1$で、$cos(\theta)$は0以上、1以下の値になります。($0 \leqq cos (\theta) \leqq 1$)
そのため、$v_G$に$cos(\theta)$を掛けるとその大きさは、$0 \leqq v_G cos(\theta) \leqq v_G$となり、元の風速$v_G$以下となります。

　次に、追い風成分の突風について説明します。
前述の通り、パラグライダーの翼は滑空で生じるラム圧で翼形状を保っています。しかし、翼の後ろから追い風の突風が吹いて、翼のインテーク部分で翼に対する空気の相対速度がゼロになると、ラム圧もゼロになり翼形状が不安定になります。そこに、別の突風や渦がかさなることで、翼が潰れる状態が生じると考えています。
　追い風の突風が、滑空方向からの空気の流れに対して角度($\theta$)を持つ場合、Mode 4を例に図2-4で図解してみます。

追い風成分が翼の滑空方向と角度($\theta$)を持つ場合、突風の滑空方向成分は$v_G cos \theta$となります。これが、翼が滑空していた速度と同等になるとラム圧がゼロになります。また、上方向からの突風の説明でも紹介しましたが、$cos \theta$は1以下の値となるため、追い風成分も滑空方向からの空気の流れ対する角度がゼロ($\theta = 0$)で$v_G$が一番厳しめの値になります。

　以上から、当「風を読む！シリーズ」では、翼の上方向からの突風も、追い風方向からの突風も、滑空速度と同等以上(風速$v\leqq v_G$)の突風が「翼が潰れる」条件の目安とします。
　パラグライダーの滑空速度＝トリム速度の場合、前回の記事でトリム速度は一律10m/sとしたので、$v_G=$10m/s以上の突風が危険の目安となります。

　他にも、前方からの突風で翼の迎角が小さくなり、その結果、インテークに対する空気速度も小さくなって前縁が押されて潰れるような機体のピッチングが「翼の潰れ」の原因となるようなメカニズムがあると思いますが、今回はこのような状況も上に示した目安に包括されると考えて進めます。

突風による 危険の種類	危険な突風の基準	危険な突風 風速の目安
煽られる	パラグライダーの沈下率以上	1m/s
翼が潰される	パラグライダーの滑空速度以上	10m/s

　パラグライダーは翼の潰れが生じたとしても、直ちに墜落するわけではなく、基本的に自然に回復するように設計されていますし、半分程度の潰れであれば、適切な操縦により飛行を続けることが可能です[1]。
　とはいえ、乱流の中では「適切な操縦」ができる保証はなく、場合によっては墜落事故につながりかねないので、「乱流」がどの辺りにあるかわかってる場合は、そこを避けて飛行するべきです。

突風の風速(平均風速からの揺らぎ)

　上の章で、乱流は「規則性がなく、大きさや強さ、そして向きも絶えず変化する」と説明し、さらに乱流の強さ=突風の風速としました。突風の風速をもう少し詳しく説明すると、平均的な風速に対して、どれくらい風速が揺れ動いているかのことです。通常、「風速3m/sの風」などと話しますが、この3m/sという値は、一定時間(例えば1時間)の平均風速です。しかし、実際の風は、一瞬強くなったり弱くなったりと常に揺れ動いています。この平均からの揺れ幅つまり「平均風速からの揺らぎ」だけなら、私の持っているパソコンでもOpenFOAMで計算・可視化することができます。以前の記事で紹介した乱流運動エネルギー$k$の計算結果を用いて、平均風速からの揺らぎを算出する方法を解説します。

乱流運動エネルギーkとは？

　まず、高校物理の復習もかねて、運動エネルギーの基本をおさらいしましょう。
質量$m$の物体が速さ$v$で動いているとき、その運動エネルギー$E$は以下の式で表されます。

$E=\dfrac{1}{2}mv^2$　(3)

　風の「揺れ動き(風速の変動)」にも同じ考えが当てはまります。
全体の平均的な風速を差し引いた、風の「揺らぎ(変動成分)」に着目します。
単位質量(1kg)の空気塊が、$x, y, z$の3方向にそれぞれ$u\prime, v\prime, w\prime$という大きさで揺らいでいるとすると、その揺らぎが持つ運動エネルギー(乱流運動エネルギー$k$)は次のようになります。

$k = \dfrac{1}{2} \left( u\prime^2 + v\prime^2 +w\prime^2 \right)$(単位: m²/s²)　(4)

ここで、話を分かりやすくにするために「3方向のゆらぎの大きさは同じ($u\prime=v\prime=w\prime$)と仮定します(この仮定を「等方性乱流の仮定」と言います)。すると、乱流の運動エネルギーは、以下のようにシンプルな形で表せます。

$k = \dfrac{3}{2} u\prime^2$(単位: m²/s²)　(5)

　当シリーズで紹介する、OpenFOAMでの計算は、実は上記の説明のように、風の一瞬一瞬の変化を追うのではなく、時間平均した風の流れを計算します。波に例えると、「今、何メートルの波が何回打ち寄せたか？」を一つずつ数えるのではなく、「その海域全体の波の荒れ具合(波浪エネルギー)」を計算しているイメージです。そのため、個々の突風を直接計算しているわけではありませんが、エネルギーの収支としてその場所に蓄えられる$k$の値(平均的な荒れ具合)は一定の値として求まります。(この計算の背景には、「レイノルズ応力」などの乱流計算モデルの理論がありますが、ここではその説明は省きます。)

風速の変動

式(5)の形を変えて$u\prime^2$を左辺にもってくると以下のようになります。

$u\prime^2 = \dfrac{2}{3}k$ (単位: m²/s²)　(6)

この$u’^2$は統計学で「分散(データのバラつき具合)」に相当します。

　ここで、分散の意味を具体的な例で考えてみます。
　たとえば、山頂で風速を朝の9時から10分ごとに5回計測し、以下の表3-1の結果になったとします。
平均風速を計算すると10 m/sです。では、「平均風速からどれくらい風速が揺れているか(プラスマイナスの揺れ幅)」の平均を求めたいときはどうすればよいか？　単純に「平均との差(偏差)」の平均をとるとプラスとマイナスが打ち消しあって合計(平均)がゼロになってしまいます。そこで、偏差を二乗してから平均をとることで、プラス側に揃えて振れ幅の大きさを計算します。この「偏差の二乗の平均」のことを分散($\sigma^2$)と呼びます。(“$\sigma$“は「シグマ」と呼びます)

項目	9:00	9:10	9:20	9:30	9:40	平均
風速 m/s	8	10	12	9	11	10
風速平均との差(偏差) m/s	-2	0	2	-1	1	0
偏差の二乗(分散$\sigma^2$ ) m2/s2	4	0	4	1	1	2

　表3-1で求めた分散$\sigma^2=2$は、一度、二乗しているため、単位がm²/s²になってしまっています。
これをもとの風速の単位 m/sに戻すために平方根を取ったものを標準偏差($\sigma$)と呼び、これが「平均的な揺らぎの大きさ」を表します。表3-1の例では$\sigma$は以下のようになります。

$\sigma = \sqrt 2 \simeq 1.41$ m/s

　話を式(6)の$u\prime^2$に戻すと、標準偏差、つまり「平均風速からどれくらい風速が揺れているか」の平均は以下のようになります。次の話につなげやすく$u\prime=\sigma$とします。

$\sigma = \sqrt{\dfrac {2}{3} k}$ (単位: m/s) (7)

　式(7)から風速の標準偏差$\sigma$がわかりましたが、一番知りたい風速の揺らぎ、つまり突風の現実的な最大値について、$\sigma$の何倍まで考慮する必要があるかを考えてみます。
　表3-1では平均風速10m/sでしたが、平均風速3m/sの日に1sごとに5時間計測したイメージグラフを図3-1に示します。実際に風速を計測したわけではなく、平均風速3m/sで乱流運動エネルギー$k$=0.5 m²/s²として乱数を使って波形を作成しています。5時間の理由は、パラグライダーのクロスカントリーでの一般的な飛行時間が5±2時間くらいのようなので、ここでは5時間と設定しています。
　図3-1の右側が各風速の発生頻度をグラフ化したもので、このグラフは正規分布と呼ばれていて、風速の揺らぎ方もおよそこの分布になることが知られており、$\sigma$の倍率が大きくなるにしたがって発生頻度は低下していきます。
　また、図3-1は風速の揺らぎのイメージを、わかり易く、等方性を仮定した1方向(風向き方向)だけの揺らぎ成分を取り出して描いたものです。実際には、すでに説明した通り、風速の揺らぎの方向は絶えず規則性無く変化しています。
　当ブログでは、風速の揺らぎの大きさは正規分布に従うと仮定して話を進めます。

各$\sigma$の超過確率(例えば4$\sigma$の超過確率は正規分布で4$\sigma$の外側の確立のこと)を表にすると以下のようになります。

$\sigma$	1	2	3	4	5
確率 ($\sigma$以内)	68.269%	95.450%	99.730%	99.994%	100.000%
超過確率 ($\sigma$以上)	31.731%	4.550%	0.270%	0.006%	0.000057%
超過確率にあたる秒数（飛行時間ごと）
1時間 (3600秒)	1142.3	163.8	9.7	0.2	0.002
2時間 (7200秒)	2284.6	327.6	19.4	0.5	0.004
3時間 (10800秒)	3427.0	491.4	29.2	0.7	0.006
4時間 (14400秒)	4569.3	655.2	38.9	0.9	0.008
5時間 (18000秒)	5711.6	819.0	48.6	1.1	0.010

表3-2を確認すると、揺らぎの大きさが4$\sigma$でも5時間の飛行を行うと1.1秒発生することになります。図3-1の場合、5時間飛行すると1.1秒だけ風速の最大揺らぎつまり突風が5.3(4$\sigma$での風速)-3(平均風速)=2.3m/sに達する可能性がありうるわけです。
　約1秒とはいえ、4$\sigma$の突風が発生する可能性があるのと、2026年9月時点、計算モデルがかなり粗いので、安全側に考えて当ブログでは風速の揺らぎは4$\sigma$まで考慮することにします。

計算結果確認方法とこの計算方法の限界

　計算結果確認方法

当ブログでの、OpenFOAMの計算結果を可視化した図の見方を説明します。

　図4-1は乱流の計算結果の例で、2, 3, 4m/sの突風が吹くエリアを分けて可視化したものです。
図中のA, Bに以下のような情報を示しています。

1. 上から順に山の名前、風の方向と風速(計算条件)、表示している解析結果の説明
   今回の例では、以下のようになっています。
   Ohira-san(日本語入力ができないので、ローマ字表記)
   
   wind: E at 6m/s
   風向が東(E)で、風速6m/sという初期条件で計算
   
   Gust 4sigma
   表示している突風(Gust)が正規分布$\sigma$の何倍かを示しています。
   今回は”4sigma”ですので、4倍です。(“$\sigma$“を直接入力できないので”sigma”としています)
   図4-1は4$\sigma$の突風2, 3, 4m/sのエリアを可視化していますので、重なっている部分も確認できるように、突風エリアを半透明にしています。
   
2. 表示している解析結果の成分と値の大きさごとの色分け凡例
   図4-1は4$\sigma$の突風分布を表しているので”Gust 4sigma”で単位は”m/s”と示しています。
   何かの理由があって、3$\sigma$の突風分布を表すときは”Gust 3sima”とします。

現時点の突風計算の限界

　図4-1は、大平山の東風、風速6m/sの計算結果です。
この計算を行った理由は、知人から以下のような依頼があったからです。

　大平山のランディングでは、東風の時にランディングの東にある山からローターが発生して着陸に影響する。
過去には東風が強い日に、これが原因で翼が潰れて重大事故も発生していている。
　大平で飛ぶ人は、東風の時にはランディングから離れた、ローターの影響がなさそうな場所に着陸するよう気を付けているが、着陸に影響のある範囲の目安を調べられないか？

　それで、大平のランディング近くで翼が潰れるであろう10m/sの突風が発生する条件を探していたのですが、東風の風速6m/sという強めの風の環境下でもランディング上空で2m/sの突風しか発生しないという計算結果となっています。この計算結果から単純な逆算をすると、大平のランディング上空で10m/sの突風が発生するには、風速30m/sの東風と台風並みの風が必要となりますが、そんな環境で飛ぶ人はさすがにいません。
　また、下の動画4-1を見ていただくと、ランディングの隣の山ではなく、二つ隣りの山からの乱流がランディング上空付近に届いている結果となっており、大平で飛んでいる方々の意見と違う結果になっている。

　大きな視点で見た場合、私がOpenFOAMを用いて行っている計算自体は大外れの結果にはなっていないと考えてはいますが、実際には翼が潰れるような事故も起きているわけで、計算結果と現実が合っていないということになります。
　このズレの原因は計算の境界条件や計算モデルに木や建物が入っていない、Meshサイズが30mと粗い、実際はウィンドウシアーなどが発生していてその影響があったか、そもそも上述の翼が潰れる条件が間違っている等々、いろいろ考えられますが、少しずつ探っていく予定です。
　まだ手探り状態ではありますが、パラグライダー飛行で注意すべき突風の強さやおよその危険エリアを可視化できるようになったので、最初の一歩は踏み出せたかなと考えています。当ブログで紹介していく乱流エリアの計算結果は、局所的な部分までは網羅できていないということを念頭に、あくまで安全飛行のための「目安」として参考にしていただければ幸いです。

参考文献

[1]パラグライダー 最新テクニック ブック　著者名 日本パラグライダー協会(監修 ) パラワールド編集部
　ISBN978-4-8022-0285-5

[2]国土地理院の地図データ
https://www.gsi.go.jp/tizu-kutyu.html
計算モデルの3D形状は国土地理院の3Dデータを用いています。
また、計算モデルにマッピングした空中写真も国土地理院のものを使用しています。

今回用いたオープンソースのソフトウェア

　図3-1は、下記のオープンソース(フリー)のソフトウェアを利用しました。これらを利用することで、企業や大学で行うような計算や、その可視化を自宅で行うことができます。とてもすごい時代ですね！

Blender

https://www.blender.org/
3Dモデリングなどを行える3DCGソフトウェア
ダウンロードした国土地理院の3Dデータ(stlファイル)を、OpenFOAM用の計算モデルへ変換するために使用

OpenFOAM

①https://www.openfoam.com/ or
②https://openfoam.org/
流体や熱移動、化学反応などをコンピューターシミュレーションできるオープンソース熱流体解析ソフトウェア
このブログで紹介する山周辺の風の流れの計算は①を用いています。

ParaView

https://www.paraview.org/
科学技術計算データを3Dで可視化することができるオープンソースのソフトウェア
OpenFOAMの計算結果の可視化に使用しています。

QGIS

https://www.qgis.org/ja/
地図の作成、分析ができるオープンソースの地理情報システム ソフトウェア
計算モデルの地表への空中写真のマッピング用データ作成に用いています。

Gnuplot

http://www.gnuplot.info/
プログラミングの様な方法で指示をする必要がありますが、とても強力なグラフ作成ソフトウェアです。
Excelでは作成が難しいグラフでも作ることができます。
今回は、図3-1を作成するのに用いました。