数理解析研究所講究録
1184 巻 2001 年 41-4941

植物の成長を考慮した砂丘の模型
西森拓 (Hiraku Nishimori)
田仲浩久 (Hirohisa Tanaka)
大阪府立大学工学研究科

Graduate School of Engineerin , Osaka Prefect$g$

1$\mathrm{u}\mathrm{r}eU\mathrm{n}i$versity

はじめに砂丘の形状 (図 1 参照) や運動を決める因子として、 様々なものが知られている。 例えば、 それぞれの砂丘地域における砂の量、 卓越風の方向の定常性、 あるいは砂丘表面を覆う植物の量などである。 その中でも、 特に、 砂丘の運動と植物の成長の関係には未だ解明されていない点が多い。 これは、 これまで蓄積されてきた個別の観測事実 [1, 2, 3, 4, 5, 6] と理論的な考察とを対応づける方法論が十分発達していないことによる。 これまでの研究の中で、 例外的に、 個別でない般論を提案したものとして、 Hack の先駆的な報告 [2] が挙げられよう。 その中でHack は、 それぞれの砂丘地域での卓越風の状態や砂の量、 砂丘表面における植物の被覆率とその地域で主に形成される砂丘の形状の関係を簡単な相図として表した。 しかしながら、 その正当性や適用範囲についての考察は十分ではない。 今

著者らがこれまで開発して来た、 植生の伴わない砂丘の動力学の計算模型を植生のある場合に拡張し、 Hack の相図およびその他の観測報告などと比較を行う。 以下、 まず、 模型を簡単に紹介し、 その後その結果と観測事実の比較検討を$\text{回_{、}}$行う。

2模型
計算には 2 次元の格子模型を用いる。 場の変数として、 次の二つを変数を用いる。

i) 局所的な砂表面の高さ $h(i, j, n)$ii) 局所的な植物密度 ($c$$i,$ $j$, n)(砂丘表面が植物で覆われている割合)である。 これらの変数は水平方向の $N\mathrm{x}N(N=100)$ 格子中の各セルに割り当てられる。 ここで $\{(i, j)|1\leq i, j\leq N\}$ と $\{n|0\leq n\}$ は、 それぞれ、 位置と時間

42

横列型砂丘

バルハン砂丘

放物線型砂丘$arrow$

風向き

図 1: 典型的な形状の砂丘。 植生が存在するとき形成される砂丘でもっとも特徴なのものは、 放物線型砂丘である。 植生に覆われていない砂丘と違い、 両翼が風上側に形成されその上に植物が成長する。 両翼の間では風による侵食が進行し窪地が形成され、 その砂が前方の尾根に移動し砂丘は前進する。 一旦緑地化された砂地の土地が再び砂漠化される過程においてこのような放物線型砂丘が形成される場合がある。
を表す。 ただし、 ( , のと名付けられた各セルは、 –個–個の砂粒の大きさ (直径$i$0.lmm のスケール) に比べて十分広く、 –方、 砂丘の大きさ (–辺以上) に比べて十分狭い領域を代表している。 また各単位時間ステップ は、 それぞれの砂粒の運動の素過程であるサルテーション (風の方向へのジャンプ) の持続時間 (1秒以下) に比べて十分長いものとする。 卓越風は–方向 ( の正の方向) に–定の強さで吹くものとする。
二つの場の変数 $h(i, j, n)$ と $c(i, j, n)$ は後述される抑制因子によって互いに影響を与える。 以上の要素に加えてその他の多くの要素が系全体のダイナミクスに関与してくるが、 今回は植生の入った砂丘の動力学のエッセンスを少しでも簡単な表式で表したいという意図から、 上述の要素に限定したミニマムな動$10\mathrm{m}$$n$$i$$a_{\beta}(i, j, n)$

力学だけで観測事実をどこまで再現できるかを検証していく。 このような、 ミニマム模型による砂丘のシミュレーションは植生のない場合について、 著者たちのグループ [8] とたちのグループによって行われて来たが、 植生のある場合についても先駆的な計算が Castro によって行われている [10]。 Castro の模型には、 植生がない限りうまく砂丘が形成されないという欠点はあるものの、 植生がある場合の系のダイナミクスはうまく単純化されており、 例えば、 植物の平均密度と隣合う横列砂丘間の間隔に関する観測事実をよく表している。 しかしながら、 この論文に続く理論面からの研究は非常に少ない。 言い替えれば、 理論模型を使って植生と砂丘の相互発展を研究するための方法論が確立されていない。
ともかく、 今回は著者の従来の砂丘模型を元にしてこれを植生のある場合に拡張する。植生がない場合から植生のある場合まで–貫した議論を行うことで、 砂$\mathrm{W}\mathrm{e}\mathrm{r}\mathrm{n}\mathrm{e}\mathrm{r}[7]$丘の運動に与える植生の効果を明確にしていく。

43$\mathrm{c}\mathrm{r}\mathrm{i}_{-}\mathrm{i}_{-\mathrm{n}}+1)$

