惑星大気研究会オンラインセミナー，神戸大学惑星科学研究センター，2010年8月26日

高分解能数値実験による地球ダストデビルの強化・維持過程の解析
竹見哲也
京都大学防災研究所
気象・水象災害研究部門

内容• ダストデビル（塵旋風）とは– 形態・構造– 形成メカニズム– 大気境界層の特徴– 環境場の特徴• ダストデビルの高分解能数値シミュレーション– 渦の強化・維持のメカニズム– 一般風の影響

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

ダストデビル（塵旋風）とは• ダストデビル：晴天時に平坦地での大気境界層内に発生する鉛直軸回りの渦– 砂漠など乾燥地で頻発（平坦地形が好条件）– 強い地面加熱による対流混合層の発達が重要（絶対不安定な接地層の存在）– 強風下では発生しにくい– 湿潤過程は関与しない

竜巻：積乱雲に伴う渦(Balme and Greeley 2006)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University古くから存在は知られていたロシアの探検家によるツァイダム盆地のダストストームの観察（Przheval’skii 1883，絵はRoborovskii）(Nalivkin 1983, Fig. 216)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

北極圏でも発生しうる
エルズミーア島（カナダ）でのダストデビル(Nalivkin 1983, Fig. 219)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

日本におけるダストデビルの発生状況
発生時期(件) (a)
98
7
654
32
10
2月

(件) (b)

３月

４月

５月

６月

７月

８月

発生時間帯

43

強さは藤田スケールでF0程度

21（気象庁「竜巻等の突風データベース」より）

09時 10時 11時 12時

13時 14時

15時

16時

ダストの舞い上がり過程

長い浮遊（<20μm）

風短い浮遊風の乱れ（20-70μm）

跳躍（70-100μm）

1.5 m

地面転動（>500μm）

跳躍（70-500μm）

ダストデビルに関する研究• 観測：Sinclair (1969, 1973), Carroll and Ryan (1970), Kaimaland Businger (1970), Hess et al. (1988), Hess and Spillane(1990), Bluestein et al. (2004), Markowski and Hannon (2006)• 高分解能モデリング– Kanak et al. (2000), Kanak (2005), Ito et al. (2010)• 対流境界層のLESにより渦の発生・強化過程について
研究：stretching（引き伸ばし）, tilting（立ち上げ）効果• 環境条件– 一般風の影響：Toigo et al. (2003), 田中ら (2008)• Toigo et al.では一般風の影響は小さい vs 田中らでは
弱風が好条件– 安定度の影響：Hess et al. (1988), Hess and Spillane (1990)• 不安定度指数 w*/u* > 5 を満たす時に塵旋風が発生– 地面粗度の影響：Lyons et al. (2008)

ダストデビルの観測的特徴
接線風速・水平風速・鉛直風速の特徴
接線風速

接線風速と大きさ（直径）の関係

上昇速度

水平風速

気圧偏差

(Balme and Greeley 2006)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

火星のダストデビル(Balme and Greeley 2006)（NASA‘s Mars Exploration Rover Spiritによる）
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

地球と火星のダストデビル
直径(m)

高さ(m)

接線風速
(m/s)

鉛直速度(m/s)

回転の向き

気温偏差(K)

気圧偏差(Pa)地球

1-150(5-30)

1-2400(2501500)

5-203-15(4)ランダム

2-8(3-4)

250-450

火星

5-1750(250)<8000(20006000)

2-93

N/A

ランダム

0.5-6

1.1-4.6(2-3)（カッコ内の数値は平均値または典型値）（Kanak 2006)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

ドップラーレーダーによる観測
Wバンド（3 mm波長）ドップラーレーダ観測：•方位角方向分解能：3-5 m•動径方向分解能： 30 m

(Bluestein et al. 2004)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

ダストデビルの多重発生

水平格子幅2.75 mでのLES(Kanak 2008)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

砂漠域の対流混合層
5月平均の仮温位
1991年5月の中国ゴビ砂漠南部のMinqinにおける20時BSTの鉛直プロファイル

上位5日分の仮温位

上位5日分の水蒸気混合比(Takemi 1999)(Takemi and Satomura 2000)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

ダストデビル時の摩擦速度と安定度
ダストデビル発生−h/ Lw∗ / u ∗550

h: 境界層厚さ
L: Monin-Obukov長(Balme et al. 2003)(Hess et al. 1988)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

晴天時の境界層日変化とダスト輸送
地上気象変化気温

ADEC観測時にライダーで観測されたゴビ砂漠南部における晴天時のダストの鉛直分布の変化。
期間：2002年4月11日－12日

風速（30分平均値・最大値）

相対湿度(Yasui et al. 2005)
4月11日 4月12日 4月13日(Takemi et al. 2005)（背景写真：情報通信研究機構 安井氏提供）

境界層過程とダスト輸送(Koch and Renno 2005)

Dust mobilization (µg/m2/s)

西サハラ

東サハラ

Dust surface conc (µg/kg)

Diurnal contribution

サヘル アラビア半島 オースト
ラリア

Synoptic contribution

