確率セルオートマトンで表現される 離散的移流現象の漸近挙動の解析
Analysis on asymptotic behavior of discrete advection phenomena expressed by
stochastic cellular automata
2020 年 3 月
延東 和茂
Kazushige ENDO
確率セルオートマトンで表現される 離散的移流現象の漸近挙動の解析
Analysis on asymptotic behavior of discrete advection phenomena expressed by
stochastic cellular automata
2020 年 3 月
早稲田大学大学院 基幹理工学研究科 数学応用数理専攻 非線形システム研究
延東 和茂
Kazushige ENDO
目 次
第1章 序論 3
1.1 背景 . . . . 3
1.2 本研究の目的 . . . . 5
1.3 先行研究. . . . 7
1.3.1 SBCAとその基本図 . . . . 7
1.3.2 4近傍系 . . . . 10
第2章 空間周期無限大における分解仮設を用いた解析 14 2.1 SECA212の基本図の理論式 . . . . 14
2.1.1 ECA212の方程式の表現 . . . . 14
2.1.2 SECA212の基本図 . . . . 15
2.2 SECA84の基本図の理論式. . . . 16
2.2.1 ECA84の方程式の表現 . . . . 16
2.2.2 SECA84の流束 . . . . 18
2.2.3 SECA84の基本図 . . . . 19
2.3 3値化されたSBCAとその基本図 . . . . 21
2.3.1 SBCAの多値化 . . . . 22
2.3.2 3値化されたSBCAの基本図 . . . . 23
2.4 考察と今後の展望 . . . . 31
第3章 空間周期有限における遷移確率行列を用いた解析 33 3.1 空間周期有限におけるSBCA . . . . 33
3.1.1 空間周期有限におけるSBCAの漸近挙動 . . . . 33
3.1.2 空間周期有限におけるSBCAの基本図 . . . . 36
3.2 SBCAの拡張と漸近挙動の理論的導出 . . . . 37
3.2.1 4近傍拡張系 . . . . 37
3.2.2 連立拡張系 . . . . 40
3.3 考察と今後の展望 . . . . 45
第4章 GKZ超幾何関数を用いた熱力学的極限の計算 47 4.1 GKZ超幾何関数の定義と性質 . . . . 47
4.2 SBCAの基本図の熱力学的極限 . . . . 49
4.2.1 SBCAの基本図のGKZ超幾何関数による表現. . . . 49
4.2.2 F1,L−1,mが満たす微分方程式 . . . . 51 4.2.3 F1,L−1,mとF0,L−2,m−1の隣接関係式 . . . . 52 4.2.4 隣接関係式の極限. . . . 53
第5章 結論と今後の展望 55
謝辞 57
参考文献 57
第 1 章 序論
1.1 背景
移流現象とは,物質,またはその物理量が流れに従って移動する現象であり,物理学,化学,生 物学,あるいは社会科学において数多く現れる現象である.例えば,水や空気などの連続体が変化 していく様子は物理における移流現象として捉えられるが,それらの現象を解析する際には,偏微 分方程式が数理モデルとして通常用いられる.偏微分方程式を用いることで,設定した初期条件,
境界条件の下で物質の物理量が微小な時間,空間変化に応じてどのように変化するかを記述するこ とができる.この偏微分方程式を解くことにより,様々な連続体の運動の予測が可能となる.しか し,現象のモデリングによって導出された偏微分方程式は非線形のものが多く,それらは一般的に は解くことが困難とされている.そこで,解析解を求める代わりに,計算機を用いた数値解法によ り,微分方程式の近似解を求めるといった手法がとられることがある.このような手法の代表的な 例として,差分法や,有限要素法と呼ばれるものがある.差分法は,極限によって定義される微分 を差分へ置き換えることにより,元の微分方程式の近似式として差分方程式を得る手法である.ま た,有限要素法は,対象を微小領域に離散化し,それぞれの微小領域の運動を形状関数と呼ばれる 簡潔な数式で表すことで,連続体の運動を数値的に解析する手法である.
このように,連続体の離散化を用いて近似計算を行う手法については盛んに研究されており,差 分法,有限要素法の他にも粒子法,格子気体法,格子ボルツマン法などが存在する[1, 2, 3].粒子 法は,連続体を多くの粒子の集合体であると解釈し,各粒子の動きを質点の運動方程式で記述する ことで,偏微分方程式を連立常微分方程式に帰着させる手法である.格子気体法,格子ボルツマン 法は,離散的な格子の上を複数の粒子が移動するようなモデルをとるもので,粒子同士の衝突や並 進など,多数の粒子の相互作用を計算することで,連続体の運動を近似的に求める手法である.
上述した手法は,連続体の近似計算としての離散的手法であるが,一方,はじめから離散的かつ 微小な粒子を想定し,それらの集合体としての物質の状態や運動を取り扱う学問として,統計力学 がある.統計力学は,微視的な物質の運動や状態が確率分布に従う変数であるとし,確率変数の期 待値や種々の統計量によって巨視的な物質の性質,エネルギー,エントロピーなどを理解すること を目的とした分野である.特に,非平衡な現象を取り扱うような非平衡統計力学の分野において,
ASEP(asymmetric simple exclusion process)と呼ばれる,移流現象を解析する離散的な確率モデ ルが存在する.ASEPは,多数の粒子が一次元の格子上を体積排除の原理に従いつつ確率的に移 流するような,多体ランダムウォークモデルである.このモデルはSpitzerによって提案され[4],
Derridaによって行列積の方法と呼ばれる解析手法が編み出された[5, 6].その後,笹本によって 直交多項式との関係が解明されるなど[7],現在でも活発に研究が進められているモデルである.
連続と離散の関係に着目した移流現象の研究においては,可積分系における試みも挙げることが できる.可積分系とは,具体的な解を求めることができるような微分方程式,差分方程式である.
例えば,可積分系の1つであるKorteweg-de Vries (KdV)方程式は,衝突しても壊れないような 複数の孤立波が,衝突後に位相のずれが起こるものの元の状態に戻るような,ソリトン解を特解と して持つ偏微分方程式である[8].一方,一次元の格子上を複数の粒子が移流していくモデルの一 つとして,箱玉系がある[9].箱玉系の時間発展において,隣接する粒子の集合を波と解釈すると,
それらの波がソリトン性を持ちながら移流する様子が観察できる.さらに,超離散化という特殊な 極限操作を通して,KdV方程式と箱玉系を代数的に関連付けることができる.このように,連続,
離散に限らず移流現象を表すモデルに関する研究が可積分系において存在している[10].
離散的な移流現象を数学的に取り扱う手法の一つとして,セルオートマトンというモデルが存在 する.セルオートマトンは,時間,空間,状態変数が全て離散的であるような時間発展系である.
1960年代に,Neumannによって,自己を複製する機構を持つ機械の運動モデルとして提唱された [11].その後,生命現象の過程を非常に単純な相互作用によって再現,計算しようと,ライフゲー ムというモデルがConwayによって開発された[12].また,徐々にセルオートマトンが力学系の一 種であると解釈され始めると,Wolframによって,一次元のセルオートマトンの漸近挙動を定性 的に4つのクラスに分類するなどの研究がなされた[13, 14].前述した各分野の研究にもセルオー トマトンは関係しており,数値計算における格子気体法,格子ボルツマン法,統計力学における ASEP,可積分系における箱玉系はセルオートマトンの一種とすることができる.また,いくつか のセルオートマトンはmax-plus方程式と呼ばれる形式によってその時間発展が表現される場合が ある.このmax-plus方程式は,前述したKdV方程式と箱玉系を超離散化によって対応付ける際に 用いられる形式であるが,超離散化に限らず,セルオートマトンをmax-plus方程式で表現し,そ れらに代数的操作を施すことで,解の振る舞いを解析するような研究も行われている[15, 16, 17].
一次元の2値3近傍セルオートマトンは,セルオートマトンにおける最も基本的な系のクラスで あり,ECA(elementary cellular automata)と呼ばれている.ECAには256の独立した時間発展 則が存在しており,そのうちの184番目の発展則を持つ系は,非自明な離散的移流現象を表す系と して注目されてきた.この系は,前述した超離散化によって,Burgers方程式と呼ばれる非線形偏 微分方程式と代数的に対応付けることができているため[18],BCA(Burgers cellular automaton) とも呼ばれることがある.この系も箱玉系同様,一次元の格子上を複数の粒子が移流する様子が観 察できるが,格子数に対する粒子数の割合(粒子密度)にある閾値が存在する.すなわち,粒子密 度が閾値を超えると,解の時間無限大における漸近挙動が劇的に変化するような,相転移が起こる ことが知られている.そのような相転移の様子を明らかにするための手法として,基本図の導出が ある.基本図は,粒子密度の変化によって,時間無限大での粒子の平均流束,具体的には,1時刻 あたりに1つの粒子がどれだけ移流したかという量の空間平均をとったものが,どれほど変化した かをグラフとしたものである.BCAは,排他律によって粒子が1方向に移流する系であり,これ を原初的ではあるが交通流を表すモデルであると解釈することができる.よって,基本図という平 均流束の粒子密度依存性を調べることによって,交通流の漸近的な渋滞生成のメカニズムが,簡易 的にだが理解できることになる.
BCA以外にも,セルオートマトンを用いて交通流を理解しようという研究は盛んにおこなわれてき た.代表的な例として,1990年代にNagel, Schreckenbergによって提案されたNagle-Schreckenberg (NS) モデルがある[19, 20].これは,設定された最高速度(1時刻で移動しうるサイト数)の下,
最も近くのサイトにいる粒子との距離に応じて粒子が加速と確率的な減速を行いながら移動し時間 発展していくような,確率的なセルオートマトンである.加速や減速の機構が組み込まれているこ とで,前述したBCAより現実に近い交通流を再現できるモデルとなっており,実際の交通流を計 測した結果とモデルの解析結果を比較するような研究もなされている.モデルにおいて最高速度が 2以上の場合,平均場近似を用いた解析が行われているが,最高速度が1の場合,いくつかのサイ トにわたる各粒子の存在の条件付確率についての仮設を用いて,理論的に基本図の導出が行われて いる.
確率セルオートマトン以外にも,離散的な確率的現象を取り扱うための手法としてマスター方程 式を用いた解析が存在する.マスター方程式とは,微小時間における離散的な状態の遷移速度を行 列方程式によって表現したマルコフ過程で,例えば,生物の集団の個体数の確率的な時間変化を表す モデルの一つとして,出生死滅過程がある[21].ある時刻tにおいて生物の個体数がi(= 0,1,2,· · ·) である確率をp(t) = (p⃗ i=0(t), pi=1(t), pi=2(t),· · ·)とすると,
d ⃗p(t) dt =p(t)⃗
−λ λ 0 0 . . .
µ −(λ+µ) λ 0 . . .
0 µ −(λ+µ) λ . . .
0 0 µ −(λ+µ) . ..
... ... ... . .. . ..
(λ, µ:正の実数)
で定義される確率過程は,出生死滅過程の一つの例である.この確率過程を基に,遷移確率行列Q が,
Q=
0 1 0 0 . . .
µ
λ+µ 0 λ
λ+µ 0 . . .
0 µ
λ+µ 0 λ
λ+µ . . .
0 0 µ
λ+µ 0 . .. ... ... ... . .. . ..
として導出される.この遷移確率行列に対してπQ=πとなるベクトルπを計算することによっ て,出生死滅過程の定常分布を求めることができ,この場合,ρ= λ
µ <1の時に定常分布が存在 し,幾何分布πi= (1−ρ)ρi(i= 0,1,2,· · ·)となることが知られている.
1.2 本研究の目的
本論文では,以上のような背景を基に,確率セルオートマトンによって表現される離散的かつ確 率的な移流現象について,その漸近挙動を理解することを目的とし,いくつかの手法を用いて研究
を行う.なお、本論文で扱うセルオートマトンすべてに周期境界条件を課す.2章では,高次保存 量を持つ系と状態変数が多値に拡張された系について解析を行い,3章では,2章とは異なるアプ ローチを提案したうえで,まず基本的な系に対してその手法を適用し,その後,多近傍拡張系,連 立系への拡張を提案し,同様の手法を用いて解析を行う.4章では,ある種の級数によって表現さ れた系の漸近挙動の熱力学的極限を計算する.
前節において紹介した研究は,移流現象を連続体,あるいは離散的な要素の集合体とみなして解 析を行うものがほとんどであるが,その場合,移流物体の生成消滅がなく系全体の質量が保存され ている前提が最も基本的な系の構成である.確率セルオートマトンで表現される移流現象の漸近 挙動の解析においても,粒子数が不変であるような粒子系に対して行われるものがほとんどであっ た.2章ではこの点に着目し,まず章の前半において,粒子数ではなく,より高次の保存量を持つ ような確率セルオートマトンに対して,その保存量密度を与えるパターンが移流していく様子を解 析することを目的とする.その際,Fuk´sによって提案された,高次保存量を持つ系と粒子系を対 応付ける変数変換を用いる[22].この変数変換によって,高次保存量を持つ確率セルオートマトン を粒子数が保存される確率セルオートマトンに帰着させ,後に先行研究で示す分解仮設を用いた手 法によって,基本図を理論的に導出する.分解仮設は,前述したNSモデルの解析に用いられた,
いくつかのサイトにわたる各粒子の存在の条件付確率についての仮設と類似するものであり,局所 的な特定のパターンの存在確率をより短い長さのパターンの存在確率の組み合わせに分解するため の仮設である[23, 24].この分解仮設の正当性を数値実験によって確認したうえで,分解仮設を平 均流束に関する平衡方程式に複数代入することで,平衡方程式を平均流束のみで閉じた形式で表現 する.その方程式を解くことで,平均流束を粒子密度と確率パラメータのみで表し,基本図の理論 的導出とする.章の後半では,状態変数のとり得る値が3値に拡張された系に対して,同様に分解 仮設を用いた基本図の理論的導出を目指す[25].それまで解析が行われてきた系と異なり,数値的 に成り立つと確認できる分解仮設の数が少なく,それを補填するための平衡方程式を複数連立させ つつ,利用できる分解仮設を用いてその連立方程式をなるべく少ない変数で表す.
2章では空間周期無限大の設定の下で系の解析を行ったのに対して,3章では空間周期有限の条 件での時間無限大における漸近挙動の解明を目的として,2章とは異なるアプローチを提案する.
これは,2章における空間周期無限大の下での局所的なパターンの存在確率の定義の曖昧さを回避 するためである.これらの存在確率は,本来空間周期無限大の設定において,時間無限大での収束 値として与えられる.しかしながら,この条件下での分解仮設の厳密な導出は困難である.した がって,存在確率の正当性を担保するためには,有限の空間サイズおよび有限時間での数値計算に 頼らざるを得ない.これでは空間,時間共の二重の意味での近似になってしまう.このことを踏ま えた上で3章ではまず,与えられた有限の長さの全空間サイトにわたる実現値の列(configuration と呼ぶ)を一つの状態とみなし,それらの集合上を確率的に遷移するマルコフ連鎖として,確率セ ルオートマトンの解析を行う.すなわち,粒子系に対して空間サイト数と粒子数を設定したうえ で,確率過程としての極限分布の一意性を仮定できるようにconfigurationの状態空間を設定し,各
configuration間の遷移確率を考慮し,前述の出生死滅過程のように遷移確率行列を導出すること
で,時間発展を離散時間のマスター方程式で表す.その方程式の固有値1に対応する固有ベクトル を求め規格化することで,系の各configurationの極限分布を導出する.また,このような解析を 異なる空間サイト数,粒子数に対しても繰り返し適用しそれらの極限分布をみることで,一般の有
限な空間サイト数,粒子数に対する極限分布の予想を立てる.また,予想により導かれた極限分布 を基に,基本図を確率変数の期待値として導出する.さらに,このような手法が適用できるような 拡張系を探索し,極限分布の予想と基本図の導出を行う.その際,得られた拡張系の極限分布は2 章において解析を行った高次保存量を持つ系のそれと等価であることを示す.
4章では,3章で導出した空間周期有限における基本図の熱力学的極限を計算する.その際,基 本図がGKZ超幾何関数と呼ばれるGelfand, Kapranov, Zelevinskyが提唱した超幾何関数で表現 されることを利用する[26].GKZ超幾何関数は,ガウスの超幾何関数の多変数拡張とも呼べるも ので,ガウスの超幾何関数と同様,GKZ超幾何関数が満たす隣接関係式や微分方程式が導出でき る.それらの性質を利用することで,隣接関係式の極限として,基本図の熱力学的極限の計算を
行う[27, 28].また,ここで導出された結果が,先行研究における結果と一致することも確認する.
すなわち,2章までの空間周期無限大における分解仮設を用いた解析結果と,3章の空間周期有限 における遷移確率行列を用いた解析結果を結びつける,熱力学的極限の手法としてのGKZ超幾何 関数を用いた計算手順を紹介する.後にも述べるが,分解仮設,遷移確率行列を用いた手法におい ては,共に時間無限大における状態を予想するための仮設が基幹的な役割を果たしており,それら の仮設は数学的証明が待たれるものである.しかし,それら独立な手法で導かれた結果がGKZ超 幾何関数を用いて対応付けられることで,それら手法を用いた解析の確かさを議論する上で有用 な傍証となりうる.なお,このGKZ超幾何関数を用いた手法は,misanthropeモデルと呼ばる離 散確率粒子系の,漸近挙動の熱力学的極限を導出する際に用いられる手法と同様のものであり,金 井,筧,金丸によって考案されたものである[28, 29].
以上に述べたように,本論文ではいくつかの確率セルオートマトンに対して,それらが表現する 移流現象の漸近挙動を解析することを目的とし,有限もしくは無限の空間周期の設定の下で異なる アプローチと仮設を用いつつ,基本図,あるいはその背後にある極限分布の導出に関する報告を 行う.
1.3 先行研究
まず先行研究として, 分解仮設を用いた確率セルオートマトンの解析手法についての説明を 行う.具体的には,前述したBCA(Burgers cellular automaton)に対して確率変数を導入した,
SBCA(stochastic Burgers cellular automaton)を扱う.
1.3.1 SBCA とその基本図
SBCAは
un+1j =unj +qjn−1−qnj, qnj = min(anj, unj,1−unj+1), anj =
1 (確率α) 0 (1−α)
(1.1)
によって与えられる確率セルオートマトンである.ただし,anj の値はj, n毎に独立に与える.こ の時,unj はj番目のサイトの時刻nでの状態変数を表しており,確率変数と{unj−1, unj, unj+1}の
値によって次の時刻の状態変数,すなわちun+1j が定まる.このような次の時刻への状態の遷移の 仕方を表にしたものを遷移表といい,SBCAの遷移表は下表のようになる.
unj−1unjunj+1 111 110 101 100 011 010 001 000
un+1j 1 1−anj anj−1 anj−1 1 1−anj 0 0 (1.2) さらに,図1.1に時間発展の例を示す.時間発展に際しては周期境界条件をとっており,1周期分
- j
? n
図1.1: SBCAの時間発展の例.α= 0.5. ■,□はそれぞれu= 1, 0のサイトを表し,横軸が右向
きにj,縦軸が下向きにnとなっている.
のサイト数をLとしている.式(1.1)の両辺を空間サイトjに関して足し合わせると,右辺のqjn−1 とqjnが消去され,X
j
un+1j =X
j
unj となり,状態変数の総和が時間不変である.u= 1のサイト を粒子が存在,u= 0のサイトを粒子が存在しないとみなした時に,系が粒子の総数を変化させず に粒子が移動するものとみなせる.このことを本論文では系の状態変数が粒子性を有すると言うこ ととする.粒子性の下で,右辺のqnj−1とqjnはそれらの符号から,1ステップ分の時間発展によっ てj番目のサイトに流入する粒子の数と流出する粒子の数であると無矛盾に解釈することができる.
さらに,qjn−1= min(anj−1, unj−1,1−unj)より,粒子が流入する条件はanj−1= 1, unj−1= 1, unj = 0 であり,流入量は1となる.流出する条件,粒子数も流入と同様であるため,サイトが1をとる状 態を粒子とみなすと,粒子のダイナミクスは,粒子が存在するサイトの右側が0である時に,確率 αで右側のサイトに移動すると解釈することができる.また,この系は前述したASEPのある特 殊な場合と等価である.
10 確率α
このような確率的な移流現象を表す系に対して,時間無限大での粒子の動きを解析することを目指 す.具体的には,時間無限大における流束の空間平均を
Q= lim
n→∞ lim
L→∞
1 L
XL j=1
qnj (1.3)
と定義する.これは,時間無限大において,1ステップ分の時間発展によって粒子が移流する総空 間サイト数の空間平均である.この平均流束Qの,状態変数の密度
ρ1= lim
L→∞
1 L
XL j=1
unj (1.4)
との関係を示す基本図を導出することが目的である.Qは図1.2に示すように,初期値によらず一 意的にuの密度ρ1および確率パラメータαに依存することを示唆している.ただしnは無限大の 極限の代わりに有限の値で観測を行い,周期長Lも同様に有限の数をとっている.Lを無限大に した設定におけるSBCAの基本図の理論式は
図1.2: SBCAの基本図(実線が理論曲線,プロットが数値計算,上から順にα= 1,0.8,0.6,0.4,0.2)
Q=1−p
1−4αρ1(1−ρ1)
2 (1.5)
で与えられるが,その導出の概要は以下の通りである[24].
まず,時刻nの解において隣り合うm個の空間サイトにわたる状態値x1x2...xm(∈ {0,1}m)の 密度を
ρnx
1x2...xm = lim
L→∞
1 L
XL j=1
#x1x2...xm (1.6)
で定義し,パターン密度と呼ぶ.ここで,#x1x2...xmはm近傍のパターンx1x2...xmの数とす る.SBCAの場合は,ρn1 については前述した粒子性より時間によらずρn1 =ρ1である.さらに,
n→ ∞における極限が一意に存在すると仮定し,これをρx1x2...xmと記す.これらのような時間 無限大におけるローカルなパターンの密度に関して,以下の分解仮設[24]
ρx1x2x3 =ρx1x2ρx2x3 ρx2
, ρx1x2x3x4 =ρx1x2x3ρx2x3x4 ρx2x3
(1.7) が成り立つと仮定する.これらの形式では,左辺における近傍数の大きいパターン密度が,右辺に おける近傍数の小さい密度の組み合わせに分解されていることがわかる.この性質が確率セルオー
トマトンにおける基本図の理論的導出に際して重要となる.しかし,この分解仮設は数学的に存在 が証明されておらず,数値実験によって近似が成り立つことが確かめられているのみである.例え ば,ρ110, ρ100, ρ101, ρ010, ρ1010に対する分解仮設は,
ρ110= ρ11ρ10
ρ1 , ρ100= ρ10ρ00
ρ0 , ρ101=ρ10ρ01
ρ0 , ρ010=ρ01ρ10
ρ1 , ρ1010= ρ101ρ010
ρ01 (1.8) となる.数値的検証は,適当な周期長Lにおいて粒子密度を変えながら十分時間発展させ,左辺 と右辺のパターンの個数を測定する方法を採用している.両辺の実測値を重ねてプロットしたもの を図1.3に示す.つぎに,ρ10について,次の時刻でパターン10を与える遷移確率を考慮すること により,n→ ∞において他のパターン密度と以下の平衡方程式が成り立つ.
ρ10= (1−α)ρ10+α(ρ110+ρ100) +α2ρ1010 (1.9)
(1.9)に(1.8)を代入すると,右辺はすべて2近傍以下のパターン密度で表現される.また,2値セ
ルオートマトンの恒等式として
ρ0= 1−ρ1, ρ00=ρ0−ρ10, ρ01=ρ0−ρ00, ρ11=ρ1−ρ10 (1.10) があり,これらを用いて右辺を整理すると,以下のρ10に関する2次方程式が得られる.
α(ρ10)2−ρ10+ρ1(1−ρ1) = 0 (1.11) この方程式を解いて0≤ρ10≤1を考慮し,2次方程式の符号を選択すると,
ρ10= 1−p
1−4αρ1(1−ρ1)
2α (1.12)
となる.Qはqnj の平均で与えられるので,(1.1)のqjnが1となる条件,すなわち,(anj, unj−1, unj) = (1,1,0)を考慮すると,
Q=αρ10 (1.13)
で与えられることがわかる.よって,Qの理論式は(1.5)になることが示される.
1.3.2 4 近傍系
独立な確率変数bnj を
bnj =
1 (確率β) 0 (1−β)
(1.14)
で定義し,流束qnj が以下で定義される3種類の4近傍の確率2値セルオートマトン
qjn= min(anj−1+unj, anj +bnj, uj−1+unj,1−unj+1), (1.15) qjn= min(max(−min(bnj, unj+1),min(anj, unj−1+unj −1)),1−unj+1), (1.16) qjn= min(max(0,min(anj, unj−1+unj −1)),1−unj+1) (1.17)
について,基本図の数値計算結果を図1.4に示す.これら基本図の理論式も[24]に報告されてお り,たとえば(1.16)については
Q=
−1−ρ1−p
(1−ρ1)2−4βρ1(1−2ρ1)
2 (ρ1≤1/2)
ρ1−p
ρ21−4α(1−ρ1)(2ρ1−1)
2 (ρ1>1/2)
(1.18)
で与えられる.
0.2 0.4 0.6 0.8 1.0 0.02
0.04 0.06 0.08 0.10
ρ010
0.2 0.4 0.6 0.8 1.0
0.02 0.04 0.06 0.08 0.10 0.12
ρ100
0.2 0.4 0.6 0.8 1.0
0.02 0.04 0.06 0.08 0.10
ρ101
0.2 0.4 0.6 0.8 1.0
0.02 0.04 0.06 0.08 0.10 0.12
ρ110
0.2 0.4 0.6 0.8 1.0
0.01 0.02 0.03 0.04
ρ1010
図1.3: 各分解仮設の数値的検証.横軸がρ1,縦軸がパターンの出現割合.黒点がρx1x2x3ρx2また はρx1x2x3x4ρx2x3,白丸がρx1x2ρx2x3またはρx1x2x3ρx2x3x4の数値.α= 0.7,L= 4000において n= 100000まで時間発展させ,n= 88000以降の数値の時間平均をとり,ρ1ごとにプロットした もの
(a)
(b)
(c)
図1.4: (a), (b), (c)はそれぞれ(1.15), (1.16), (1.17)の基本図(数値計算)
第 2 章 空間周期無限大における分解仮設を用 いた解析
前章では,SBCAをはじめいくつかの粒子性を持つ確率セルオートマトンについて,その漸近 挙動における基本図の理論式を導出する手法について説明した.その際,平衡方程式として導出さ れた粒子の流束の発生条件を短い近傍のパターン密度で記述するために,分解仮設は非常に重要な 役割を果たしている.また,先行研究として紹介したこれらの解析は,すべて2値の粒子性を持つ 確率セルオートマトンに対してのものであった.
しかしながら,平均流束の密度依存性のような重要な力学的指標の解析は系が2値であったり,
粒子性を持つことに限定する必要はない[22, 30, 31, 32, 33].例えば2値セルオートマトンにおい てパターン01と10の密度の和σ01+σ10が保存量となるものがECAに存在し,それらを変数変 換を通して粒子系に帰着させることができる.また,max-plus表現されたセルオートマトンに着 目し,状態変数を多値化した系の性質について解析された例もある[16].
この章では,そのような2値粒子系以外の確率的な系について,その漸近挙動の理論的な解析を 行う.その際,先ほど紹介した分解仮設を用いた手法を用いる.具体的には,まず2次保存量を持 ついくつかのセルオートマトンをmax-plus方程式で表現し,それに保存量を保ったまま確率変数 を導入する.その後,超離散コール-ホップ変換と呼ばれる変数変換を用いて[22],それらの系を 粒子性を持つ確率セルオートマトンに帰着させる.そのセルオートマトンに対して分解仮設を用い て基本図の導出を行う.また,多値化された系についても,同様に分解仮設を用いて基本図を導出 する.
2.1 SECA212 の基本図の理論式
2.1.1 ECA212 の方程式の表現
まず,以下で与えられる時間発展方程式を考える.
vn+1j =vnj ⊕min(vnj−1⊕vnj,1−vjn⊕vj+1n ) (2.1) ここで⊕は
x⊕y=x+y mod 2 (x, y∈ {0,1}) (2.2) で定義される2進演算であり,x∈ {0,1}に対して以下の公式が成り立つ.
x⊕0 = 0⊕x=x, x⊕1 = 1⊕x= 1−x, x⊕x= 0
(2.1)は以下のルール表で定義することもできる.ルール表を見ると,時刻n+ 1における状態変数を 数列とすると11010100となる.これを10進数に直すと212となることから,このルール表で定義さ れる系をECA212と表記する.なお,ECA212は文献[32]で取り扱われているECA142とreflection において等価である.この時,reflectionにおいて等価とは,vjn+1=f(vnj−1, vjn, vnj+1)によってECA の時間発展のルールを定義するとき,2つのルールf1とf2がf1(vjn−1, vnj, vj+1n ) =f2(vnj+1, vjn, vnj−1) となることである.
vjn−1vnjvnj+1 111 110 101 100 011 010 001 000
vjn+1 1 1 0 1 0 1 0 0 (2.3)
変数変換
unj =vjn−1⊕vnj (2.4)
を考えると,(2.1)より次のuに関する時間発展方程式が得られる.
un+1j =unj ⊕min(unj−1,1−unj)⊕min(unj,1−unj+1) (2.5) さらに,この方程式の右辺は2つの⊕をそれぞれ演算+,−に置き換えても同じ結果を与える.
un+1j =unj + min(unj−1,1−unj)−min(unj,1−unj+1) (2.6) この方程式はECA184に他ならず,保存量ρ1を有している.したがって(2.4)を通じてECA212 (2.1)は2次保存量
σ01+σ10= 1 L
XL j=1
vnj ⊕vnj+1 (2.7)
を有している[22, 32, 33].混乱を防ぐため,uとvの系の密度をそれぞれρ,σで記す.なお,周 期境界あるいは境界条件 lim
|j|→∞unj = 0を与えた無限領域においては,σ01=σ10であるので,σ01, σ10がともに保存量となる.我々は以降でσ01=σ10を前提とする.
2.1.2 SECA212 の基本図
ECA184は確率変数を導入して1次保存量を壊さずに確率化が行え,SBCA (1.1)が得られた.
この仕組みを利用すれば変換(2.4)を通じてECA212の確率化をSBCAと関連付けられる.その 図式は以下の通りである.
un+1j =unj + min(anj−1, unj−1,1−unj)−min(anj, unj,1−unj+1) (SBCA)
↕ (2値なので同値)
un+1j =unj ⊕min(anj−1, unj−1,1−unj)⊕min(anj, unj,1−unj+1) (SBCA)
↕ unj =vnj−1⊕vnj (変数変換)
vjn+1=vjn⊕min(anj, vnj−1⊕vnj,1−vjn⊕vj+1n ) (SECA212)
(2.8)
この確率化によってもSECA212はσ01+σ10が保存量となることは明らかである.なお,SECA212 のルール表は以下のようになる.
vjn−1vjnvj+1n 111 110 101 100 011 010 001 000
vn+1j 1 1 0 anj 1−anj 1 0 0 (2.9)
またSBCAの平均流束Qは,変換によってSECA212では
Q= lim
n→∞ lim
L→∞
1 L
XL j=1
min(anj, uj−1⊕vjn,1−vnj ⊕vnj+1) (2.10)
となる.unj = 1は(vnj−1, vjn) = (0,1)もしくは(1,0)を与えるので,このQはSECA212において 連続する1個以上の0と1の境目の平均流束を表している.そしてSECA212における基本図は,
図1.2で示されたSBCAの基本図の横軸を2σ01 (=σ01+σ10)に,縦軸を(2.10)に置き換えたも ので与えられ,その理論式は
Q= 1−p
1−8ασ01(1−2σ01)
2 (2.11)
となる.
2.2 SECA84 の基本図の理論式
2次保存量を持つECAはECA212以外にも存在する.この章では,そのひとつであるECA84 について述べ,確率化を行い,基本図の理論式を導出する.
2.2.1 ECA84 の方程式の表現
ECA84のルール表は
vjn−1vnjvnj+1 111 110 101 100 011 010 001 000
vjn+1 0 1 0 1 0 1 0 0 (2.12)
で与えられる.この表を4近傍のルール表に書き改めると
vnj−1vnjvj+1n vj+2n 1111 1110 1101 1100 1011 1010 1001 1000 vn+1j vn+1j+1 00 01 10 11 00 01 10 10
0111 0110 0101 0100 0011 0010 0001 0000
00 01 10 11 00 01 00 00
(2.13)
となる.これよりvのある時刻nにおける密度σについて
σn+101 =σn1110+σ1010n +σn0110+σ0010n =σ10n σn+110 =σn1101+σ1001n +σn1000+σ0101n =σ10n
(2.14)
となり,
σ01n+1=σ10n+1=σn10 (2.15) となる.周期境界あるいは境界条件 lim
|j|→∞unj = 0を与えた無限領域においては,常にσ10n =σ01n であるので,両者はそれぞれ保存量σ01,σ10となる.
ルール表(2.12)を方程式で表現すると,たとえば
vjn+1=vjn⊕min(vjn−1⊕vjn+vj+1n ,1−vjn⊕vnj+1) (2.16) となる.この方程式に(2.4)と同じ変換
unj =vjn−1⊕vnj (2.17)
を施すと, lim
|j|→∞vjn= 0の境界条件の場合は,境界条件 lim
|j|→∞unj = 0の下での un+1j =unj ⊕min(unj−1+ ∞⊕
i=j+1uni,1−unj)⊕min(unj + ∞⊕
i=j+2uni,1−unj+1) (2.18) という無限近傍のuに関する時間発展方程式が得られる.ただし,
k2 i=k⊕1
fk =fk1⊕fk1+1⊕ · · · ⊕fk2 (2.19) と定義する.(2.18)はさらに
un+1j =unj + min(unj−1+ ∞⊕
i=j+1
uni,1−unj)−min(unj + ∞⊕
i=j+2
uni,1−unj+1) (2.20) と保存形で書き換えてもu∈ {0,1}ならば右辺の値は変わらない.よってuの密度ρ1が保存量と なる.
しかしながら,上式の2つの⊕uiの項は,周期境界条件のときにうまく定義できない(not well- defined).そこで,(2.20)のu= 1を粒子とみなし,時間発展方程式をEuler表現からLagrange 表現に書き換える.まず,(vnj−1, vjn) = (0,1)を記号Lで,(1,0)をRで表すこととする.すると,
以下に示すようにu= 1の値を持つサイトはLとRの2種類の記号で表される.
vjn: . . .00011110001010011000. . .
unj : . . .000L000R00LRLR0L0R00. . . (2.21) LとRは,vnj の空間パターンにおける1の列011. . .110の左端と右端に対応する記号である.ま た,LとRを粒子とみなしたとき,uのサイトj,j+ 1間の流束は(2.16)より
qnj = min(vjn−1⊕vjn+vj+1n ,1−vnj ⊕vnj+1) (2.22) で与えられる.ルール表で書き改めると
vjn−1vnjvnj+1 111 110 101 100 011 010 001 000
qnj 1 0 0 1 1 0 0 0 (2.23)
となり,上の例では
vjn: . . .00011110001010011000. . . unj : . . .000L000R00LRLR0L0R00. . . qjn: . . .00011101000001010100. . .
(2.24)
となる.一般に,時刻nでの左からk番目のL,Rのサイトをそれぞれxnk, yknと表すと,時間発 展のLagrange表現は
xn+1k =ykn−1, yn+1k = min(ykn+ 1, xnk−1) (2.25) となる.すなわち,LとRは以下のルールで表される粒子モデルとなる.
• LはRが右に隣接していればその場にとどまり,そうでない場合は自分の右の最も近いRの 左隣のサイトに動く.
• Rは右にLが隣接していればその場にとどまり,そうでない場合は右に1サイトだけ動く.
なお,このルールより,n≥1ではLとその右のRの配置はLRかL0Rのどちらかになる.
2.2.2 SECA84 の流束
ECA84の保存量σ01(=σ10)を壊さずに外部変数(anj ∈ {0,1})を以下のように導入する.
vn+1j =vnj ⊕min(min(anj, vnj−1⊕vjn) +vj+1n ,1−vjn⊕vj+1n ) (2.26) 二進ルール表の形で時間発展則を表すと以下のようになる.
vjn−1vnjvnj+1 111 110 101 100 011 010 001 000
vjn+1 0 1 0 anj 0 1 0 0 (2.27)
この表より,前節のL, Rの移動ルールは
• LはRが右に隣接していればその場にとどまり,そうでない場合は自分の右の最 も近いRの左隣のサイトに動く.
• Rは右にLが隣接していればその場にとどまり,そうでない場合はanj = 1なら右 に1サイトだけ動き,anj = 0ならその場にとどまる.
(2.28)
となる. lim
|j|→∞vnj = 0の境界条件の下では(2.17)の変数変換により un+1j =unj+min(min(anj−1, unj−1)+ ∞⊕
i=j+1
uni,1−unj)−min(min(anj, unj)+ ∞⊕
i=j+2
uni,1−unj+1) (2.29) となる.前節と同様にuの時間発展則をLとRの粒子系の移動ルールとみなすと,Lagrange表 現で
xn+1k =ynk −1, ykn+1= min(ynk +anyn
k, xnk −1) (2.30)