〔原著論文〕
Linear instability of the boundary layer flow along an archery arrow Toshinari HASEGAWA, The University of Erectro-Communications Takeshi MIYAZAKI, The University of Erectro-Communications Kentaro TANAKA, Tokyo Denki University Naoya TAKAHASHI, Tokyo Denki University Reima IWATSU, Tokyo Denki University (Received 8 November 2018; revised 4 February, 2019; accepted 24 April, 2019; published 25 June 2019) The linear instability of the boundary layer flow formed along an archery arrow is investigated based on the parallel flow approximation. The present study covers a Reynolds number range 1.0 × 10 4 ≤ Re ≤ 1. 2.0 × 10 4 (diameter of the arrow shaft being used as the length scale) , whose highest value corresponds to a typical arrow velocity shot from a recurve bow. Two types of points, i.e. streamlined point and bulge, are attached to the arrow head and the axisymmetric basic flows are numerically computed by a stream function-vorticity method.  The linear disturbances are Fourier decomposed in the axial and azimuthal directions. A Chebyshev collocation method is used in the radial discretization. The most quickly growing mode has the azimuthal wave number m = 2 and the second and third modes are with m = 1 and m = 3, respectively. We find that the axisymmetric mode with m = 0 never grows in the studied Reynolds number range.  These results are common to two types of points.  The estimated N -value is less than 2.5 even at Re = 2.0 × 10 4 , indicating that the boundary layer may remain laminar, if the arrow aligns with the local wind, as was observed in the previous MSBS wind tunnel tests of our group. It is also indicated that the laminar-turbulent transition found in our free flight experiments might be a result of bypass transition process caused by a small angle of attack of the arrow inflight. (KEY WORDS) : Archery arrow, Boundary layer flow, Linear instability, Parallel flow approximation, e N method 1絡言オリンピックでのアーチェリーは70 m先の的を狙う得点競技であり，メダルを目指すためには直径12.2 cmの中心円を射貫くことが求められる。競技者の放つ矢の速度Uはおよそ60 m/sで，矢の直径d = 2a（aは半径）を代表長としたRe数（Re = Ud/v：vは空気の動粘性率）は約2.0 × 10 4となる。飛翔中の矢に働く空気抵抗は大きく重力の0.5倍を超えるが，その半分以上は矢のシャフトに働く摩擦抵抗からの寄与で，残りは矢羽に働く空気抵抗である．矢側面の境界層流れが層流である場合と乱流である場合では空気抵抗が有意に異なり,鉛直方向に20 cm程度の的ずれが生じる．このような矢の飛翔時の乱流遷移の詳細なメカニズムを解明して、それを制御することはスポーツ流体力学の観点からも重要である.矢の側面の境界層流れに関する実験的な研究として, Park et al. ¹）の水槽実験が知られている．椎型鏃をつけたアーチェリー矢のシャフト側面の境界層を可視化し, Re = 2.74 × 10 4では先端から乱流境界層になることを示した.一方，矢の空力特性の定量的な測定は，JAXAの60 cm磁力支持天秤装置付き風洞（MSBS風洞）を用いた実験や飛翔実験によって行われている．鈴木ら²ʲはクロスポウ（細長比57）の矢の抗力係数CD，揚力係数CL,ピッチングモーメント係数CMを求めた.測定された10 4 ≤ Re ≤ 5 × 10 4の範囲では，境界層の状態は鏃形状に依存することが示された．流線形鏃の場合には層流境界層であり、椎型鏃および凹凸型鏃を装着した場合は境界層は乱流状態となった.アーチェリー競技者の使用するEaston社製のA/C/E矢の空力特性はMiyazaki et al.  3,4)の支持干渉のない実験により詳細に調べられている．JAXAのMSBS風洞を使用した測定では，鏃形状に関わらず矢を気流に対して平行に置く場合にはRe ≤ 1.3 × 10 4で境界層は居流となる。一方，椎型鏃を装着した矢を気流に対して0.75度傾けた場合には，Re = 1.1 × 10 4付近で乱流適移することが示された．また，圧縮空気によって矢を発射する飛翔実験結果は，椎型鏃を装着した場合は1.3 × 10 4 < Re < 1.8 × 10 4の領域でCDが居流境界層に相当する値（CD ～ 1.5）と乱流境界層に相当する値(CD ～ 2.7)に二値化し，流線形態の場合ではこの乱流遷移領域がより高Re数領域(1.3×10 4 ≤ Re ≤ 2.4×10 4 )に広がった.  Miyazaki et al. 4)はCDの二値化の原因が矢の迎角のばらつきにあることを指摘し，矢の迎角が気流に対して±0.6度以上つくと境界層が乱流遷移すると指摘した．実際，迎角がつかないように放たれた椎型鏃の矢のCDは，遷移領域のRe = 1.7 × 10 4でも層流値（CD = 1.34）となることが確認されている．これらの実験結果を理論的に説明するための第一段階として,境界層流れの線形安定性を調べることが本論文の目的である。円柱側面の境界層流れの線形安定性についてはTutty et al. 51による平行流近似に基づいた局所的な線形安定性解析が知られている。彼らは最も低い臨界Re数( Re c = 2120）を与える不安定モードが周回方向波数m = 1のモードであることを示した．また，m = 0, 2, 3のモードの臨界Re数はそれぞれRec = 24878,12140,20204であることが示された。アーチェリー矢のRe数領域ではこれらの不安定モードが成長する可能性がある。また,飛翔中の矢は矢羽の効果によって100 rps程度で回転しており，無次元回転数（Sp = πdf/U ～ 0.03：fは回転数）は小さいものの，その影響も調べる必要がある． Muralidhar et al. 6.7)は対称軸周りに回転する円柱側面の境界層流れに対して，平行流近似に基づく線形安定性解析を行い，線形不安定性成長率は非回転時に比べて顕著に大きくなることを示している．また，Herrada et al.8)は非平行流効果を放物化された安定性方程式（PSE）に基づいて調べている。しかしながら，これらの理論的研究では非常に細長い（細長比105）円柱の下流側が主な解析対象とされたため，先端ではBlasius境界層が仮定されており，鏃形状の影響は考慮されていない．したがって、細長比が120程度のアーチェリー矢側面に発達する境界層流れに焦点を絞った線形安定性解析を実行しなおす必要がある.本論文では、2つの鏃形状（流線形鏃と椎型鏃）を反映した軸対称基本流を精密に計算し，その線形安定性解析を平行流近似に基づいて行う．基本流として，鏃を装着した矢羽無しの矢（ベアーシャフト）を想定し、シャフトは回転しないもの（Sp = 0）として計算を行う．線形安定性の解析手法は，Muralidhar et al. 6)の高精度な手法に従う.軸方向(波数α)と周回方向(波数m）についてスペクトル分解し，半径方向の離散化にはChebyshev選点法を使用する．得られた解析結果に基づいてMiyazaki et al. 3,4)の実験結果を考察する.特に,不安定性をもたらす撹乱モードについてN値を評価し,乱流遷移との関係を議論する.
2　軸対称基本流の数値計算
流線形鏃および椎型鏃を装着したベアーシャフトを過ぎる基本流の計算を行う．シャフト部分は円柱とみなし，半径をaとして長さが240α（細長比120）になるように設定する。以下ではaを長さスケール，一様流速Ｕを速度スケールとして基礎方程式を無次元化する.また，矢筈（矢後端）は半球とする．座標系は矢軸方向をz，半径方向をr，周回方向をθとして，鏃の付け根を原点にとる。軸対称基本流成分を(Uz,Uz,O)とするとき，Stokesの流れ関数と周回方向の満度成分，
を導入すれば，軸対称基本流の基礎方程式はゆ，ωの関係式および渦度輸送方程式で表される。
式(3)右辺のRe数は直径を代表長さにしている．上式(2,3)を時間発展させて，定常解を求める．計算に使用するのは時間４次，空間４次精度の差分法コード9で，空間方向には一般座標変換を行い，メトリクス成分および計算平面での導関数はともに4次の中心/片側差分で近似する．式(2)は３レベルのマルチグリッド法で解き，式(3)はCarpenter - Kennedyによる5段4位の低容量ルンゲ・クッタ法10で時間積分する.物理平面(z,r)に対して計算平面を(ε,η)としてO-型格子を生成する.矢(鏃とシャフトおよび矢筈)表面に沿う方向をを,矢表面でを線と直交し表面から離れる方向をηとおけば，矢表面はη = 1，遠方境界はη = N n , ,鏃前方および矢筈後方の対称軸はそれぞれ= 1, ξ = Nεで表される．ここで，遠方境界はシャフト中点より10aだけ後方の点を中心とする半径240αの半円とした。このとき，境界条件は矢表面で流速が0，遠方で一様流，軸上で対称条件となり，ゆとωに対して計算平面上において以下の式で与えられる。
格子生成にはEisemanによる4平面補問法!!)を使用した．格子点数はNε × N π = 1201 × 301とする．鏃の形状モデルはCADデータに自然スプラインとVinokurの集中関数¹²を用いて作成する。矢表面におけるη方向の最小格子幅は2.0 × 10−3，ξ方向の最小格子幅は2.0 × 10−²とする．図1に，実寸の流線形鏃と椎型鏃の概形（単位: mm）および，作成した計算格子を示す.鏃の，半径による無次元長さは流線形鏃が11，椎型鏃が4.8となる。また，計算格子はξ方向を12本毎に，η方向を４本毎に描画している.初期条件は一様流とする。時間刻みは△t = 1.0×10−3に設定する.  Re数はRe = 1.0 × 10 4 , 1.2 × 10 4 , 1.4 × 10 4 , 1.6 × 10 4 , 1.8 × 10 4 , 2.0 × 10 4 とした．定常解と判定する基準は，すべての格子点において時間ステップnと次ステップn + 1での∪ ,の相対変化量の絶対値が1.0 × 10−4以下になった場合とした．０ ≤ t ≤ 1200の間でz = 50, 100, 150, 200, r − 1 = 0.2におけるUz, Uzの時刻歴はもに関してそれぞれ単調に減少/増加しながら一定値に漸近する挙動を示した。定常解に到達した時刻は，Re = 2.0 × 10 4において流線形鏃，椎型鏃ともにｔ〜 1200であった.図2に，流線形鏃の場合に対する，軸方向の速度成分U zをz = 0,50,100,150,200において示す. 同じ位置におけるU,の速度分布ではz = 0において表面に向かう弱い流れがあり，z = 50においてr − 1 − 0.2で表面から離れるごく弱い流れがある。しかし，z ≥ 100においてし,の値はO(10~³)以下の大きさとなる.図３に，椎型鏃の場合に対するU,の分布を図２と同じ位置において示す.図３のｚ＝０からわかるように，椎型鏃では鏃の付け根に逆流領域が生じる．可視化図（省略）からは，鏃表面で剥離し，シャフト表面で再付着する扁平な剥離泡がみてとれる．剥離泡の長さは，ここで調べたRe数の範囲ではほぼReに比例して長くなり，Re = 1.0 × 10 4において−0.8 ≤ z ≤ 0.6, Re = 2.0 × 10 4において−1.0 ≤ z ≤ 1.6が逆流領域となっている．U，の分布は，z = 0以外の位置において流線形鏃とほとんど類似の図となる．逆流領域が発生する影響と椎型鏃が流線形鏃よりも短いために，椎型鏃の境界層厚さは0 < z ≥ 60で流線形鏃のそれより薄くなるが，60 < zでは両者の差はほとんどなくなる．軸方向の圧力勾配は鏃の部分を除いてほぼ０である.
3　平行流近似に基づく線形安定性解析
基本流計算で得られた主流(U 2 , U 2 , 0 , P)に,微小撹乱( u z , u r , u o , p )を加えて，その時間発展を調べる.平行流近似ではひ,の軸方向依存性を無視し，ひ,はひ,に比べて十分に小さいとして無視する．平行流近似での主流は半径方向のみに依存するため，線形撹乱を軸方向と周回方向でスペクトル分解できる。速度と圧力の撹乱u, pは以下のように記述される.
ここで，αは実数であり撹乱の軸方向の波数，mは整数で周回方向の波数を表す。ωは複素周波数であり，虚部ω;は撹乱の成長率を表す。ω;の正負によって撹乱の振幅が時間発展するかどうかが決定される．ω； > 0のとき，撹乱が成長し，主流は線形不安定である．ω；＜０ならば線形安定な流れであり，撹乱は減衰する。また, ωの実部ω,は線形撹乱の角振動数を表す.式（7）をNavier-Stokes方程式の速度・圧力成分に加えて線形化すると以下の式を得る。
ここで，プライムはrに関する微分を表す．速度場の撹乱は境界条件（12）を満たす.
式（8）-（11）は半径方向rを独立変数とする撹乱の常微分方程式となり，ωに対する固有値問題に帰着される.固有値問題を数値的に解くために，ζ = 2(z/Re) 1/2 , σ = (r − 1)/ζの座標変換を導入する。 σ方向の離散化にはChebyshev選点法を使用し，線形撹乱方程式をベ  が成長する．特に，z = 20付近でm = 2のモードが線クトル行列形式で取り扱う．選点数Nに対し，方向形増幅し，上流側の線形不安定性に最も寄与する。また，Re = 2.0 × 10 4 （図４（f））ではm = 3のモードがの点とChebyshev点は以下の式の関係がある.
ôはσ方向のChebyshev点の分配を決定する実数で,シャフト付近での点の分配の粗密を決定する．本論文ではδ = 15とした。式（8）-（11）から，圧力項qを消去し，境界条件（12）を適用し，3 N − 6次の正方行列の固有値，固有ベクトルを計算する．固有値計算はMathWorks社の数値解析ソフトウェアMATLABのeig関数を使用したQR法で行う．本研究では選点数Nを64とした．また，周回方向波数mについてはTutty et al. 5)の先行研究の臨界Re数から，不安定性を生じる可能性があるm = 0, 1, 2, 3の場合を調べた．事前に，計算コードの検証のためにBlasius境界層の固有値を計算し，６桁の精度を確認した．
4　線形安定性解析結果
線形安定性の解析結果として，局所最大成長率を図4に示す.軸方向距離zを，縦軸は撹乱の成長率ω;を表す。また，図中の□，＋はそれぞれ流線形鏃，椎型鏃を装着した場合の計算結果で，ω;が最大となるαを探した．Re数が増加するにつれて局所成長率は大きくなる．不安定モードのm = 1（黒），2（青）ではRe数の増加につれて中立安定点が上流側へ移動していく. m = 3（緑）のモードはRe = 2.0 × 10 4（図4（f））でのみ，z = 20付近で成長する．その成長率は，m = 1, 2のモードよりも小さく，z = 45で再び線形安定となる.また，シャフト中間部を過ぎるとm = 3のモードのω iは再び増加する.鏃による局所成長率の違いはどのRe数でも，シャフト前半部で現れる。上流側ではいずれのモードも，流線形鏃を装着した場合の方が，わずかであるがより不安定となる．流線形鏃では椎型鏃に比べて境界層がより厚くなるためである。境界層の厚さの差が大きいz ＜ 60（図2,  3）の範囲では，成長率の違いが見られる．ｍ = 2のモードはz ≥ 70で，m = 3のモードはz ≥ 30で椎型鏃の場合の方が，成長率がわずかに大きくなる．解析したいずれのモードもz ≥ 70では鏃による不安定性の違いがほとんど見られなくなる．これは，基本流である境界層流れに対する鏃の影響がほとんど現れないためである.アーチェリー選手の放つ矢のRe数 (1.6×10 4 ≤ Re ≤ 2.0 × 10 4 ,図4 （d） - （f））ではm = 1,2,3のモードz = 20付近で線形増幅するが，m = 1,2のモードに比べるとその増幅率は小さい。 m = 1のモードは局所成長率がzとともに単調増加しており，解析したどのRe数でも最下流側で成長するモードである。したがって,今回調査したRe数全域で矢側面の境界層流れは，非軸対称モードに対して線形不安定な流れである．一方で， m = 0（赤）の軸対称撹乱は全く成長しない。軸対称モードm = 0はa → ∞では2次元撹乱に対応する．この事実は平板境界層の線形不安定性において2次元TS波が主要な役割を果たしていることと対照的である.図5に，流線形鏃の場合における局所最大成長率を与える波数α maxを示す．線形不安定となる場合のαmaxのみを示している．黒実線はm = 1を，青破線はm = 2を表し，それぞれの色の曲線で上からRe = 2.0, 1.8, ..., 1.0×10 4を表す.  m = 2はRe = 1.4×10 4まで示す．鏃による定性的な違いはないので流線形鏃の結果のみを示す. α maxは境界層の厚さによって変化する.下流に向かうほど境界層が厚くなるので，αmaxは小さくなる。また，Re数が増加すると，αmaxが増加する.これはRe数の増加につれて，境界層が薄くなるためである.  m = 1,2のモードを比較すると，m = 2のモードのα maxが小さくなる．後端での撹乱の波長λ = 2π/αは，シャフト半径の8～30倍程度となる．Re = 2.0×10 4の中立安定点でのm = 1, 2のモードの波長は，それぞれ半径の3.6倍，2.5倍程度となった。図6に，流線形鏃の場合のRe = 2.0 × 10 4における,成長率が最大となるz = 200(m = 1),66(m = 2)での固有関数を示す。横軸は半径方向距離r − 1を表し，シャフト表面を0としている．図6の(a)， (b)， (c)はそれぞれz, r, θ方向の固有関数である. それぞれの図で固有関数の各成分の絶対値を実線で表す。流線形鏃の場合,成長率が最大となる場合の軸方向波数αは，m = 1のモードでα = 0.76,  m = 2のモードではα = 1.22であった．また，鏃による違いはほとんど見られず，椎型鏃の場合はm = 1ではz = 200, α = 0.78,  m = 2ではz = 70,  α = 1.24となった。固有関数の各成分の半径方向の広がりを排除面積¹³ﾝから得られる排除厚さに相当するδ(z)と比較する.
式（15）よりδ(66) = 0.143であり，δ(200) = 0.230である．固有関数のピークはシャフト近傍に現れ，不安定モードは境界層内で成長する。モードによる違いはピークの場所に現れる．ｍ = 1のモードに比べて，ｍ = 2
のモードの固有関数はさらに表面付近に集中している.
5実験結果との比較と議論
線形安定性解析結果（図4）を，我々のグループの実酸結果（Miyazaki et al.3,4)と比較する。図７に，椎型鏃を装着したA/C/E矢のCD-Re数依存性を示す. 図7中の●はベアーシャフトで迎角0度の場合の風洞実験結果である．このCD値は，理論計算14,15)の結果に重なる層流値を示した．◆は矢羽をつけた場合の風洞実験結果で，CD値は矢羽に働く抗力分だけ増加したことを示す。迎角を0度にした場合，Re数増加に伴いCD値が滑らかに減少していることから層流状態を示す。+は迎角を0.75度付けた場合で，Re = 1.0 × 10 4からCD値が跳ね上がる. これはRe = 1.0 × 10 4を超えたところで境界層が乱流遷移し，抗力が増加したことを示す. ◇は国立スポーツ科学センター（JISS）での圧縮空気による飛翔実験結果である．Re = 1.0 × 10 4に設定して矢を飛ばした場合は迎角によらずCD値は層流を示した. 1.3 × 10 4 ≤ Re ≤ 1.8 × 10 4の範囲ではCD値の平均は層流値と乱流値の間をとり，標準偏差は大きくなった.このRe数領域ではCDの層流値と乱流値は混在し，迎角のばらつきによるものであると考えられている．一方，Re = 1.7 × 10 4に設定し，打ち出しに使用するノズルの先端を10度カットして迎角を制御した場合(☆) , CD値は唇流値をとる結果が得られている．×は競技者が放った矢の解析結果であり，CD値は乱流値を示している．これは矢の迎角が大きいこと，および弓から放たれた際の矢のたわみ振動が乱流遷移を誘発するためである。競技者の放つ矢のCD値のばらつきは非常に小さく，完全に乱流遷移していることを示唆している.図8に，流線形鏃を装着したA/C/E矢のCD-Re数依存性を示す.  MSBS風洞実験結果については，椎型鏃の場合の結果と一致した．一方，飛翔実験結果（△）は，椎型鏃の場合との差が見られた.  Re = 9.0 × 10³で矢を飛ばした場合は，CD値は層流値を示した．1.3 × 10 4 ≤ Re ≤ 2.0 × 10 4ではCDの平均値は層流値と乱流値の間をとり，標準偏差が大きくなる．この遷移領域は椎型鏃の場合に比べてより高Re数領域に広がることがわかる. Re = 2.4 × 10 4の椎型鏃の場合はCD値の平均は乱流値(CD = 2.6±0.03)をとり，標準偏差が小さくなっていることから，境界層は完全に乱流状態である．一方，流線形鏃の場合（C D = 2.2 ± 0.4）は標準偏差が椎型鏃の場合と比べて大きい．このことからRe = 2.4 × 10 4でも乱流遷移が完結していないと考えられる.以上より，迎角が非常に小さい場合は鏃の形状によらず，飛翔実験においても居流境界層が実現する．線形安定解析結果によればm = 1,2  (Re = 2.0 × 10 4ではm = 3も）のモードが成長するはずであるが乱流遷移をもたらすことはないと言える.線形不安定性領域が存在するにも関わらず，層流を観測した要因として次の２点が考えられる. (i)線形撹乱の成長率が小さいこと, (ii)使用したMSBS風洞の乱流強度が小さく(0.05 %)，境界層への撹乱の導入が弱かったこと，である．まず，群速度cg = dωr/dαを計算し，線形撹乱の波束の振る舞いを考察する．図５に示す波数αの撹乱の波束は速度c gで伝撒する。図9に,流線形鏃の矢の場合について，線形成長する撹乱モードごとにcgを示す。群速度は一様流速の0.3～0.4倍程度でRe数が大きくなるほど小さくなり，撹乱の滞留時間が長くなる．m = 1のモードの群速度は下流に向かうほど小さくなるが，負にならない．これはすべての撹乱波束が下流に向かって流れることを示す。特に中立安定点でも，撹乱が上流の安定領域に遡ることはない．したがって，撹乱は二次不安定16)を起こし乱流遷移をもたらす前に下流に流される可能性が高い.次にe^法に基づいて乱流遷移の可能性をより定量的に検討する．e^法は境界層の乱流遷移を線形安定性解析によって予測する半経験的な手法である177.時間増幅率からGaster変換18)により空間発展を評価する。ここでは式（16）を使用して線形安定性解析結果から包絡的なN値を計算する.
局所成長率と群速度をスプライン補問により関数近似し、不安定性の現れる範囲で式（16）を数値積分し，N値を計算する.図10に，2つの鏃に対して，モードごとに計算したN値のRe数依存性を示す.  Re < 1.6 × 10 4ではm = 1のモードのN値が大きくなり，Re > 1.6×10 4ではm = 2のモードがN値が大きくなる.  Re = 1.0×10 4における，m = 2のモードは成長しないので計算をしていない．椎型鏃の場合よりも流線形鏃の場合のN値が大きいのは上流側の線形不安定性の成長率がわずかながら大きいことによる. N値と乱流強度からMSBS風洞実験における乱流遷移の可能性を評価する。平板境界層の場合について， Mack¹9により乱流遷移に至るN値の評価式(17)が与えられている.
式（17）の一様流中の残留乱れTuへの依存性は，2次元TS波の線形成長を経由した乱流遷移に対するものである。円筒面に対する同様の評価式が見いだせないので，本研究では，参考として式（17）を使用して遷移の可能性を評価する。MSBS風洞の乱流強度は0.05 %であるからTu = 5.0×10−4となり，遷移位置のN値を計算するとN = 9.81となる．図10に示されるN値はいずれもこの値を下回る。さらに，Vinod and Govindarajan 201はm = 1, 2, 3の不安定モードが有限振幅（主流の2 %）に成長した状態の二次不安定性を調べた．彼らは円柱側面に発達する境界層は，平板境界層よりも不安定化しにくいことを示している．したがって，本研究でのN値は乱流遷移をもたらす程大きくないと判断される. MSBS風洞実験で迎角を0度に設定した場合や，飛翔実験において迎角を制御して矢を発射した場合にCpが層流値となったことは，撹乱が乱流遷移をもたらす程成長しなかったためと考えられる.その一方で，気流に対してわずかに迎角のついた矢の側面における境界層流れは乱流化することが，MSBS風洞実験と飛翔実験から分かる。図７と図８を比べると，流線形鏃を装着した場合の遷移領域は椎型鏃を装着した場合よりも高Re数領域に広がっている.図10に示すN値は流線形鏃を装着した場合の方が少し大きいので，この飛翔実験結果を説明することはできない．鏃形状の基本流への影響は原点（鏃の付け根）付近で顕著なので，矢先端（z = 1,2）における境界層流れの不安定性も調べる必要がある．ただし，椎型鏃の場合の基本流は剥離を伴うため平行流近似が適用できず，定量的に信頼できる結果を与えるためにはbiglobal²¹)安定性解析を行う必要がある。ここでは，平行流近似に基づいて線形不安定性を調べ，定性的な目安とする．図11にz = 1における解析結果を示す。実線は流線形鏃の場合の，鎖線は椎型鏃の場合の成長率のα依存性を表し，不安定性は椎型鏃の場合にのみ現れる．椎型鏃の場合の成長率のピークは上からm = 0（赤）, 1（黒）, 2 (青)，３（緑）のモードであり，流線形鏃の場合は上からm = 3, 2, 1,0のモードである．椎型鏃の場合の不安定性成長率はO(1)であり，N値への寄与も4.6程度に達する。その他のモードの不安定性成長率も同程度の値となる．流線形鏃の場合には，原点付近でも線形安定である．この先端付近における安定性特性の差が，鏃形状による乱流遷移領域の差につながっているようである.飛翔実験で待られたCD値はどちらの鏃の場合でも層流値と乱流値に二値化していることから，気流に対するわずかな迎角が剥離に影響を与えたり境界層内にストリーク構造を誘起して，シャフト下流部における線形不安定性を経由しない乱流遷移をもたらす可能性もある.迎角がついた場合，主にm = 1のモードの撹乱が矢の先端から導入されるはずである。この撹乱は線形過程では，矢の後方でのみ成長する。しかし，飛翔実験でCD値が二値化することは，矢の先端から乱流遷移が引き起こされていることを意味する．bypass遷移²²のように線形不安定性を経ない遷移の可能性も含めて，より詳細
6結言
2種類の鏃をつけた非回転アーチェリー矢の軸対称流を計算し，平行流近似に基づく線形安定性解析を実際のアーチェリー競技で取りうるRe数領域で行った. 1.0 × 10 4 ≤ Re ≤ 2.0 × 10 4の範囲ではm = 1のモードの線形不安定性が現れ，1.2 × 10 4 ≤ Re ≤ 2.0 × 10 4の範囲でm = 2のモードの線形不安定性が現れた。m = 3のモードの不安定性はRe = 2.0 × 10 4で現れた。m = 0のモードの不安定性は，解析したRe数の範囲では現れなかった．境界層や線形安定性の鏃による違いはz ～ 60程度の箇所まで現れ，わずかであるが流線形鏃を装着した場合の方が不安定性は強くなった．包絡的にN値を計算した結果、撹乱は乱流遷移をもたらすほど成長しないことが示された．このことは迎角が0度の場合のMSBS風洞実験や，選手領域のRe数での飛翔実験で迎角を制御した場合に層流を観測したことと整合する。その一方で，矢側面の境界層流れの乱流遷移は線形不安定性を経由しないbypass遷移²²である可能性がある.椎型鏃の継ぎ目近傍では剥離泡による強い線形不安定性が現れた．この領域での平行流近似に基づく解析は不適切であり，biglobal 21な線形不安定性解析を行うことが必要である．矢羽による回転の影響を考慮した詳細な安定性解析を行い，迎角がついた矢側面の境界層流れの乱流遷移過程を解明しなければならない。実験データの整理に関して，電気通信大学の安藤，芹野両君に協力をいただき，ここに感謝します．本研究は科学研究費補助金基盤研究(C) (一般) 18K03943 (2018 ～ 2020年度)の支援の下に行われたものである. 1) Park, J. L., Hodge, M. R., Al-Mulla, S., Sherry, M., and Sheridan, J.: Air flow around the ponit of an arrow, J. Sports Engineering and Technology, 227 (2013) 64– 69. 2) 鈴木一史,桝井和典,向山桂太,宮器武,澤田秀夫:矢の空力特性-境界層遷移に対する先端形状の影響-,ながれ, 29 (2010) 287–296. 3)  Miyazaki, T., Mukaiyama, K., Komori, Y., Okawa, K., and Taguchi, S.: Aerodynamic properties of an archery arrow, Sports Engineering, 16 (2013) 43–54. 4) Miyazaki, T., Matsumoto, T., Ando, R., Ortiz, J., and Sugiura, H.: Indeterminacy of drag exerted on an ar- row in free flight: arrow attitude and laminar-turbulent transition, Eur. J. Phys., 38:064001 (2017). 5) Tutty, O. R., Price, W. G., and Parsons, A. T.: Bound- ary layer flow on a long thin cylinder, Phys. Fluids, 14 (2002) 628–637. 6) Muralidhar, S. D., Pier, B., and Scott, J. F.: Instability of flow around a rotating, semi-infinite cylinder, Phys. Rev. Fluids, 1:053602 (2016). 7) Muralidhar, S. D., Pier, B., Scott, J. F., and Govin- darajan, R.:  Flow  around  a  rotating,  semi-infinite cylinder in an axial stream, Proc. R. Soc. London A, 472:20150850 (2016). 8) Herrada, M. A., Pino, C. D., and Fernandez-Feria, R.: Stability of the boundary layer flow on a long thin ro- tating cylinder, Phys. Fluids , 20:034105 (2008). 9) Iwatsu, R.:  Application of fourth-order explicit schemes to the computation of axisymmetric swirling flows of Boussinesq fluid in containers with arbitrary cross sections, Num. Heat Trans. B , 56 (2009) 191– 210. 10) Carpenter, M.K. and Kennedy, C.A.: Fourth-order 2N-storage Runge-Kutta schemes, NASA TM , 053602 (1994). 11) Eiseman, P. R.: Mesh generation using algebraic tech- niques, NASA CP , 2166 (1980) Numerical Grid Gen- eration Techniques , 73–120. 12)  Vinokur, M.: On one-dimensional stretching functions for finite-difference calculations, J. Comput. Phys., 50 (1983) 215–234. 13) Tutty, O. R.: Flow along a long thin cylinder, J. Fluid. Mech., 602 (2008) 1–37. 14) Seban, R. and A., Bond, R.: Skin-friction and heat- transfer characteristics of a laminar boundary layer on a cylinder in axial compressible flow, J. Aerosp. Sci., 18 (1951) 671–675. 15) Kelly, H. R.: A note on the laminar boundary layer on a cylinder in axial compressible flow, J. Aerosp. Sci., 21 (1954) 634–637. 16)  Herbert, T.: Secondary Instability of Boundary Layers, Ann. Rev. Fluid Mech. , 20 (1988) 487–526. 17) 徳川直子,高木正平,上田良稲,井門敦志:矩形翼境界層の自然遷移に対する外乱の影響,ながれ, 24 (2005) 629–639. 18) Gaster, M.: A note on the relation between temporally- increasing and spatially-increasing disturbances in hy- drodynamic stability, J. Fluid Mech. , 14 (1962) 222– 224. 19)  Mack, L. M.: Transition Prediction and Linear Stabil- ity Theory, AGARD Conf. Proc. , 224 (1977) 1-1-1-22. 20) Vinod, N. and Govindarajan, R.: Secondary instabil- ity in incompressible axisymmetric boundary layers: Effect of transverse curvature, J. Fluids Engineering, Transactions of the ASME , 134:024503 (2012). 21) Theofilis, V.:  Advances in global linear instability analysis of nonparallel and three-dimensional flows, Prog. Aerosp. Sci. , 39 (2003) 249–315. 22) Schmid, P. J. and Henningson D. S.:  Stability and Transition in Shear Flows , Applied Mathematical Sci- ences No. 142 (Springer, 2001) 402.