図 2: a)(局所的) 植物密度の時間発展のルールを表す写像の模型:$c(i, j, n)arrow$) がないとき, 植物密度は飽和(I) もし表面の高さの時間変化:b(i,まで単調に増加する。 (II) 一定以上の速さで表面の高さが上昇または下密度降すると表面の植物は埋没したりあるいは根の露出によってダメージを被り植物$c(i, j, n+1)$$j,$ $n$。

$c_{\max}$

密度は最低で$c_{\min}$まで現象する。

1. まず局所的な植物密度 $c(i, j, n)$ の時間発展則として, 次のような観測事実を考慮する。 すなわち、 砂の移動がないとき、 植物密度は時間とともに成長し、 やがて、 ある飽和値に達する。 逆に、 砂の移動によって砂丘の表面が急に高くなった
り、 あるいは急に下がることで植物の埋没や根の露出が起こり、 これによって植物密度は極端に減る。 これらを反映させたミニマム模型を構成するために、 植物の成長に関して、 次のような区間線形写像を考える (図 )$2\mathrm{a}$$c(i, j, n+1)$$=$$A(c(i, j, n)-b(i, j, n))+c_{mim}$$c(i, j, n+1)$$=$$c_{\max}$$((c_{\max}-c_{\min})/A+b<c)$$c(i, j, n+1)$$=$$c_{mim}$$(c<b)$$(b\leq c\leq(c_{\max\min}-C)/A+b)$

ここで $b\equiv b(i, j, n)\equiv|h(i, j, n)-h(i, j, n-1)|,$ $c\equiv c(i, j, n)$ 。また、 $A$ は植物の成長率を決める定数。

2. 砂丘表面の高さ $h(i, j, n)$ の時間発展については, 離散化された砂の保存則,$h(i, j, n+1)-h(i, j, n)=Q_{in}(i, j, n)-Qout(i, j, n)$${ }$(1)が成り立つ。 ここで, $Q_{in}(i, j, n)$ は時刻 と $n+1$ の間にセル $(i, j)$ に飛び込んで来る砂の総量。 –方 $Q_{out}(i, j, n)$ は同じ時間帯にセル ( , のから出て行く砂の$n$$i$総量である。 $Q_{in}(i, j, n)$ と $Q_{out}(i, j, n)$ には、 下で説明する 2 種類の砂の移動形態、 i) サルテーション ii) クリープの両者による移動を含んでいる。

[サルテーション]:サルテーションとは卓越風の揚力および他の砂の衝突によって引き起こされる砂のジャンプが連鎖的に起こる過程である。 サルテーションは風44

図 2: b) 式 2 で表された、 サルテーションによる砂の流れ

風方向への) 局所的傾き$(=q_{S}a\iota)$。砂丘の (卓越の関数として表される。 これは頂上付近での$\nablah<0$
非常に強い風と風下側での風速の急激な低下という観測事実 $[8, 11]$ を定性的に再現したものである。 c) 砂の流れに対する抑制因子 は植物密度 の関数として与えられる (式 (3))。ある閾値
以上の植物密度のもとで、 砂の流れは極端に減少する。
$\nabla h$$\nabla\sim 0$