東アジア(Luo et al. 2004)

対流混合層の構造とダストデビルの発生機構(Hess and Spillane 1990)(Kanak 2000)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University一般風の影響高度別の渦度の最大値

U=0.5

U=5.0

(Toigo et al. 2003)

U=15.0(Ito et al. 2010)

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

本研究で調べること• 対流混合層中に発生した多数の渦のうち強い渦に発展するものは数が限られている。強い渦に強化され維持されるのはなぜか？
• 一般風があるとダストデビルの発生にもなんらかの影響をもたらす：– 弱風時が好条件（観測、数値実験：Ito et al. 2010）– 弱風時のみならず強風時でも発生する（数値実験：Toigoet al. 2003）
高分解能数値実験により• 渦の強化・維持過程• 一般風の影響を調べる
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

既往の高分解能モデル研究の諸元

B=Boussinesq, SC=supercompressible, FC=fully compressible
FS=free slip, NS=no slip

(Kanak et al. 2006)

格子点数

領域(km)

格子間隔(m)

境界条件

一般風(m/s)B

90 x 90 x 60
200 x 200 x 130

4.5 x 4.5 x 34 x 4 x 2.650
20NS0.5, 5.0, 15.00.5

B

210 x 210 x 1304.2 x 4.2 x 2.6

20NS

0.5, 5.0, 15.0

方程式系
田中ら
(2008)
Ito et al.
(2010)

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

数値実験• 数値モデル：Weather Research and Forecasting (WRF) model– Large-eddy simulation (LES)モードで利用– 乱流モデルはSmagorinsky-Lilly• 問題設定の簡略化：– 水平一様な基準状態– 境界層スケールの計算領域– 側面境界は周期条件– 物理過程は乱流以外は考慮しない• 実験1：水平格子幅3 m；一般風なし– 渦の強化・維持過程の解析• 実験2：水平格子幅10 m；一般風U0=0, 5, 10, 15 (m/s)– 一般風に対する感度の解析
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

1.5km

実験1の設定1km領域

水平1km・鉛直1.5km

水平解像度3m

鉛直解像度

下層3mで上層は粗い
Stretch grid

境界条件

水平：周期境界
下端：バルク法

地表面の熱Flux 0.24(K m/s)≒290W/ m 2
初期温位(左図) 高度800mまで300K
以高 10K/kmで上昇

図：初期・3600s後の鉛直温位分布

乱流モデル

Smagorinsky-Lillyモデル

一般風なし

雲物理なし

対流混合層の構造
時刻5100 sでの高度50 mにおける鉛直速度の水平断面

乱流統計量
Heat flux(m/s)

VarianceSkewness

対流混合層内の鉛直渦対流セル

Y(m)

(m/s)

X(m)
図: t=5100sの高度50mにおける鉛直速度(色付)と鉛直渦度(等値線：絶対値0.20(1/s)以上)の水平断面図。

圧力偏差と最大風速半径10秒間隔で180time step間（計算開始後60-90分）に発生していた渦を抽出した結果

図：横軸が中心からの距離で縦軸が圧力偏差を表す。

図：横軸が最大風速半径を格子幅で、縦軸が渦の発生数を総数(225)で規格化したもの

・発生していた渦の大多数は弱くて小さい渦。ごくまれに強大な渦が発生。

渦の融合5050s

5060s

5070s

5080s

高度1.5mの鉛直流（色付）と渦度（等値線：間隔は0.15 s-1）と風速ベクトル

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University閾値|p’| > 10 Pa−1|ζ| > 0.15 s

渦の強度変化閾値を超える強い渦の中心における(上)圧力偏差と(下)鉛直渦度の時間発展。
各色線が一つ一つの渦を表し, 矢印は渦の融合が起こった時間を示す。

渦の強化過程：立ち上げと引き伸ばしω = (ξ ,η , ζ ) u = (u, v, w)z(m) ∂u ∂v Dζ∂w∂w≈ξ+η− ζ  + ∂x∂yDt
 ∂x ∂y 

鉛直流（色付）と立ち上がり項（等値線：間隔は0.03 s-2 ）x(m)
渦度（影）と水平発散（等値線：間隔0.1 s-1）の鉛直断面1.5km

実験22km領域

水平2km・鉛直1.5km

水平解像度

10m

鉛直解像度

下層10mで上層は粗い
Stretch grid

境界条件

水平：周期境界
下端：バルク法(熱・運動量)

地表面温度

320K (Ryan 1972)

初期温位

高度800mまで300K
以高20K/kmで上昇

一般流 U00; 5; 10; 15 m/s

乱流モデル

Smagorinsky-Lillyモデル

雲物理なし

図： t=3600sにおける温位の鉛直分布

平均プロファイル

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University一般風の違いによる対流構造の変化
U0=0U0=5U0=10U0=15

高度50 mでの鉛直速度（カラー）および渦度の鉛直成分（実線）

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

気圧偏差と渦度の最大値の違い(Toigo et al. 2003)

渦の構造（U0=5）

T=4800sZ = 15my(m)

z(m)