$a$

$c$

$c_{\mathrm{c}r}$

による砂の移動形態の中で特に重要なものである。 今、 サルテーションによって(的) にある量 $m$ の砂が時刻 とに $n+1$ の間に距離 だけ風下までジャンプしたものとする。 このときサルテーションによる流量 qsal $(i, j, n)\equiv l\mathrm{x}m$ を$n$$l$

qsal $(i,j, n)=q_{\mathrm{o}(an}th(\nabla h(i, j, n))+1.0)(tanh(-\nabla h(i, j, n))+1.0+\alpha)$ .

(2)と模型化する。 ここで $\nabla_{i}h(i, j, n)\equiv h(i, j, n)-h(i-1, j, n)$ を意味する。 また、は風上側における砂の流量を特徴づける定数である。上の、一見奇妙な表式は、次のような卓越風の性質 $[8, 11]($ 図を反映したものである。 すなわち、 1.
風上ではほぼ–定の風が吹く。 2. 頂上付近で風の鋭いピークがある、 3. 頂上付近$\alpha$$2(\mathrm{b}))_{\text{、}}$

における流れの剥離により、 風下面では風の力が極端に衰える。

[クリープ]: クリープとは表面に沿った転がりや滑動による砂の移動形態である。
クリープ $q_{creep}(i, j, n)$ は風方向への–定の運動に加えて、 砂表面の起伏によって影響を受ける。
$q_{C\Gamma}eep(i, j, n)=q_{creep}^{0}-D\nabla h(i, j, n)$(3)

ここで、 $q_{creep}0$ は砂丘表面の傾斜によらず、 風下向かって砂が–様に移動することを表す定数、は重力による、 一様な流れへの補正。 これには砂丘風下の急斜面での砂のなだれも含まれる。
$- \mathrm{D}\nabla h(i, j, n)$

3. さて、 植物の成長によって砂の運動がどのように影響を受けるのか。 重大な効果として、 砂丘表面における植物の被覆率が上昇するにつれて表面付近の卓越風の応力が急激に低下してくることが挙げられよう。 結果として、 砂の移動量が急激に低下する。 言い替えると、 植物の密度の増加にともなって、 砂の動きが抑制45される。砂の動きが急激に抑制され始める限界植物密度を$a_{\beta}(i, j, n)=1+\beta(tanh(c_{c}r-C)-1)$ .

$c_{cr}$として、抑制因子を(4)$(0\leq\beta\leq 1/2)$と表す (図) 。ここで は最大の抑制効果を決定づける定数である。 の値は、サルテーションとクリープとで異なった値をとる。 上の抑制因子を考慮して、 サ$2(\mathrm{c})$$\beta$$\beta$レテーション/クリープによる砂の流量を表すと、 $q(i, j, n)=a_{\beta}(i, j, n)2/q(i, j, n)$となる。 ここで $q(i, j, n)$ と $q’(i, j, n)$ は, それぞれ、 植生がある場合と無い場合の$)$砂の流量である。

以上が砂丘と植物の成長を表現する基本ルールである。

3

数値計算

Hack による相図をシミュレーションで再現することをまず第– の目標として次の 2 つの制御パラメータを導入する。
i) 風上の砂源における砂の量。具体的には、最も風上に位置する 2 列のセル:{(i, ) $|i=$$1\leq j\leq N\}$ における平均の砂高 $<h_{SO}urCe>$ である。 もっと詳しく言えば、$1,2$$j$$2<h_{SO}urCe>$と$0$の間の–様乱数を各セルに 1 ステップ毎に割り当てた。