x(m)
鉛直流（色付）と圧力偏差（10Paごと）

x(m)
Y=1180の鉛直断面（緑: v 1(m/s)ごと)

渦の構造（U0=15）

T=4800s

Z=15my(m)

z(m)

x(m)
鉛直流（色付）と圧力偏差（10Paごと）

x(m)
Y=1930の鉛直断面（緑線はv: 1(m/s)ごと)

一般風の違いによる渦の強度変化

閾値|p’| > 25 Pa −1s|ζ| > 0.15(a) U0 = 0(b) U0 = 5(c) U0 = 10(d) U0 = 15閾値を超える強い渦の中心における圧力偏差の時間発展の図。
各色線が一つ一つの渦を示している。

一般風ごとに違いが見られた要因考えられる要因
・組織的な（セル状）対流構造が壊れてくる
風速が強くなると・・・

・風速変動が大きくなり、拡散性が強まる

・地表面からの熱フラックスが増大する

Negativerole
Positiverole

安定度との関係(w* = gH 0 hθ)

1−1 3：対流速度スケール
0u* = C D | U |k =1 ：摩擦速度

まとめと今後の展開• 渦の強化・維持においては、複数の渦の融合により強い渦が形成されること、立ち上げ・引き伸ばしによる強い渦の維持、が重要なメカニズム• 一般風の影響は、弱から中程度の風速の場合に最も多数・最も強い渦が形成される、ということに表れる• 渦の強さは混合層の深さによってどう変化するか？
• 空気密度の違いはどう影響するか？

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University参考論文• Ohno, H., and T. Takemi, 2010a: Mechanisms for intensification andmaintenance of numerically simulated dust devils. Atmospheric ScienceLetter, Vol. 11, pp. 27-32., doi: 10.1002/asl.24.・・・渦の強化・維持過程• Ohno, H., and T. Takemi, 2010b: Numerical study for the effects of mean
wind on the intensity and evolution of dust devils. Scientific Online Letterson the Atmosphere, Vol. 6A, pp. 5-8, doi:10.2151/sola.6A-002.・・・一般風の影響• 大野 洋，竹見哲也，2009：数値計算パラメータによって変化する塵旋風の構造の感度実験．京都大学防災研究所年報，第52号B, pp. 465-470.
(available online at http://www.dpri.kyoto-u.ac.jp/nenpo/no52/ronbunB/
a52b0p49.pdf )・・・解像度依存性

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

参考：積乱雲による渦の形成
風の鉛直シアー（dvh/dz）が水平方向に軸を持った渦度場を表現している。
簡単のため、風向変化がない鉛直シアーの場を考えると、立ち上げ効果による渦度変化は次式で表せる：dς dU ∂w=
dtdz ∂y(5)

この式によれば、•上昇流の南側（∂w/∂y>0）では正の渦度（低気圧性循環）が生成され、•北側（∂w/∂y<0）では負の渦度（高気圧性循環）が生成される。
さらに、立ち上げられた渦度は上昇流による引き伸ばし効果により強化される。
(Klemp 1987, cf. Rotunno 1981)
Tetsuya Takemi, Disaster Prevention Research Institute, Kyoto University

CBL（対流混合層）T = 5100sz = 300mz = 200mz = 100mz = 50mz = 10m
X(m)渦の選出の方法
閾値
圧力偏差・・・ 10(Pa)
渦度の絶対値 ・・・ 0.15(1/s)(i-1, j-1)

(i, j+1)

1)閾値を超え、なお且つ周囲よりも圧力偏差が大きい地点を仮の中心とした。
2)格子間を40等分し、ランキン渦の圧力分布からのずれを最小二乗法で評 (i-1, j)価することにより真の中心を決定する。
3)各格子点の風速を動径方向と回転方向に分解して、周囲の中から最大風速半径を決める。
4)これらの操作を180time step繰り返す。

Realcenter(i+1, j)20 lines(i-1, j-1)

(i+1, j+1)

Temporalcenter(i, j-1)
20 lines(i+1, j-1)

渦度収支（U0=5）

∂w
Dζ
∂
∂
R− (∇ ⋅ u )ζ + (∇T × ∇P ) + K∆ζ= (η + ξ ) w + ζ∂x∂y∂zP
Dt
Total

Tilting

Stretching Divergence

Solenoidal

渦中心を中央にして、9格子×9格子の領域で各項を平均した時系列

Diffusion

渦度収支（U0=15）
Z=5m

Z=15m

Z=33m

Z=50m

水平断面：セル状構造（高度約50m）
20m

15m

5m

3m

10m

色付：鉛直流等値線：圧力偏差 5(Pa),
渦度 0.15(1/s)
以上の領域

風速分布
横軸：中心からの距離20m

15m

5m

3m

縦軸：速度の回転成分

10m

圧力分布
横軸：中心からの距離

縦軸：圧力偏差

20m

15m

5m

3m

10m

最大風速半径の頻度分布
横軸：最大風速半径/格子幅

縦軸：全渦数に対する割合

TOTAL:31

TOTAL:86

TOTAL:512

TOTAL:1157

TOTAL:218