ii) 卓越風の強さ。 これと、 サルテーションによる砂の移動量の間には互いに単調を卓越風の強さを指標として与えた。
増加の関係があるので、 式 (2) におけるまた、 Hack の相図では植物の密度も相図の軸の–つとなっているが、 これっは、上の二つの制御変数を決定すれば、 従属的に決定されるので今回の模型では制御$q0$変数としては取り扱わない。
シミュレーションでは (上で定義した砂源以外の場所で) 小さなゆらぎを除いて、 ほぼ平らな初期地形から始めて、 系にどのような形状の砂丘が出現するか観察した。境界条件は卓越風に垂直方向には周期境界条件を適用し、風下側の境界は自由境界とした。 また、実際の砂丘地域に対応して、卓越風によって動き得る乾燥砂層の厚みを有限とし、それ以下は硬い層もしくは水分を含んだ層と見なす。具体的にはある–定の高度以下 (ここでは、砂源を除いた初期の高さ $h(i,$ $n+1)=0$ )に砂の面が堀下げられることを禁止した。
$j,$4

結果

シミュレーションで得られた結果は次の通りである。 まず、 卓越風の風力がごく弱いときは、 砂源における砂の供給量の大きさによらず、 明確な形状の砂丘は現われない (図) 。卓越風の風力が–定の強さを超えた場合、 以下に記すよう$3(\mathrm{a})$に大きく分けて 2 種類の砂丘群が系の中に形成された。

46

図 3: 制御変数空間 (砂源での砂の量-風の強さ、 図 4 参照) 中のいくつかの点において現われる典型的な砂丘のパターン。 各図の左側は植物密度の空間分布 (より濃いところはより植物密度が高く白いところは植物は生えていない) を表し、 右側は砂表面の高さ (より濃いところはより高い場所) を表している。砂は格子の最も左の 2 列から供給され、 風は左から右方向に–様に吹いている。 典型的には次の 4 つのパターンがある。 (a):卓越風の風力が弱いとき、 明確な形状の砂丘は形成されない。 同様のケースは砂源での砂の供給が少ない場合にも当てはまる。 (b):風力が強く、 かつ、 砂源での砂の供給が非常に豊富な場合、 系には、主に横列砂丘ができる。 (C):風力と砂源での砂の量が (a) と (b) の中間であるときに放物線型砂丘が系に発生する。 その両翼はバルハン砂丘とは反対に風上側にのび、 両翼には植物が生えている。と類似の環境で砂の供給量が (c) より小量の場合延伸放物線型砂丘が形成される。 形はヘアピンに近く、 各翼だけに着目するならば線形砂丘に類似している。
$100\cross 1\mathrm{o}\mathrm{o}$$(\mathrm{d}):(\mathrm{c})$47i) 比較的強い卓越風 および十分に多量の砂の供給があるとき横列砂丘、 バルハン砂丘、 もしくはその両者が系の中で形成された。 (図).
$.3(\mathrm{b})$ii)–方、 砂の供給量、 もしくは卓越風の強さが i) の場合より少し劣る場合放物線型砂丘が系の中で多く形成される (図 3(c)(d))。放物線型砂丘は、 その両翼を風上側に向け、 両翼の上には植物が生えている。 両翼に挟まれた凹部では、 低い位置にある硬い地層 (もしくは水分を含む層) に至るまで砂の侵食が行われる。 侵食された部分の砂は砂丘の先端部に供給される。 両翼の長さは、 卓越風の強さとともに単調に長くなる傾向を持っている。 (放物線砂丘が見える制御パラメータの範囲の中で) 比較的砂の供給量が少ない場合は、 ヘアピン型のいわゆる延伸放物型砂丘 (elongate parabolic). 先端部を除いて、) $[4]$ が形成される (fig
これは線形砂丘と非常に似た形状をしている。
$.2(\mathrm{d})$$\mathrm{d}\mathrm{u}\mathrm{n}\mathrm{e}\mathrm{s}$以上の結果は図 4 にまとめられている。 これは 3 つの軸をもった Hack の相図とは直接比較できないが、 ここで選択された二つの制御パラメータと出現する砂丘の形状の関係に関して、 Hack の観測結果に定性的にはよい対応を示している。 とくに、砂源にある砂の量や風の強さに応じて、 1. 砂のみから成る横列砂丘 (transverseや、 3 放. 表面の–部に植生のともなう放物線型砂丘 parabolic物型砂丘がさらに発達した延伸放物線型砂丘 (elongate parabolic dune) が形成され、 4. 砂の量がそれらに達しない場合や卓越風の風力が弱すぎる場合、 植生がほぼ全面を覆うという、 4 つの観測事実が我々の簡単な模型計算でも得られることは重要な意味を持っているように思われる。 今回の模型では、 砂丘風下側でのすべり面での砂の運動をクリープの –部と見なしたり植物の成長を線形写像で表すなど系の動力学を極端に簡単化している。 しかしながら、 (植生を含まなくても十分複雑な砂丘の運動に) 植生の影響まで含める試みはほとんどなく、 今回のよ$($$\mathrm{d}\mathrm{u}\mathrm{n}\mathrm{e})_{\text{、}}2$$\mathrm{d}\mathrm{u}\mathrm{n}\mathrm{e})_{\text{、}}$

うな手法による解析は砂丘の緑化過程とその逆の砂漠化過程の理解に–定の貢献ができると思われる。 今後、 この模型により定量的な側面を加えるために植物の成長率や被覆による卓越風への影響など必要に応じて系の詳細を模型に付加して行く予定である。 また、 ルールの簡単さゆえ、 今回用いた動力学を偏微分方程式系に置き換え、 放物線型砂丘の成因や植物との共存 (もしくは排他的な) 関係について数理的な解析をすることも今後の課題となる。

参考文献[1] Cook,R.,Warren,A.,Goudie,A.
press: London

1993.

DesertGeomorphology.UCL[2] Hack,J.T. 1941. Dunes of the western Navajo Country: Geographical Review

31:240-263.

48

図 4:. 色々な制御変数空間での位置と、 各位置で形成される主な砂丘の形の関係を表した相図。 図中のアルファベットの各文字の位置では典型的な形の砂丘が形成される。各文字は図 3 での各回に付したアルファベットの名前と対応している。
はそれぞれ、 放物線型砂丘、 横列型 (またはバルハン) 砂丘をあらわ
記号
し、は砂丘が形成されない状況を意味している。 の中でも*をともなうもの
は、 放物線型砂丘の中でも特に両翼が延びた延伸放物線型砂丘を表している。 口の位置では典型的な形状の砂丘して分類することが困難であった。
$\triangle_{\text{、}}\bullet_{\text{、}}$$\triangle$$\cross$[3] Thomas,D.S.G.,Tsoar,H. 1990. The geomorphological role of vegetation indesert dune systems: Vegetation and.
$\mathrm{E}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{i}_{\mathrm{o}\mathrm{n}}:471-489$[4] Pye,K. 1982. Morphological development of coastal dunes in a humid tropical environment, cape bedford and cape flattery, North Queensland. Geografiska $Annaler64\mathrm{A}:21$ 3-227.
[5] Tsoar,H.,Moller,J-T.1982.
Unwin: Boston: 75-97.

AeorianGeomorphology.

Allenand[6] Halsey,L.A., Catto,N.R., and Rutter,N.W. 1990. Sedimentology and development of parabolic dunes, Grande Prairie dune field, Alberta. CanadianJoumal of Earth Sciences27: 1762-1772.
[7] Werner,B.T. 1995. Eolian dunes: Computer simulations and attractor interpretation. Geology 23: 1107-1110.

[8] Nishimori, H., Yamasaki, M., Anderson, K. H. 1998. A simple model forthe various pattern dynamics of dunes. Intemational Journal of Modern
Physicsl 2:257-272.

49[9] Wasson,R.J.,Hyde,R.1983. Factors determining desert dune type Nature304: 337-339.
[10] de Castro,F.1995. Computer simulation of the dynamics of a dune system: Ecological Modelling78: 205-217.

[11] Rassmussen,K.1989. Proc of Symposium coastal sand dunes. Royal Societyof Edinburg. Edinburgbf B96:129-147.