Comments
Description
Transcript
非線形力学系の研究
非線形非平衡物理 坂口英継 これまで、大自由度の非線形ダイナミクス、具体的にはカオス結合系、振動子 集団の同期現象、パターン形成などの計算機シミュレーションを行ってきた。 われわれの計算結果を中心に非線形非平衡物理学の研究分野の一端を紹介する。 目次 1. 非線形ダイナミクスとカオス 2. 非平衡系とランジュバンモデル 3. イジングモデルとカオス結合系 4. 相転移とギンツブルグーランダウ方程式 5. 振動子集団の引き込み相転移 6. 散逸構造と位相ダイナミクス 7. 時空カオス 8. スパイラルパターン 9. 局在パターンとソリトン 10. 成長パターン 11. 進化とネットワーク形成モデル １．非線形ダ イナミクスとカオス a. カオスとロジスティック写像 物理学では 、観測した結果を数値化していくつかの法則としてまとめる。観測 量の時間変化は微分方程式（ 偏微分方程式）の形でまとめられることが多い。古 典力学のニュートン方程式、電磁気学のマクスウェル方程式、量子力学のシュレ ディンガー方程式、流体力学のナビエーストークス方程式などである。マクスウェ ル方程式とシュレデ ィンガー方程式は線形方程式の形をし ており、重ね合わせの 原理が成り立つ。線形方程式の解法は昔から数理的に詳し く研究されてきた。非 線形方程式は一般に解くことが難しい。天体力学の３体問題、流体力学の乱流問 題など に非線形性の難しさがよく表われている。非線形方程式は一般に解けない ので 、数値計算でその振る舞いが研究され 、その過程でソリトンやカオスといっ た非線形ダ イナミクスの中心課題が浮かび 上がってきた。 カオスは時間発展の法則が簡単にもかかわらず、その時間発展は不規則で、長 時間先の値の予想が困難な現象をいう。天気予報の長期予想が困難な原因は大気 のダ イナミクスがカオス的だからという説明がされることがある。カオスは３体 問題に関連してポアンカレのころから一部には知られていたようであるが 、流体 のモデル方程式（ローレンツ方程式）の数値計算によりローレンツによりその実 体がわかってきた 1]。 散逸力学系では時間が十分経過するとアトラクタと呼ばれる軌道に吸引され ていく。アトラクタの定性的変化を扱う分岐理論の観点からカオスへの遷移過程 が研究されてきた。典型的には 、定常状態の不安定化により、リミットサイクル と呼ばれる時間周期解が生じ 、リミットサイクルがさらに何回か分岐を起こし て カオスになる。典型的なカオスへの分岐過程に 2n 分岐のカスケード がある。こ れは、リミットサイクルの周期が倍々になり、周期が無限になって非周期すなわ ちカオスになるという過程である。よく知られている 2n 分岐を起こすモデルに ロジスティック写像がある。ロジスティック写像は Xn+1 = AXn (1 ; Xn ) (1) で表される。このモデルはある昆虫の年ごとの個体数の変動を表すモデル方程式 である。A は個体数の増加率を表すマルサス係数に対応する。最初、メイによっ て、2n 分岐やそのカオス化が研究された 2]。昆虫の個体数が周期的に変動する ことは実際に観察されている。インフルエンザや害虫の流行が何年かごとに起こ ることに対応している。その後この 2n 分岐の普遍性やスケーリング則がファイ ゲンバウムなど により研究された。図 1(a) は A を変えた時の 2n 分岐を表してい る。われわれは 、このロジスティック写像を２つ結合した系（結合ロジスティック 写像）の計算を行った 3]。 Xn+1 = AXn (1 ; Xn ) + D(Yn ; Xn ) Yn+1 = AYn (1 ; Yn ) + D(Xn ; Yn ) (2) 結合項は拡散を表す。この系の特徴は各変数 Xn Yn があるパラメータ領域で 2n 周 期を示すのであるが、D = 0 の極限では、その出現する順番がずれていてもよい。す なわち、X が x1 x2 x3 x2 の順であらわれるとき、Y が xm xm+1 xm;1 とずれた順番であらわれてもよいということである。すなわち 、2n 個の順番が n 1 1 0.9 0.8 0.8 Y X 0.7 0.6 0.6 0.5 0.4 0.4 0.3 2.6 2.8 3 3.2 3.4 3.6 0.2 0.2 A 0.4 0.6 X 0.8 1 Figure 1: (a) ロジスティック写像の 2n 分岐 (b) 結合ロジスティック写像の４重 安定解。４つの記号がそれぞれのアトラクタをあらわす。 ずれた解がこの 2 次元写像系のアトラクタとし て存在する。D = 0 なら 、これら の解はまったく等価であるが D = 0 だと、この縮退が解け 2n 個の異なったアト ラクタになる。このような、単純な機構で 2n 個の多重安定な解が生じ 、2n 個の 異なったアトラクタは 、それぞれ異なった分岐過程を経てカオス化していく。図 1(b) は A = 3:5 D = 0:02 のときの４つのアトラクタを記号を変えて描いている。 自由度をあげると非線形系のふるまいがますます複雑になることを示す簡単な例 である。 6 b. カオス拡散 カオスや乱流の特徴にその混合性がある。ものを溶かすときに 、スプーンで コップの中の水をかき混ぜるのは 、乱流の混合性を利用している。流れにのった粒 子のカオス的拡散をラグランジアンカオスと呼ばれる。Gollub 等は振動するベナル対流において、対流セルをまたいで粒子が拡散する様子を実験と計算機シミュ レーションで示した 4]。粒子の慣性を無視すると、粒子は流れと同じ速度で運動す る。その運動方程式は 、流れ関数を (x z t) = A= sin x + B sin(!t) sin(z ) とすると、 f dx = @ dz = ; @ dt @z dt @x g (3) この方程式は 、ハミルトン系の方程式で 、面積要素を保存する保存力学系になっ ている。保存力学系には 、アトラクタはなく、初期条件でカオスになったり , 安定 周期軌道になったりする。われわれは 、粒子の慣性をストークス近似で取り入れ 2 0.5 0.51 40 0.505 20 x(t) z 1 0.5 -20 0.495 0 0 0.5 x 1 0.49 0.49 0 -40 0.495 0.5 0.505 0.51 0 5000 t 10000 Figure 2: (a) A = 1 B = 0 = 20 k = 2 壁に吸引される。(b) A = 1 B = 0:001 = 20 対流セルの隅 (x = 0 z = 0) 付近の拡大図。カオスアトラクタを 示す。(c) x の時間発展。カオス的拡散を表す。 たモデル dvx = ( @ ; v ) dvz = (; @ ; v ) (4) dt @z x dt @z z を計算し た 5]。vx vy は粒子の速度をあらわす。 は質量に逆比例する係数で、 質量 0 の極限で元のモデルに戻る。慣性があるために 、流速と粒子の速度にわず かにずれが生じ 、このモデル方程式は散逸方程式になる。系の振る舞いは初期条 件によらずそのアトラクタで決まる。定常な対流条件では対流セルの壁に漸近す る 、無限周期の軌道（ ヘテロクリニック軌道 ）がアト ラクタになる。粒子は流 体を閉じ 込めている上下の壁の隅の方に溜まってくる。この軌道はわずかな摂動 （ B = 0) でカオス化する。カオス化すると , 対流セルを渡り歩くようになり、粒 子はカオス的に拡散する。無限小の摂動でカオス化する変わった力学系の例でも あり、カオス拡散の簡単なモデルでもある。 6 ２．非平衡系とランジュバンモデル a. ランジュバン方程式と乗法雑音 熱平衡統計力学では 、ある温度 T であるエネルギ− E をもつ状態が実現され る確率は 、ボルツマン分布 P (E ) exp( E=kT ) で与えられる。一般の非平衡 系では定常でも、その分布は一般には分かっていない。非平衡系の一つの大きな 研究分野に、熱平衡へ至る緩和過程や線形応答の研究がある。Green-Kubo 公式 や揺動散逸定理など がその大きな成果である 6]。 現象論的な緩和過程の確率モ デルにランジュバン方程式がある。一番簡単な線形ランジュバン方程式は / ; dv = ;v + (t) + f (t) (5) dt である。(t) は h(t)(t0 )i = 2T (t ; t0 ) を満たす白色ガウス雑音で、f (t) は外力 を表す。f = 0 の時の v の定常分布は 、P (v) / expf;v2 =(2T )g となる。これは 速度 v に関するマクスウェル分布をあらわし 、運動エネルギー E = v2 =2 に対す 3 1 0 -0.1 0.1 P(y) -0.2 0.01 -0.3 -0.4 0.001 -0.5 0.0001 0.001 0.01 -0.6 0.1 1 10 100 0 0.5 1 1.5 2 2.5 3 y Figure 3: (a) = 1 T = 1=2 D = 3=4 に対する定常分布。両対数グラフ (b) = 1 D = 1=8 and F = 0:1 の時間相関関数のフーリエ変換および 線形応答係数 から求めた 00(!) るボルツマン分布でもある。正弦的力 f (t) = F sin !t に対する線形応答係数とゆ らぎ の間に揺動散逸関係と呼ばれる比例関係が成り立つ。外力に敏感な系ほど よ く揺らぐ というもっともらしい関係を定量的に表した関係である。その逆を擬人 的に言うと、頑固な人は考えは揺らがないが他人の意見も聞かないというわけで ある。 この、揺動散逸関係はより一般的な線形ランジュバン方程式でも成り立つ。われ われは 、線形成長率（ マルサス係数）がランダムに揺らぐ 線形ランジュバン方程 式の計算を行った 7]。 dy = (; + (t))y(t) + (t) + f (t) (6) dt (t) は h (t) (t0 )i = 2T (t ; t0 ) を満たす別の白色ガウス雑音である。成長率がゆ らぐ 系は生態系や株価の変動など の社会経済現象ではそれほど 不自然ではない。 f (t) = 0 の時、やはりフォッカープランク方程式の定常解から P0(y) / (T + Dy2 )1( +D)=(2D) (7) 3(a) にこの定常分布の例を両対数グラフで 示す。この定常分布はエネルギ ー E (y) = 1=2( + D)y2 に対する q 次（ q = 1 + 2D=( + D) ）のツァリス分布 P0 (E ) / 1=f1 + (q ; 1)E=T g1=(q;1) の形をし ている。ツァリス統計はボルツマ ツァリス統計の次数 q の意味は ン統計の一般化として最近議論されている 8]。 よく分からないことが多いが 、今の場合は２種類の雑音の強度できまる。D = 0 になると元のボルツマン分布にもど る。この分布は y が大きいところではベキ分 布に近づく。成長率のゆらぎでしばらくの間成長が持続すると、指数関数的に y 図 の値が大きくなるためである。この系でも正弦波外力に対する線形応答の虚数部 4 分が 00(!) = ( ; D;)!2 + !2 : (8) 00 G(!) = ;2 !(!) hy2 i (9) となり、外力がないときの y(t) の自己相関関数のフーリエ変換 G(!) との間に揺 動散逸関係 st が成り立つ。 st は定常 ! = 0 的な力に対する応答係数、 y2 は外力がないとき の y のゆらぎ（ 二乗平均）この関係式をシミュレーションで確かめたのが図 3(b) である。D = 0 とすると 、従来の揺動散逸関係になる。 h i ラチェットモデル ｂ. 熱力学の法則によれば 、一定温度の熱平衡系から 熱エネルギーを取り出して 仕事に変換することはできない。高温源から熱エネルギ ーを取り出して、その一 部を低温源に捨てて、残りを仕事に変換するのが熱機関である。理想的熱機関が カルノーサイクルであるが 、ファインマンはその教科書「ファイン物理学」で熱 機関（ カルノーサイクル ）のわかりやすい例とし てラチェットモデル（ 爪車と歯 止の系）を提案している 9]。このモデルは仮想的なモデルとし て提案されたもの であるが 、近年、分子モータ−（ ミオシン 、キネシンなど ）の運動のメカニズム と関連しているのでないかと考えられ 、いろいろなラチェットモデルが研究され ている 10]。いずれも熱ゆらぎの中で一方向に運動が生成される機構を考えてい る。われわれは元々のファインマンの本に書かれているファ インマンラチェット をそのままシミュレ ーションするようなモデルを提案した 11]。 １次元の周期的ポテンシャルの中の粒子の運動を考え、荷重 L が懸かってい るとする。さらに熱ゆらぎに対応するランダムノイズを受けているとし て、ラン ジュバン方程式をたてる。 @x = ; @U (x) ; L + (t) (10) @t @x (t) は h (t) (t0 )i = 2T (t ; t0 ) を満たす別の白色ガウス雑音である。ポテンシャ ル U (x) は図 4(a) のような鋸歯状ポテンシャルを仮定する。ポテンシャルのピ ー クが 2w で x = x0 から登り坂になっている。x = 2n（ポテンシャルに 2w の ギャップがある位置）に、ゲート（関所）を置き、粒子は左から右には自由に行け るが 、左から右への移動は確率 P0 = e;2w=T0 でしか許されないものとする。T と T0 が２つの熱源の温度に対応する。フォッカープランク方程式の定常解から 、 定常分布が求まる。フォッカープランク方程式の定常解は図 4(a) で x0 = 0 の時は P (x) = f=(w + L) + e;(w+L)x=T (11) R f と は境界条件と正規化条件 dxP (x) = 1 から求まる。この確率分布を図 4(b) に示す。確率の流れがない、すなわち詳細釣り合いが成り立つ条件は 、 2(w + L)=T = 2w=T0 5 (12) 0.1 0.5 0.4 0.3 v P(x) ’0’ ’T’ ’1.5T’ 0.05 0 0.2 -0.05 0.1 0 0 1 2 3 x 4 5 6 -0.1 -0.5 -0.3 -0.1 0.1 L Figure 4: (a) ポテンシャル U (x) (b) 定常分布 P0(x) (c) w = 0:6 T = 0:5 0,0.5,0.75 の時、荷重 L と平均速度。 で温度 T0 が である。この状態では 、温度 T の熱浴から 2(w + L) の熱エネルギ ーを受け取 り、L のおもりを 2 だけ持ち上げる。仕事効率は、 = 2(2wL+ L) = T ;T T0 すなわち、2 つの温度 T と T0 の熱浴のする最大仕事効率（カルノーサイクル ）に なっている。この条件を満たさない場合は、一般に f = 0 で 、粒子は平均速度、 2f でゆらぎを伴いながら一方向に運動する。この平均速度は一般的には、ファ インマンの本にある v e;2(w+L)=T e;2w=T0 とは異なり、ポテンシャルの詳 細な形に依存する。図４ (c) は x0 = 0 w = 0:6 T = 0:5 で温度 T0 が 0,0.5,0.75 の時、荷重 L を変えたときの粒子の平均速度をあらわす。 筋肉の分子モーター、アクチンミオシン 系が上のようなゆらぎ が重要な役割 を果たすという確率モデル説は日本の大沢先生たちによって提唱され 、柳田先生 等による１分子計測によりある程度確かめられたとされているが、ハクスレー等 のボ ートこぎ 説にはじまる力学モデル説との間の論争に決着がついたかど うかは よくわからない 12]。図 4(c) で T0 > T で L < 0 の場合を考えると、x0 で強く熱 励起され 、後は 、勾配に従って自然に負の方向に進む。負の方向に進みながらだ んだん上の方に上昇し 、おもりを持ち上げる。このパラメータ領域では 、ゆらぎ を含みながらもほぼ力学的に進むので 、確率モデルと力学モデルの折衷案になっ ているように思うのだがど うであろうか。 上のラチェットモデルは非平衡系の簡単なランジュバンモデルでもある。非平 衡系に対して一般的に成り立つ法則はあまり多く知られていないが 、その一つに Gallavotti-Cohen の揺動定理がある 13]。これは 、 時間 の間にエントロピ ーが だけ増大する確率と、 減少する確率の比 P ()=P ( ) が が大きくな ると ; 6 / ; ; ; P () P (;) = exp( ) (13) となるという関係である。非平衡系ではエントロピーは増大するが有限時間では ゆらぎによって減少してもおかしくない。その確率を見積もった関係式と解釈で 6 0.3 0.5 0.0025 500 0.002 400 0.0015 300 200 0.001 100 0.0005 0 -0.005 0 0 0.005 0 0.001 0.002 Figure 5: (a) エントロピー生成率 の分布 (b) lnfP ()=P (;)g= と の関係 きる。われわれの系では粒子がド リフト 運動をし て x = 0 から 2n まで移動した とすると、温度 T の熱浴から Q = 2n(w + L) の熱を吸収し 、温度 T0 の熱浴に Q0 = 2nw の熱を放出したことになる。このときの両熱浴のエントロピー増大は = 2n((w + L)=T w=T0) となる。ド リフトの確率分布をシミュレーション で求めたのが図 5(a) でこれから ln P ( )=P ( ) = と の関係をグラフ化し たのが図 5(b) である。確かに 、傾きがほぼ１の直線にのっており、揺動定理が成 り立っていることが分かる 14]。 ; f ; g ３ イジングモデルとカオス結合系 イジングモデルは磁性体の簡単なモデルであり、特に二次元イジングモデル はオンサーガ−により厳密解が求められ 、相転移の存在や相転移点近傍のスケー ル則など の相転移の統計力学において重要な位置を占めてきた 15]。 モンテカルロ法など で イジングモデルをシミュレーションするときには 通常 乱数を発生させ、確率的にスピンを反転させる。イジングモデルのエネルギ ーは X si sj (14) si は i 番目の格子点のスピンで 、上向き（ si = 1) または下向き（ si = ;1) をと る。和は最近接格子点で取り、J は相互作用の強さを表す。このエネルギーに対 するボルツマン分布 P / exp(;E=T ) が平衡状態で実現するように遷移確率を考 える。二次元正方格子を考え、スピンの配位が fmij g をとる確率 pn(fmij g) に 対するマスター方程式は X pn(fmij g) = pn;1(fm0ij g)w(fm0ij g ! fmij g): (15) E=J mij 0 と書ける。この式はスピン変数が遷移確率 w( m0ij mij ) にしたがって確 率的に変化することを表している。たとえば 、遷移確率とし て w= Y ij f g!f g 1=21 + tanhfK=4 mij (m0i+1j + m0i;1j + m0ij +1 + m0ij ;1)g]: (16) 7 を選べば 、詳細釣り合いが成立し 、平衡分布 peq (fmij g) / expfK=8 X ij mij (mi+1j + mi;1j + mij +1 + mij ;1)g (17) が得られ る。この分布は２次元イジングモデルの K=8 = J=T とおけば 熱平衡 分布と等価になる。２次元イジングモデルは厳密解が知られている。無限系では P K = Kc = 2ln(1 + 2) で相転移がおこり、１格子あたりの平均磁化 Sij =N は p m = 1; 1=8 1 K > Kc sinh4 K=2 (18) で与えられる。 われわれは 、この二次元イジングモデルを決定論的な結合カオス力学系で実 現できることを示した 16]。要素の力学モデルとして、ベルヌイ写像 Xn+1 = 1 +2 a (Xn + 1) ; 1 Xn+1 = 1 ;2 a (Xn ; 1) + 1 (19) を用いる。図 6(a) にベルヌイ写像を示す。a は 0 < a < 1 のパラメータ。ベルヌイ 写像はカオスを示し 、スピン変数 sn を Sn = sgn(Xn+1 Xn ) 、すなわち Xn < a に対して Sn = 1 、Xn > a に対して Sn = 1 と定義すると、スピン変数の平均 値は a で時間相関が０となることが分かる。これはベルヌイ写像系と裏と表の出 る確率が異なるコイン投げとが等価であることを示している。 このベルヌイ写像系をもとに 結合写像格子系を構成する。パラ メータ aij n が １ステップ前の最近接格子点のスピン変数によって決まるとする。 各正方格子点 上にベル ヌイ写像を置き、そのパラメータ a を隣接格子点のスピン変数によって 次のように決まるとする。 ; ; K i+1j i;1j ij +1 i:j ;1 aij n = tanhf 4 (Sn;1 + Sn;1 + Sn;1 + Sn;1 )g (20) こうすると、スピン変数の時間発展が上にのべたマスター方程式そのものになる とことが分かる。従って、この結合写像系を時間発展は決定論的であるが、スピン 変数のみ見ると確率的なイジングモデルと等価な振る舞いをする。実際、100 100 の正方格子で結合写像系を計算した、結果を図６に示す。図 6(b) はある時刻のス ピン変数の空間分布を示している。K = 1:7 と相転移点の近くなので、スピンの 長距離相関があらわれ 、クラスター構造が見えている。図 6(c) は平均磁化と K の関係を示したもので 、有限サイズ効果で相転移がなまっているが 、理論に近い 平均磁化が得られている。 ベル ヌイ写像は時間反転できない写像であるが 、パイコネ変換の変形版であ る２次元写像 Xn+1 = 1 +2 a (Xn + 1) ; 1 8 (a) (b) (d) 1 Y 0.5 0 -0.5 -1 -1 -0.5 0 X 0.5 1 Figure 6: (a) ベルヌイ写像。(b) K = 1:7 におけるスピン変パターン。 (c) 結合変形パイコネ変換の (Xn Yn) 平均磁化と K のグラフ。実線は理論値。 (d) のプロット。 Yn+1 = 21 (Yn + 1) for ; 1 < Xn < a Xn+1 = 1 ;2 a (Xn ; 1) + 1 (21) Yn+1 = 12 (Yn ; 1) for a < Xn < 1 を用いると。時間反転できる写像になる 17]。この写像を使った結合写像格子系 を同様に作ることができる。 2 (X ij + 1) ; 1 n 1 + aij n Ynij+1 = 12 (Ynij + 1) for ; 1 < Xnij < aij n 2 (X ij ; 1) + 1 Xnij+1 = n 1 ; aij n ij Ynij+1 = 12 (Ynij ; 1) for aij n < Xn < 1 Xnij+1 = こうすると 、スピン変数の過去のデ ータは 9 Y (22) 変数に蓄積されていく。この過程 を考えるため、y = (1 + Y )=2 とおく。y は 0 1] 内の実数である。y を２進数 と表せる。y の時間発展はその２進数表現では 、 展開すると、y = 0:s1s2 s3 0:s1 s2 s3 0:s0 s1 s2 s3 となる。つまり、過去のスピン情報は（ スピン＋１を １，スピン−１を０に変換して）２進数表現の下の方の桁に保存されていく。時間 発展は桁を１つずつ下げて、現在のスピン状態を小数点第１桁に加えることに対 応する。逆に言うと、ある時刻ｎの Ynij を詳細に見れば 、過去のスピンの時系列 がわかるということである。図 6(d) は (Xnij Ynij ) のスナップショットを表してお り、X 方向はほぼ一様であるが Y 方向には微細な構造がある。この構造は、相転 移後なのでスピンが上を向く確率と下を向く確率が異なるため、y = 0:s1s2 s3 で表すと、si が１をとる確率が 1/2 からずれているためである。確率分布にフラ クタル性（ マルチフラクタル性）が現れているといってもよい。このモデルは作 り方から分かるように、２次元正方格子モデル以外でも一般的なイジングモデル の計算に使うことができる。ただ 、最近は乱数発生が容易に行えるので 、通常の モンテカルロシミュレーションと比べて優位かど うかは分からない。カオスとい う決定論に従う系でも統計力学系が構成できることを示せたのはおもし ろいこと ではないだろうか。時間反転可能なモデルではある時刻の状態が分かれば過去も 未来も完全に決定されているのに確率的なふるまいをするのである。 ! ４．相転移ダ イナミクスとギンツブルグ−ランダウ方程式 a. ギンツブルグ−ランダウモデルとド メインの成長と制御 イジングモデルのような実数の秩序変数をもつ系の相転移はギンツブルグ− ランダウモデルで記述できる。系の自由エネルギ−を磁化のような秩序変数 u で ベキ展開する。系の対称性から x の偶数ベキのみ現れる。 Z F = ;au2 =2 + bu4=4 + D=2(ru)2dx (23) 相転移のダ イナミクスは 、非保存系では du = ; @F = (au ; bu3 + Dr2 u) dt @u (24) となる。この式を時間に依存するギンツブルグ−ランダウモデルという。相転移 温度以上では a < 0 で u = 0 の一様な解に近づく。相転移温度以下では a > 0 で 定常解 u = (a=b)2 に近づく。a を相転移点からのずれで表すと、a = a0 (Tc T ) なので、秩序変数は (Tc T )1=2 で増大する。この 1/2 乗則は平均場近似の指数 と等価である。系の対称性から、秩序変数は正と負の２つの安定点を持つ。力学 系の言葉では双安定系である。ほぼ一様な系から、急冷して相転移以下にすると、 u が正のド メインと負のド メインができその競合が起こる。ド メインのサイズが 大きくなると、空間次元が２次元以上では表面張力で支配されるようになり、ド メインの成長は時間 t の 1/2 乗となることが多い。このような相転移のダ イナミ クスはこれまで詳し く研究されてきた 18]。 われわれは、ギンツブルグ−ランダウ系のような双安定系で大域フィード バッ クによってド メインの大きさがコントロールできることを示した 19]。非線形系 の制御はカオス制御など 、応用を視野に入れて、最近いくつか研究がある。１次 ; ; 10 (a) (b) 1 400 0.5 t u 300 200 0 -0.5 100 -1 0 0 50 Figure 7: (a) 100 x 150 200 0 50 100 x 150 200 ２ド メイン共存系 (b) ３ド メイン共存系 元系での、大域フィード バックモデルは次のようになる。 @u = u ; u3 ; K (hui ; c) + @ 2 u (25) @t @x2 RL 右辺の第３項がフィード バック項で 、K > 0 とする。hui = 1=L 0 udx は u の平 均値を表し 、c はある定数である。u = 1 の双安定系なので 、u = 1 と u = ;1 の２つのド メインに分かれる。もし 、hui が c より大きいと、u = 1 のド メイン と u = ;1 のド メインを分かつ界面が u = ;1 のド メインを大きくする方向に動 き、hui が小さくなるようになる。逆に 、hui が c より小さいと、界面が u = 1 ド メインが成長する方に移動し 、hui が大きくなる。いずれにしても、hui は c に近 づく。その結果、u = 1 のド メインのサイズは (1 c)L=2 に落ち着く。逆に言 えば 、u = 1 のド メインサイズを L にしたければ 、c = 2 ; 1 と置けばよい。 このように、任意のド メインサイズの共存系を大域フィード バック系で実現でき る。ド メインサイズをコントロールできるわけである。 上の系では双安定系を用いたが３重安定系の場合は 、２つの大域フィード バッ ク項を導入すれば 、３つの相を好きなサイズ比にすることができる。 たとえば 、 ５次のギンツブルグ−ランダウモデル @u = ;au + u3 ; u5 + e + @ 2u (26) @t @x2 を用いると、３重安定点が得られる。特に 、３重点 a = 3=16 e = 0 では 、３つ p の安定点 u = 0 3=2 の自由エネルギ− F = au2=2 ; u4=4 + u6=6 + eu が等し くなる。この、３つの安定状態のド メインサイズをコントロールするために、２ つの大域フィード バック項を導入する。 @u = u3 ; u5 ; (hui ; c ) ; (hu2 i ; c )u + @ 2 u 1 2 @t @x2 (27) もし 、この２重フィード バックがうまくいけば、安定定常状態では、 u = c1 u2 = c2 + 3=16 が得られる。u = 3=2 のド メイン長を l 、u = 0 のド メイン長を l0 p 11 hi h i (c) (a) 80 100 350 300 60 250 200 40 50 150 100 20 50 0 0 0 10 20 30 40 50 60 70 0 20 40 60 80 0 x Figure 8: (a) p 100 200 x ２ド メイン共存系 (b) ３ド メイン共存系 とすると、 u = 3=2(l+ l; )=L = c1 u2 = 3=4(l+ + l; )=L = 3=16 + c2 。逆 にいうと、こうなるような c1 c2 を選べば 、望むサイズ比の３つのド メインの安 定共存系が得られる。図 7(a) は２つのド メインの安定共存状態を、図 7(b) は３ つのド メインの安定共存状態を表す。 hi b. ; h i 階層的バンチングのモデル 合金の相分離過程など では 、相転移転移以下で 、一様混合相から２相に分離 する。その際、金属元素の量は保存するので 、相転移のダ イナミクスは 、保存系 のギンツブルグ−ランダウ方程式（ Cahn-Hilliard 方程式）になる。 du = ; r2 (au ; bu3 + Dr2 u) dt (28) ２相の体積分率が大きく違うときは、液滴状の小さなド メインが生成され 、それ らが成長するプロセスが見られる。多くの小さなド メイン間である種の競合が起 こる。２次元以上では表面張力効果により、より小さなド メインが蒸発し 、その 成分が拡散により大きなド メイン付近に輸送され 、大きなド メインがさらに大き く成長する。Lifshitz と Slyozov はこの競合成長過程を詳し く解析し 、ド メイン 数は 1=t で減少することを見いだした 20]。 表面物理の分野でも似たような現象が見られる。(1,1,1) 面のような小さな指 数の面からわずかにずれた面を微斜面という。(1,1,1) 面からずれた微斜面は本来 エネルギー的に損なので 、(1,1,1) 面にもど ろうとする。しかし 、界面の平均角度 は保存されるので 、保存系のダ イナミクスになる。微斜面が (1,1,1) 面と大きな角 度の界面に相分離する傾向がある。特に 、Si(1,1,1) 面に近い微斜面に外から表面 に沿って電流を流してやると、微斜面が不安定化し 、微斜面の相分離が起こるこ とが実験的に知られている。ミクロに見ると、原子の層の段差（ ステップ ）が束 ねられて大きなステップ群が生成され 、それが大きな角度の斜面になる。ステッ プ間隔が長い部分は平らな界面に対応する。平らな面と大きな角度の斜面への相 分離は 、ステップから見ると、ほぼ等間隔に並んでいたステップがバンチングす 12 300 400 る不安定化と見ることができる。実験によると、平らな面の長さが t1=2 で増大す ることが観測されている 21]。このような、時間のベキ則で進むバンチングを階 層的バンチングと呼ぶ。われわれは 、簡単な階層的バンチングの１次元モデルを 作った 22]。ステップ位置を yi とし 、ステップの運動方程式として dyi = ;f (y ; y ) + f (y ; y ) (29) i+1 i i i;1 dt f (l) はステップ間の力で、簡単のため f (l) = l + a=l3 とする。l 十分近いときは f (l) 1=l3 の反発力が働く。この反発力のため、ステップは融合することはない。 l が l が大きいときの力は計算を簡単にする l に比例するとしたが 、実際は必ずし もそうなっていない。またこのモデルでは一方向向きのド リフト 運動は考慮して いない。この l が大きいときの l に比例する項のため、等間隔にステップが並ん だ解 yi = l0 i は不安定化する。直接数値計算すると、まず隣ど うしがバンチし あ い、ステップのペアができる。ペアど うしがさらにバンチしあって大きなクラス タになる。さらに 、クラスタど うしが束ねられといったプロセスがつぎつぎ 階層 的に起こることが分かった。図 8(a) にその時間発展を示す。この系では 、平らな 面の幅に対応するステップ 間隔は図 8(b) にあるように、t に比例して増大するこ とが分かった。 ステップ間隔 li = yi+1 yi に関する式を連続体近似すると、 ; @l = ; @ 2 ; @U (l) + D @ 2 l @t @x2 @l @x2 (30) と、Cahn-Hilliard タイプの式になる。しかし 、ポテンシャル U (l) = 1=2l2 + a=2l;2 は、双安定ポテンシャルではなく、l(t) の時間発展はギンツブルグ ーラン ダウモデルとはかなり違ったものとなる。ギンツブルグ−ランダウ系では１次元 では双安定のド メインができ、ド メイン間の競合が非常にゆっくりと進む。今の 系では 、l(t) が小さいバンチし たクラスタの状態と、l(t) が大きい平らな面に分 かれる。l(t) の大きなところはパルス状に見え、その幅はパルスごとにあまり変 わらず、その高さが様々に分布する。図 8(c) にそのパルス群の競合的時間発展を 示す。 小さなパルスは、より小さくなり、その成分が拡散し 、大きなパルスに吸 収される。この競合成長過程の駆動力は表面張力ではないが 、Lifshitz-Slyozov が 解析したようなプロセスが１次元系でも起こっていると解釈できる。この連続モ デルでも、パルスの高さが t に比例して増大する。 ; ５．振動子集団の引き込み相転移 a. 位相振動子モデル 非線形ダ イナミクスの典型例に 、リミットサイクル振動がある。定常状態が ホップ分岐すると、リミットサイクル振動があらわれる。力学の講義で習う単振 子など とは異なり、時間が十分たつと、初期値によらず一定の振幅、周期のリミッ トサイクルに近づく。散逸系のアトラクタとし ての振動である。ベルの振動、電 気回路の発振や心臓の拍動など もリミットサイクル振動である。リミットサイク ル振動の特徴的な性質とし て引き込み現象（ 同期現象）がある。２つのわずかな 周期の異なるリミットサイクルを結合させると、その周期がそろうという現象で 13 Figure 9: (a) = 0 の振動数分布 (b) = =4 の振動数分布 Figure 10: ジョセフソン接合系 ある。最初、ホイヘンスが２つの時計の同期現象を発見したといわれる。心筋の ペースメーカー細胞は同期して振動し 、拍動のペースを作っている。脳波も脳の 神経細胞の弱い同期現象によって生じたものといえる。強く同期発振するとてん かんのような制御できない状態になり具合が悪いが 、脳の中の視覚野や海馬と呼 ばれる部分では同期現象が脳での情報処理にも役立っているのではないかという 説もある。 多数のリミットサイクルの同期現象は Winfree によって最初研究され 、後に 蔵本先生がシンプルなモデル（ 位相モデル ）を提案した 23]。そのモデルは N di = ! + K X sin(j ; i ) (31) i N dt j =1 である。K は結合強度で 、各振動子の固有振動数 !i はランダ ムにばらついてい るとする。平均場結合の位相モデルは 、熱平衡系の秩序無秩序相転移のような相 転移を示す。実際、磁性体の平均場近似と似た方法で解が求まる。K < Kc 以下 では 、各振動子はばらばらに振動する。K > Kc では 、何割かの振動子がその平 P 均振動数 !0 で同期して振動する。秩序変数を = (1=N ) exp(ij ) とすると、 14 (a) 20 10 10 0 0 y y 20 -10 -10 -20 -20 (b) -10 0 x 10 Figure 11: 20 -20 -20 -10 0 x 10 20 結合レスラー系 jj は K > Kc で０から (K ; Kc )1=2 で増加する。これは、朱に交われば赤くな るということわざ のように、集団内にいるとその集団のペースに巻きこまれてし まうことを表している。上の位相モデルとその変形版はその後詳し く研究されて いる。 上のモデルでは 、平均振動数に引き込まれる傾向があるが 、実際の同期現象 では 、高い振動数に同期されやすい。われわれはその効果を取り入れたモデルを 考えた 24]。 N di = ! + K X sin(j ; i ; ) i N dt j =1 (32) 位相の変調項 のため、同期した結果の振動数は 、平均振動数より高くなる。こ のモデルも 、平均場近似で解ける。図 9(a) は が０、すなわち元のモデルでの 同期後の振動数分布である。振動数がそろった振動子が多数生じ るので 、振動数 分布にδ関数状の部分が生じ る。その同期振動数は元の振動数の平均値になって いる。図 9(b) は = =4 の系の振動数分布で 、平均振動数より高い振動数に引 き込まれていることが分かる。 ジョセフソン 素子も外から 電流を流すと 、臨界電流以上で発振する。これは 交流ジョセフソン効果と呼ばれる。ジョセフソン素子を LCR 回路の中で多数直 列につなぐ と、この発振も同期する。図 10 のような一様なジョセフソン 素子を 多数直列につないだ系は電圧の標準を決定するためにも利用されている。ランダ ムなジョセフソン素子の直列結合系は 、ジョセフソン素子に対する RSJ モデルを 用いて、Wiesenfeld らによって近似的に解析された。実効的な容量（キャパシタ ンス）が無視できるとき、上の を含んだ位相モデルで近似できることが示され た 25]。われわれはキャパシタンス成分を入れた数値計算も行い、(32) 式の の 値がキャパシタンスの効果で正負の符号を含めて大きく変化することを示し た。 また、キャパシタンスをいれると、単独の素子でも発振と定常状態が双安定にな るので 、多数ジョセフソン素子を結合した系では同期の秩序変数が同じパラメー タでも連続的に変わりうることを示した 26]。 15 Figure 12: 振動子格子 (a) １次元系 (b) ２次元系 b. カオス振動子集団の同期 位相振動子は単独では正弦波的な振動をする。均一な系では 完全に同期するが 、 自然振動数の不均一さや熱ゆらぎが入ると、同期非同期転移が起こる。カオス振 動子の場合は単独でも不規則運動をするので 、均一な系でも、結合を強くしてい くと、ばらばらに振動している状態から同期するようになる。これをカオス同期 という。強く結合させると、完全にカオス運動がそろう現象が見られるが 27] 、弱 い結合の場合、振動の振幅は一致しないが位相がほぼそろうことがある。これを カオスの位相同期という 28]。レスラー系のようなパワースペクトルにかなり鋭 いピ ークがあるカオス系でよく見られ る。われわれは同一のレ スラー系を多数、 平均場的に結合した系の数値計算を行い、結合を強くしていくと、秩序無秩序相 転移のようなふるまいが見られることを示した 29]。図 11(a) は結合が弱いとき で、ある時刻での N 個のレスラー系の (Xi Yi) のプロットである。ばらばらに運 動しているのでレスラーアトラクタ上に均一に分布している。図 11(b) は結合が 強いときで、位相同期しているために 、位相分布は狭くなり、(Xi Yi) がレスラー P アトラクタのある扇形領域内に分布している。この状態では平均値 (1=N ) Xj はきれいな振動を示す。 c. 振動媒質の同期現象 16 Figure 13: 振動子媒質の非同期転移 (a) 相図 (b) 非同期時の振動の様子 位相振動子を格子上で局所的に結合させた有限次元系での引き込み現象を調 べた 30]。モデル方程式は di = ! + K X sin( ; ) i j i dt j (33) となる。和は最近接格子点で取る。１次元系、２次元系では同期したクラスタ− が局所的にできる。結合強度 K を大きくし ていくと、その同期クラスタ−が図 12 のように成長する。図 12(a) は１次元系での K を変化させて時の振動数のプ ロファイルを描いてある。同期クラスター内では振動数がそろっており、隣のク ラスタ−との間に振動数の飛びがあるのが見て取れる。２次元系では同期した隣 り合う振動子をボンド でつないだ図を描いてある。K とともに 、パーコレ ーショ ンのように同期クラスターが成長する。パーコレ ーションは２次元系で臨界現象 を起こすが、この振動子系では２次元系で引き込み相転移を生じないようである。 熱平衡系の相転移との類推で考えると、系の次元をさらにあげ るとある臨界次元 以上で、平均場モデルのような引き込み相転移が生じ るものと考えられる。 格子上のように振動子が離散的に離れていると、同期非同期現象は無理なく おこるが 、偏微分方程式系のような連続系では同期非同期転移は起こりうるのか という疑問が生じ る。非同期というのは振動数の不連続であり、空間連続性とお りあう必要がある。われわれは、非一様なギンツブルグ−ランダウ方程式を用い てこの非同期化を調べた 31]。モデルは @W = (1 + i!(x))W ; jW j2W + D @ 2 W (34) @t @x2 である。0 < x < 2 の領域で考え、W は W = Rei(t) のような複素変数。自 然振動数とし て 0 < x < 1 で !(x) = !0 、1 < x < 2 で !(x) = ;!0 とす 17 1 150000 0.01 100000 0.0001 50000 0 1 10 100 1000 10000 0 Figure 14: (a) バネ−ブロックモデル (b) (c) スリップ 発生の時間発展 臨界時のスリップの大きさ分布 る。結合がなく、D = 0 の時は 、0 < x < 1 で W = ei!0 t 、1 < x < 2 で W = e;i!0 t となり x = 1 で完全に不連続になる。D が０でないと、!0 と D の相 図が図 13(a) のようになる。D も !0 も大きいと、連続性と、非同期性のおり合 いがつかず W = 0 となる。いわゆる、振動停止である。!0 が大きい状態で D を 小さくしてゆくと、W = 0 の解がホップ分岐し 振動を始める。このときの振動は W (x) cos(!t) + iA(cos !t + ! sin !t) cos(x=2) のようになる。x < 1 では右回 り、x > 1 では左回りの楕円軌道になる。x = 1 に近づくにつれ 、楕円がだんだん つぶれていき、x = 1 で水平に振動する。この時の振動の様子 (ReW (x) ImW (x)) を図 13(b) に示す。x = 1 では周期的に原点（ 位相が定義できない点）を通過し 、 W の空間連続性と振動数の非同期性が両立する解になっている。 d. スティックスリップ 振動と自己組織化臨界現象 リミットサイクル振動は、ホップ分岐近傍では正弦波振動に近い形をするが 、 一般には変動が小さいところと急に変動が大きくなるところが繰り返す弛緩振動 型の振動が多い。お寺にあるしし おどしのような振動は典型的な弛緩振動である。 ある閾値まではほとんど 動かず、閾値を越えると急に変化し 元の状態にもど る。 神経細胞の発火振動もこのタイプである。一様な滑りが不安定化して発振すると きには、滑り状態と静止状態を繰り返すスティックスリップ振動が生じる。チョー クを黒板にこすったときのギ−といういやな音が出ることがあるが 、あの振動も そうである。地震も断層が滑って起こるので、スティックスリップ振動を起こす。 地震が周期的に起こる傾向があるのはこのせいである。地震の規模に関し ては 、 その発生頻度がエネルギ ーのベキ乗則に従うというグーテンベルグリヒター則が 18 知られている。ベキ分布ではガウス分布の標準偏差のような特に決まったスケー ルがないことが特徴である。熱平衡統計力学では 、空間相関の長さのスケールが 発散する、２次相転移点や臨界点でベキ分布が観測される。ゆっくりと駆動され る非平衡系では 、特に臨界条件を外から与えてやらなくても、自然に臨界状態に 向かう傾向があることが Bak 等により予測され 、この現象は自己組織臨界現象と 呼ばれている。Bak 等の砂山モデルでは 、砂山が崩れて起こるなだれの大きさが ベキ分布に従う 32]。地震のグ ーテンベルグ−リヒター則が自己組織臨界現象と 見なしてよいかは議論があるようである 33]。Langer 等はたくさんのブロックを バネで１次元的につないだモデルを数値計算し 、一度にスリップするブロックの 数の分布を調べ、その数が比較的小さいときはベキ分布に従うが、大きさスケー ルのスリップは比較的大きさがそろっていることを示した 34]。 われわれも、すこし 変形したバネーブロック系の計算を行った 35]。図 14(a) にバネ−ブロックモデルの模式図を示す。摩擦面の状態が流動体的であるか固体 的であるかに関するオーダーパラメータを導入し 、オーダ ーパラメータがずり速 度で変化するモデルを用いた。具体的には d2 xi = K (x ; 2x + x ) ; F ( v ) i+1 i i;1 i i dt2 di = (1 ; ) ; v + D( ; 2 + ) i = 1 2 N (35) i i i i i+1 i i;1 dt i が摩擦面の流動状態を表すオ−ダ ーパラメータである。N = 1 の１ブロック系 では 、規則的なスティックスリップ振動をする。N が大きな系では 、パラメータ をかえていくと、ずるずると一様にずれていく状態から、スティックスティック 状態に相転移的に変化することが分かった。相転移点で、スリップの大きさ分布 がきれいなベキ分布になり、相転移点からずれていくと、Langer 等が見いだした ように、小さなスケールではベキ的だが大きなスリップはかなり周期的に同程度 の大きさで起こることが分かった。すなわち、このモデルでは自己組織臨界現象 は起こらず、通常の相転移系に近い振る舞いが見られたということである。もっ とも、この系は力学系なので、そのランダムなふるまいは散逸力学系のカオスに 起因するものと考えられる。図 14(b) に臨界点でのスリップの大きさ分布を示す。 図 14(c) はスリップの発生の時間経過を示す。Bak 等のモデルで自己組織臨界現 象が見られたのはモデルで仮定されている粒子数の保存則によるところが大きい のではないかと思う。 ６．散逸構造と位相ダ イナミクス a. 振幅方程式 空間的なひろがりを持った系で、熱平衡から離れると、熱平衡状態と質的に異 なった状態に転移する。これを非平衡系の散逸構造と呼ぶ 36]。上下に温度差を 付けたときに生じ るベナール対流はその例である。温度差が小さいときは熱伝導 で熱は伝わるが 、対流状態では 、流体に流れが生じ 流れにのって熱が伝わる。化 学反応の１つである BZ 反応と呼ばれるマロン酸の酸化反応では 、化学反応が時 間的に振動したり、波動が伝搬したりする。一様定常状態から 、温度差など の外 部パラメータを変化させたとき、まず最初に起こる不安定性には、不安定化する 波数が０か有限かと、振動数が０か有限かの組み合わせで大きく４つに分けられ 19 Figure 15: 向きの異なる進行波の周期パターン jWr j と Wl j (b) ２次元杉綾模様 (a) １次元系 点線と実線は る。１つ目は、波数 k = 0 で振動数 ! = 0 の場合で 、異なった一様定常状態に転 移する。ギンツブルグ−ランダウ方程式で表されるような熱的相転移もこの範疇 に入る。２つ目は 、波数が有限で 、振動数が０の場合で 、ベナール対流のように 定常な空間パターンが生じ る。３つ目は 、波数が０で 、振動数が有限の場合で、 一様振動状態が生じ る。空間的な広がりのあるホップ 分岐である。４つ目は、波 数、振動数ともに有限の場合で 、波動パターンが形成される。 転移点近傍では、新しい状態を特徴づける変数が小さいので 、様々な近似法で 秩序変数の振幅に対する方程式が導ける。この方程式を振幅方程式と呼ぶ。いづ れも、４章でのべた 、ギンツブルグ −ランダウ方程式によく似た形をし ている。 一様定常解への転移に関する振幅方程式は 、ギンツブルグ−ランダウ方程式その ものである。 @W = W ; W 3 + r2W @t (36) 定常パターンへの転移に関する振幅方程式は 、パターンの対称性によって異なる が 、一番単純な熱対流のようなロール状のパターンへの転移に関しては 、NewellWhitehead 方程式が得られる。簡単のため、ロール方向のみ考えると、複素変数 W に関するギンツブルグ−ランダウ方程式になる。 @W = W ; jW j2W + @ 2 W @t @x2 (37) 一様振動状態への転移に関しては 、係数も複素数になる複素ギンツブルグ−ラン ダウ方程式が得られる。 @W = W ; (1 + ic )jW j2W + (1 + ic ) @ 2 W 2 1 @x2 @t (38) 波動状態への転移に関しては 、右に進行する波動と左に進行する波動が生じ 、そ の間の結合も考えると、結合複素ギンツブルグ−ランダウ方程式が得られる。 @Wr = W ; f(1 + ic )jW j2 + (g + ig )jW j2gW + (1 + ic ) @ 2 Wr r 2 r 1 2 l r 1 @x2 @t @Wl = W ; f(1 + ic )jW j2 + (g + ig )jW j2gW + (1 + ic ) @ 2 Wl (39) l 2 l 1 2 r l 1 @x2 @t 20 これらに関しては 、Cross-Newell のレビュ−が詳しい 37]。波動パターンに関し ては 、空間微分項を無視すれば 、g1 > 1 の時、進行波解 Wr = 0 Wl = 0 または Wr = 0 Wl = 0 となり、g1 < 1 の時、定在波動解 Wr = Wl が得られることが分 かる。しかし 、空間微分項の影響を考慮すれば 、進行波・定在波間の転移もかな り複雑になる。われわれは 、ある場合には 、右進行波と左進行波が空間的に周期 的に（ 交互に ）あらわれる場合があることを示した 39]。液晶のような異方的な 系では 、向きの違う斜め進行波解が交互にあらわれ 、いわゆる杉綾模様があらわ れる。液晶の電気対流ではきれいな杉綾模様が実際に観測されているが 42] 、模 様のできるメカニズムは上のものとはすこし 違うようである。図 17 は (39) 式を 計算して、空間的に左右の波が交互に優勢になった例である。 6 6 b. 位相方程式 一様定常状態への転移以外の場合では 、振幅方程式には 複素変数が使われて おり、解とし て W (x t) = exp(i) のような W = 1 から位相がずれた解を持ちう る。これは 、空間パターンの場合は位置のずれ 、一様振動の場合は時間のずれに 対応する解である。これを一般化して、位相状態がわずかにゆっくりと変動する 摂動に関して、位相方程式が導出される。線形位相方程式は 、@=@t = D@ 2 =@ 2 のような拡散方程式の形になる。D が正の場合は、そのような摂動はゆっくり緩 和し 、規則的な空間パターンや時間振動にもど る。しかし 、場合によっては D が 負になることがある。位相変調が増大する位相不安定性である。不安定性が進行 すると、線形方程式では系の発展を記述できなくなり、非線形位相方程式が重要 になる。蔵本先生は一様振動解の位相変調に対し て Kuramoto-Sivashinsky 方程 式を導出した 23]。 @ = ; @ 2 ; @ 4 + @ 2 @t @x2 @x4 @x (40) @ = ; @ 2 ; @ 4 + @ @ @t @x2 @x4 @x @x2 (41) Kuramoto-Sivashinsky 方程式は時空カオスを示す方程式として詳し く研究され てきた。 同様に、定常パターンに対する位相変調に対しても非線形位相方程式が導出 できる 38]。 この位相不安定性はエクハウス不安定性と呼ばれ 、不安定性が進むと、空間周期 構造の周期が詰まったところと開いたところができてくる。しかし 、この非線形位 相方程式を数値計算すると、不安定性は有限にとど まらず発散することが分かっ た 40]。実際の実験でも、エクハウス不安定性は有限でとど まらず、最終的には 空間周期構造が違う状態に移行する。いずれにせよこの場合位相方程式の記述は 破綻することがわかる。 Kuramoto-Sivashinsky 方程式の場合はうまく解が有限にとど まった。しかし 、 もとの一様振動解の振幅方程式である複素ギンツブルグ−ランダウ方程式を数値 計算すると、位相不安定性が弱いときは Kuramoto-Sivashinsky 方程式で近似で きる位相乱流状態があらわれるが 、不安定性が強くなると、振幅 W が０近くま で変動する振幅乱流状態になる。この状態では Kuramoto-Sivashinsky 方程式で の位相記述が破綻しているようである 40]。 j j 21 5000 4000 t 3000 2000 1000 0 0 Figure 16: (a) 不安定バーガ ース方程式の時間発展 モデルの時間発展 50 100 150 200 250 300 350 x (b) ショックの位置に関する われわれは波動パターンに関する非線形位相方程式も導出した 41]。この場 合、右進行波に対する位相 r と左進行波に対する位相 l の結合方程式が得られ る。この結合方程式は液膜流の不安定性を記述するベニー方程式の結合系のよう な形をし ている。一方向進行波の時は 、結合項がなくなり、ベニ−方程式そのも のになる。 ７．時空カオス a. Kuramoto-Sivashinsky 方程式 不安定バーガース方程式 少数自由度のカオスの研究が非常に盛んに行われた後、空間自由度を含んだ時空 カオスの研究が行われるようになった。しかし 、自由度が大きい時空カオスの本 質はなかなかはっきりし ない。Kuramoto-Sivashinsky 方程式は の x 微分を用 いた形で研究されることも多い。 @u = ; @ 2 u ; @ 4 u + u @u @t @x2 @x4 @x (42) Kuramoto-Sivashinsky 方程式は時空カオスを生成する一番簡単な偏微分方程式 系とし て知られている。様々な系で理論的にこの式が導出されているが 、この方 程式に従う現象の実験例はあまり多くないようである。Sivashinsky が最初導出 した燃焼波面の不安定性がその実験例の１つかもしれない。波数 k 1= 2 の変 動の線形成長率が最大なのでその付近のゆらぎ（セル構造）が比較的よくあらわ れる。セル構造の生成消滅がこの系の特徴ともいえるが 、セル構造自体それほど はっきりししたものではない。非線形項の形は流体力学でよく出てくる移流項の 形をしている。そこで 、時空カオスのト イモデルとし て、Kuramoto-Sivashinsky 方程式が持っている長波長安定性とガ イレ イ対称性を保持した不安定バーガース 方程式の計算をし た 43]。衝撃波に対する流体力学でよく知られたバーガ ース方 程式に長波長成分が不安定になる項を付け加えたものである。バーガ ース方程式 の時間発展と同様、ショックがたくさん生成される。ふつうのバーガース方程式 22 p 3 10 2 8 u 1 6 0 4 -1 2 -2 -3 0 0 100 200 x 300 400 0 0.05 0.1 0.15 k 0.2 0.25 0.3 Figure 17: (a)KS 方程式のショック構造。点線はある時刻での u(x t) 、実線はそ (b) フーリエ振幅の自己相関関数から求めた有効粘性係数 の長時間平均 ではショックは互いに融合しあって数が減っていき、ショック間の間隔が大きくな る。長波長の不安定項を導入してあるので 、ショック間の間隔が広くなりすぎ る と、新たなショックが自発的に生成される。この結果ショックの生成と融合を繰り 返す時空カオスがあらわれる。ショックの幅と、ショック間隔が大きく異なるの で、ショックの生成融合という描像が明瞭に見える。このショックの位置に関して 簡単な方程式をたて 、ショックの生成消滅過程を力学的に表すこともできる。図 16(a) は不安定バーガース方程式の時間発展を表す。図 16(b) は、ショックの位置 に関する簡単化した方程式を数値計算した結果である。２つのショックが近づく と１つに融合し 、離れすぎ るとその中間に新たなショックをできるとし て計算し ている。 Kuramoto-Sivashinsky 方程式の大スケールのダ イナミクス Kuramoto-Sivashinsky 方程式 (KS 方程式）の時空カオスは統計的にも詳し く研 究されている。上でのべたようなセル構造よりもずっと大きい大スケールのふる まいに着目すると、統計的に KPZ 方程式 @ = @ 2 + @ 2 + (x t) @t @x2 @x または 、微分系では noisy Burgers 方程式 @u = @ 2 u + u @u + @ (x t) @t @x2 @x @x (43) (44) と同じようなふるまいをするのではないかという推測がある 44]。ここで (x t) はガウシアンホワイト ノイズである。KPZ 方程式はランダムな界面の成長モデル とし てよく研究されているモデルであり、noisy Burgers 方程式もゆらぎのある流 体力学系としてよく研究されている。繰り込み群を用いた研究により、その長波 23 Figure 18: 流れ場との結合系 ショット (b) 時空カオス (a) ジグザグ振動 ２つの図は違う時刻のスナップ 長のゆらぎがあるスケーリング則に従うとされている 45]。ノイズ項があるとは いえ、大スケールでバーガ ース方程式にしたが うとすれば 、境界で段差を付けた 系では大スケールでの tanh 型のショック構造があらわれるのではないかと予想し て、微分型の Kuramoto-Sivashinsky 方程式を境界条件 u(0) = U0 u(L) = U0 で数値計算してみた 47]。段差 U0 が小さいときは、時空カオスにかき消されるた めか 、u(x t) を長時間平均するとほぼ０になった。U0 を大きくすると、u(x t) の 長時間平均がかなりきれいに U0 tanh (x L=2) で近似されることが分かった。 このショック構造の幅から有効粘性 eff が見積もることができる。数値計算の結 果から 、有効粘性が 、ショックの幅とともにすこしずつ大きくなり、かなり幅が 大きいときは有効粘性が１０程度になることが分かった。この値は 、従来言われ ていた有効粘性の値ともあっている 46]。元の周期境界条件での KS 方程式の長 波長成分のゆらぎ の時間相関係数を計算して、有効粘性を見積もってみると、や はり波長とともにすこしずつ大きくなり、１０程度に近づくことが分かった 48]。 この時空カオス系でも段差を付けるという外力に対する応答と、外力がないとき のゆらぎの性質が対応しているという、線形応答関係が成り立っていることを示 唆しているように思う。 ; ; b. 複合場の時空カオス 熱対流など の定常パターンに対する振幅方程式はリアプ ノフ関数を持ち、定常 解に向かう性質がある。温度差を大きくするとロールのつなぎ 換えや欠陥の生成 消滅を含んだ時空カオスが実験では見られるが 、振幅方程式では時空カオス化し ない。熱対流の実験でカオス化が起こりやすいのは粘性が小さい流体に対してで ある。これらの系では 、大きなスケールの流れが誘起されやすく、流れとロール 構造のダ イナミクスが結合して時空カオスが起こるというモデル方程式が Siggia 等によって提案された。このモデルの数値計算でも欠陥の生成をともなう乱流が 確認された 49]。緩和が遅いゆっくりとした場が、振幅方程式のようなやはりゆっ くりと動く方程式と結合すると、長波長でのカップ リングが強くなり、カオス化 が起こるものと考えられる。液晶の電気対流は甲斐先生等により実験的に詳し く 研究され 、電圧を上げると、Fluctuating Wiliiams Domain とよばれる欠陥乱流 24 200 100 0 0 100 200 Figure 19: 配向場との結合系の時空カオス が生じ ることが知られている 50]。液晶の電気対流は理論的には Kramer 等によ り詳し く研究されている 51]。この液晶対流系に対しても 、流れ場と振幅方程式 の複雑な結合方程式が提案された 52]。われわれも、そのモデルの簡略モデルを 数値計算し 、ロール構造の横方向へのジグザグ振動がおこり、それがカオス化し て欠陥カオスが生じ ることを見いだした 53]。ジグザグ振動生成の簡単なモデル 方程式も出した。定量的な解析では流れ場がこの時空カオス化の原因とも言い切 れないようであるが 、実験の様子からは横方向の振動的なグライド 運動が時空カ オスの基礎的な運動モード のように見える。図 18(a) はジグザグ振動を表してお り、図 18(b) は時空カオス状態をあらわす。 通常の液晶の電気対流では 、液晶分子を上下の壁に水平に配置させる。これ を垂直配向させたものをホメオトロピック配向という。ホメオトロピック系では、 電圧を上げていくとまず、液晶の分子配向が一部傾く。傾く向きは任意である。さ らに 、電圧を上げると、電気対流が生じる。この時、分子配向の向きとロール形 成のダ イナミクスがカップルする。傾きの方向が任意なため、向きを変える自由 度がやはり遅い変数になる。このため、対流が起こるとほぼ同時にゆっくりと変 動する時空カオスが観測される。甲斐先生たちはこれを Soft mode turbulence と 呼んでいる 54]。遅いモード が結合しあうとカオス化しやすい例だと考えられる。 われわれは分子配向の向きの自由度を複素変数のギンツブルグ−ランダウ方程式 で表し 、対流モード を Swift-Hohenberg 方程式でモデル化し数値計算した 55]。 @W = W ; jW j2W + r2W ; @u + i @u @u W + @u W @t @x @y @x x @y y @u = u ; u3 ; (1 + r2)2 u ; W 2 @ 2 u + 2W W @ 2 u + W 2 @ 2 u x y @x@y x @x2 y @y2 @t (45) 結合がなければ 、W の向きも、ロールの向きの任意の向きをとれ る。結合係数 がともに正なら 、ロールの向きは W と垂直方向を取る傾向があり、逆に、W の向きはロールに平行方向になる傾向がある。これでは互いに好ましい方向をと れない。ロールが x 軸を向いていれば 、W は x 方向になる。W は x 方向になれ 25 150 150 100 100 50 50 0 0 50 100 0 150 0 150 150 100 100 50 50 50 100 150 100 150 0 0 50 100 50 150 Figure 20: スパイラルコアの位置を操作して対消滅する過程 ば 、ロールは y 方向にずれていく。この向きに関するある種のフラストレーショ ンの結果時空カオスが生じ る。図 19 はこの数値計算であらわれた時空カオス状 態である。様々な向きのド メインが生じ 、ゆっくりと変動している。 ８． スパイラルパターン スパイラルパターンすなわち渦巻きパターンは静止しているのに回転や拡大 縮小など の動きを感じさせる図形である。１回転して、半径方向にスケール変換 すると元に戻るという対称性に起因しているものと思われる。渦巻きに吸い込ま れることを想像するからかもしれない。非平衡系のパターン形成の中では、BZ 反応のスパイラルパターンがよく知られている 56]。ある種の粘菌もスパイラル パターンを形成する。一様振動解に関する振幅方程式である複素ギンツブルグ -ラ ンダウ方程式にもスパイラル解があらわれる。中心点の周りに１周すると位相が 2 まわる位相欠陥があると、中心から位相波が伝播し , スパイラルパターンがで きる。結晶成長でも 、らせん転位の周りにスパイラル状に結晶成長が進むことが あり、このらせん成長は結晶成長の重要な成長様式の１つと考えられている。 BZ 反応では多くのスパイラルが生成されるが、実験ではほとんど そのコアの 位置は動かないように見える。我々は複素ギンツブルグーランダウ方程式の位相 欠陥を２つ用意し 、その運動を調べた。複素係数が実数の場合、スパイラルを形 成せず、右回りと左回りの欠陥対は互いに引き合うが、同方向の欠陥対は反発し あう。一般の複素係数の場合は 、逆向きの欠陥対の距離が近いと引き合い消滅す るが 、ある度程度引き離すと反発するようになることがわかった 57]。通常の BZ 反応ではスパイラルパターンとターゲットパターンが共存していると, スパイラ ルパターンの振動のほうが速いので、次第にスパイラル領域が増大していくこと 26 Figure 21: スパイラルカオスの除去過程 が知られている。複素ギンツブルグ ーランダウ系でも、自然振動数が局所的に高 い領域を作ってやるとそこを中心にターゲットパターンができる。自然振動数の 十分大きくすると, スパイラルパタ-ンの振動数より大きくすることができ、逆に ターゲットパターンの領域が増大する。ターゲット領域がスパイラルのコアに近 づくとスパイラルコアが動き出し 、最終的には計算領域から出て行ってしまう。 反応拡散系のモデル方程式として FitzHugh-Nagumo モデルがある。神経軸 策中のパルス伝播を説明するモデル方程式とし て提案された方程式である。1 次 元では興奮パルスが一定の速度で伝播する。２次元系ではパルス伝播がスパイラ ル状になることができる。われわれはこのモデル方程式 @u = u ; u3=3 + r2 u ; v + a(x y) @t @v = (u ; v + ) @t (46) を用いて、スパイラルコアをピンセットでつまんで移動させるようなことができ ることをシ ミュレーションで示した 58]。u の式の中にある a の大きさを変える と, 興奮性を変えることができる。スパイラルコアの付近で a を局所的（スポット 状）に低くすると , スパイラルコアはその領域にトラップされるようになる。その 状態で 、スポットを移動させると、スパイラルコアはスポットに引きずられて移 動する。図 20 は２つの逆回転のスパイラルを最初離しておく。上の方法で２つの スパイラルコアを移動させて最終的に対消滅させた。スパイラルの位置をコント ロールできるということである。光感受性をもつ BZ 反応ではスポット 光を使っ てスパイラルを位置をコントロールできるのではないだろうかと考えている。 ペースメーカー細胞以外の多くの心筋細胞も興奮性媒質と考えられる。神経細 胞の興奮と異なる特徴として、1 回の興奮の持続時間が興奮の頻度（ 心拍数 ）に よってかなり変化することがあげられる。この点をモデル化した方程式に AlievPanlov モデルがある 59]。このモデルでは 、スパイラルが不安定化し , たくさ んのスパイラルコアが自発的に生成される時空カオスが生じ ることが知られてい る。われわれもこの系の数値計算し 、スパイラルの不安定化がパルスの脈動不安 定性（パルスの幅すなわち興奮の持続時間が振動する不安定性）によることを示 27 2 3 1 2 0 -1 1 -2 0 6 8 10 x 12 14 0 5 10 x 15 20 Figure 22: (a) パラメータ振動系の局在パターン (b) 空間周期パターン (c) ２ 次元局在パターン した。さらに 、この脈動の周期付近の振動数で強制振動外力を加えると, この時 空カオス状態を除去できることを示した 60]。図 21 は初期にスパイラルカオス状 態にあった系が周期外力で静止状態になった例を示している。心臓における心拍 の異常は不整脈と呼ばれるが 、このスパイラルカオスのような空間的にも不規則 で、時間的にも不規則で速い心拍異常は心室細動と呼ばれ 、緊急処置を必要とす る重篤な不整脈である。人の心室細動とスパイラルカオスとの関係は必ずし も明 らかではないが 、時空カオスの制御が現実に問題になる例である。 ９．局在パターンとソリトン a. 局在散逸構造 空間的に局在したパターンが散逸構造とし てあらわれ ることがある。水とア ルコールのような２種類の流体の混合流体の熱対流では、進行波動が散逸構造と してあらわれることがある。定常状態から進行波状態の中間に局在進行波パター ンが実験的に観測されている 61]。また、粉体を鉛直方向にゆするとさまざ まな 空間周期パターンが得られるが 、やはり定常状態との中間に局在振動状態が得ら れている 62]。 複素ギンツブルグ ーランダウ方程式の係数が 、純虚数になるとエネルギーを 保存する非線形シュレデ ィンガー方程式になる。非線形シュレデ ィンガー方程式 はソリトン解という空間的に局在したパルス解をもつ。ソリトンは保存系のパル ス解であり、パルスの振幅、速度が任意の値を取ることができる。複素ギンツブ ルグ -ランダウ系の係数が虚数に近いとき（ 実部が虚部に比べて小さいとき ）ソリ トン的な解がアトラクタとし て得られることが期待される。Fauve たちは５次の 複素ギンツブルグ -ランダウ系で数値計算により、アトラクタとしての一定振幅の パルス解を見つけた 63]。われわれは同様に５次の複素 Swift-Hohenberg 系で一 定速度で伝播するパルス解および定在パルス解を見出した 64]。２成分流体の対 流の局在進行波パターンに対応する解といえるかもしれない。 また、３次の複素ギンツブルグーランダウ方程式にパラメータ振動項（ W の 28 複素共役 W ）を付け加えた式 @W = (1 + i!)W ; (1 + ic )jW j2W + (1 + ic ) @ 2 W + W (47) 2 1 @x2 @t を数値計算して、空間周期解および 局在振動解を得た 65]。特に 、c1 = c2 = 0 の 場合は sech 型の厳密解が得られている。パラメータ振動系の局在解は粉体振動層 の局在パターンに関連しているように思われる。図 22 にこの系の空間周期解と 局在解、および 、その２次元系の局在パターンを示す。 液晶の電気対流で電圧を大きくしていくと 、グリッド パターンと呼ばれる格 子状のパターンがあられる。このパターンはあるパラメータ領域で振動する。そ のとき, 局在ターゲットパターンが得られることを佐野先生等が実験で見出した 66]。空間周期構造が振動不安定化するときの振幅方程式が Coullet 等によって導 出されている 67]。その振幅方程式 @W = (1 + ic )W ; (1 + ic ) j W j2 W + (1 + ic ) @ 2 W + ( + i ) @ W 0 2 1 @x2 1 2 @x @t @ = @ j W j2 +i @W W ; W @W + @ 2 (48) @t @x @x @x @x2 を数値計算し 、局在ターゲットパターンに対応するパターンを得た。このパター ンに対応する sech 型の厳密解も得た 68]。 これらの、局在散逸構造は拡散項のため、進行速度も一定値に近づく。すなわ ち速度が０（動かない）か 、一定の進行速度を持っている。われわれは、2 つの 接近したコア部分をもつ光ファイバー中の散逸ソリトンを想定して、空間微分の 係数の実部が０になる結合複素ギンツブルグ -ランダウ方程式を数値計算し , 振幅 はある一定値になるがその進行速度が任意の値を取る局在パルス解を得た 69]。 b. ソリトン 非線形シュレデ ィンガー方程式はソリトン方程式の代表例とし て詳し く研究 されてきた。特に 、光ファイバー中のソリトンは非線形シュレデ ィンガー方程式 で記述でき、光ソリトンを使った高密度光通信の可能性が詳し く研究されてきて いる。近年のボ ーズーアインシュタイン凝縮体の研究の高まりの中 、この BEC 凝縮体でもソリトンが実験的に見出された 70]。BEC 凝縮体は磁場やレ ーザー光 で閉じこめられているので 、この物質波の従う方程式はポテンシャル中の非線形 シュレデ ィンガー方程式になる。 1 @ 2 ; jj2 + U (x) i @ = @t 2 @x2 (49) ただし 、簡単のためプランク定数、質量を１にし 、有効原子間力をあらわす非線 形項の係数も１にしてある。われわれは、ある狭い幅の領域だけで U0 > 0 を持つ ポテンシャルの壁でのソリトンの散乱の数値計算を行った 71]。ポテンシャルの 壁の高さ U0 が大きいときは 、非線形項の効果が小さく、量子力学的な振る舞い が見られるが 、U0 を小さくすると, ソリトンのまま、分裂することなく透過ある いは反射することを見出した。さらに透過と反射の境となる臨界運動エネルギー 29 (a) (b) 150 150 t t 100 100 50 50 0 0 0 10 20 x 30 40 50 0 10 20 x 30 40 50 Figure 23: ソリトンのポテンシャル散乱 (a) 反射 (b) 透過（トンネル現象） が壁の高さ U0 より小さくなることを見出した。運動エネルギ ーが U0 以下なの に、ほぼ 全部の粒子が壁を通り抜けるトンネル効果を表している。たとえていう と, 通ってはいけない赤信号も手をつないでみんなで渡ると心理的な壁が低く なるようなものである。図 23 にソリトンが反射と透過している様子を示す。 周期的ポテンシャル中のソリトンの研究は近年盛んに行われている 72]。周期 的ポテンシャル中のシュレデ ィンガー方程式の基本解はブロッホ態になる。非線 形性がきくと、ブロッホ状態のソリトンが生成される。特に、バンドギャップ付 近ではブロッホ状態の有効質量が負になる。有効質量が負のソリトンをさらに空 間的にゆっくり変化する調和ポテンシャルの中で運動させた 73]。数値計算の結 果, 上に凸な反閉じ 込めポテンシャル中でソリトンは調和振動するが 、下に凸な 本来粒子を閉じ 込めるはずの調和ポテンシャル中では、そこから加速度的に逃げ ることがわかった。質量が負の粒子なのでこの振る舞いはもっともらし いもので あるが 、質量が負の粒子が現実に BEC 系など で観測されれば面白いと思う。 １０． 成長パターン 非線形非平衡系の物理とし て、樹枝状結晶や DLA（拡散律速凝集体）の研究 が盛んに行われてきた。熱平衡条件から離れた、過冷却された液体が結晶化する ときの形態形成の問題は、数学的には動く界面に対するステファン問題の一例で あり、難しい問題とされてきた。針状結晶の成長速度とその先端曲率の関係は 、 1980 年代に放物界面の先端曲率の選択問題とし て解かれた 74]。先端曲率半径 と成長速度 v の間に v2 が一定値になるという関係が過冷却度によらず成り立つ というのが最大の成果である。DLA はやはり 1980 年代にフラクタルになる物理 系とし て詳し く研究された 75]。拡散律速凝集系では 、そのフラクタル次元が系 の種類によらずある一定値（ 2 次元では 1.7 程度）になる。実際に、電解析出で 30 450 400 350 300 250 200 150 100 50 0 0 100 200 300 400 500 600 Figure 24: (a) 結合写像モデルの６回対称樹枝状形態 (b) 線上で成長したデン ド ライト 300 250 200 150 100 50 0 0 100 200 300 400 500 600 700 Figure 25: フェイズフィールド モデルによる先端振動する樹枝状パターン できる金属葉、オイルに空気などを押し込んでできるヴィスコスフィンガリング、 バクテリアコロニーなど さまざ まな系で DLA 状のパターンが得られている。 われわれは 、結合写像系や phase eld model と呼ばれている偏微分方程式系 でこれらの成長パターンを作った。結合写像格子モデルは界面サイトで xn+1 (i) = xn(i) + c1 (Tc ; un(i)) un+1 (i) = un(i) + c2 (Tc ; un(i)) (50) の形の写像で時間発展する 76]。i は格子点の位置を表し 、u が温度場、x が結晶 化度を表すオーダーパラメータである。c2 に比例に比例する項は結晶化に際して の潜熱の発生を表している。この写像と熱拡散方程式の離散系を交互に時間発展 させ、成長パターンをつくった。系の異方性は写像を置く格子の種類で決まる。 このモデルは写像系なので 、パソコンなど で簡単に短時間でさまざ まなパターン （ 正方形のコンパクトパターン 、DLA 状のパターン 、DBM のような密に分岐し たパターン 、樹枝状パターンなど ）を作ることができる。図 24 にこのモデルで 作った、雪のような 6 回対称パターン 、および線上から成長させた競合するデン ド ライパターンを示す。 上の結合写像系は簡単にシミュレションできるが、パラメータの意味や定量性 に欠ける面がある。phase eld model はそのパラメータと界面張力など の本来の 31 Figure 26: (a) 競合進化のモデル。性質 z の分布が z が大きい方へ進化 (b) 進 化速度と競合度 (a1 ) の関係 の結晶成長モデルのパラメータとの対応がはっきりしている。特に、Karma 等に よって改良された phase eld model を使ってさまざまな計算がなされている 77]。 われわれも, このモデルを使って、樹枝状パターンの先端振動や針状形態が２つ結 合して安定化するダブロンパターンの研究を行った。樹枝状形態の横枝の生成機 構は ノイズが増幅されたものであるという説が有力であるが 、本庄先生たちは、 塩化アンモニウムの樹枝状結晶の実験で成長方向が (1 0 0) から (1 1 0) に変化 する近傍で、かなり規則的な振動成長を見出した 79]。この振動の結果、横枝も かなり規則的にでる。条件によっては定常成長が不安定化し 、振動成長がおこる 可能性があるのかもしれない。この成長方向が変化する現象はフェイズフィール ド モデルでも再現できる。図 25 はその成長方向が変わる付近で見られた規則的 先端振動パターンである 78]。先端付近の温度場がリミットサイクル的な振動を 示していることも確認している。また溝があるダブロン構造の先端振動形態も見 いだしている １１．進化とネットワーク形成のモデル 非線形力学の研究は生物系や経済など のより複雑な系にも展開されてきてお り、一般に複雑系と総称されている。物理学会の統計力学の分野でも物理学とは いえないようなさまざ まなテーマで講演がなされている。さまざ まなアプローチ で複雑系の研究がされているようであるが 、われわれは非線形ダ イナミクスの不 安定現象の一例とし てあるト イモデルの計算をおこなってみた。 a. 進化のモデル 生物の進化は突然変異と生存競争の結果起こるというのがダーウィンの進化 論である 80]。進化のメカニズムの１つに「 赤の女王」進化と呼ばれる進化形態 があるという説がある。たとえていうと, 足の速い獲物は敵から容易に逃げられ 適応力が高いので生きのび る。逆に 、足の速い捕食者も獲物を取りやすいので適 応力が高く生存競争に勝つ。互い足の速いものが競争に勝つので、足の速さはど んどん進化していく。止まっていることができず常に動いていることが強いられ る。 「 鏡の国のアリス」の中の「 赤の女王」の言葉になぞらえてこう呼ばれる。こ 32 0.07 0.025 0.06 0.02 0.05 0.015 0.04 0.01 0.03 0.02 0.005 0.01 0 0 0 50 100 i 150 200 0 50 100 i 150 200 Figure 27: (a) 中間規模の階層構造 (b) 全域にわたる階層構造 のような定向進化が実際に起こったかど うかは議論があるようである。われわれ は、レプリケータモデルとよばれる力学モデルを使って計算機中ではこの過程が 起こりうることを示した 81]。レプリケータモデルは dxi = x (f ; c) + D X(x ; x ) (51) i i j i dt と書ける。xi は i 番目の種のポピュレーション、fi はその適応度、D は突然変異 で j から i に変わる割合を示す。c はポピュレーションの総和が一定になる制限 項でこれがあるので互いに競合する。D = 0 なら 、最大適応度を持ったものだけ が生きのびる。 「赤の女王」型の進化を考えるので 、競合する２種を考え、それぞ れの種内に足の速さ z のような量が異なるものが分布しているとする。すなわち @y1 (z ) = y (z )(f (z ) ; c) + D @ 2 y1 (z ) 1 1 @t @z 2 @y2 (z ) = y (z )(f (z ) ; c) + D @ 2 y2 (z ) (52) 2 2 @t @z 2 被食者 (獲物 ）の種 1 に関し ては捕食者 2 の z よりも大きい z を持つものが 適 f (z ) が大きくなるように f (z ) の形を仮定する。たとえば 、f1 (z ) = a1 + R応度 z y (z )dz 。捕食者の適応度も同様に仮定する。実際に計算して見ると、図 26(a) ;1 2 にあるように、z の分布はパルス状になり一定速度で z が大きいほ うに移動して いくことがわかった。この型の進化は 、捕食被食関係にある種間だけでなく、互 いに同等な競合関係にある種間でも、互いに相手を出し 抜こうとし て一定速度で 進化が起こる。このモデルでは互いの力関係が接近しているほど 競合が強く進化 速度が速くなることがわかった。図 26(b) にそのその力関係 (相互作用がないと きの適応度 a1 。ただし a2 = 2 とし ている）と進化速度の関係を示している。力 関係が同等な 2 種 a1 = a2 の時、最大の進化速度が得られている。互いに競争す れば進化するという競争至上主義のト イモデルとなっている。 b. 階層構造生成のモデル 33 最近、インターネット、送電線のネットワークからウエッブサイトの引用被引 用関係や俳優の共演関係まで様々なネットワークで、そのリンク数の分布をとる とベキ分布を示すことが多いという報告があり、スケールフリーネットワークと 呼ばれている 82]。ネットワーク構造にある種のフラクタル性がみられるといえ る。確率モデルを用いてスケールフリーネットワークの研究がされていることが 多いようである。 生態系には食う食われ るの関係（ 食物連鎖）にもとづいたピラミッド 状の階 層構造がある。森林には太陽光を求めての競合にもとづいた、高木、低木、下草 の階層構造がある。近代国家の代表的な階層構造は軍隊と官僚組織と言われる。 階層構造が生まれる要因は一般には複雑で多様だと考えられるが、われわれは系 のダ イナミカルな不安定性によって、一様状態から階層的構造が生成される簡単 な力学モデルを考えた 83]。一様状態が不安定化し 空間周期パターンが形成され るのと同様なメカニズムで階層構造が生成されるト イモデルをつくる。進化する 生態系を念頭に置いて次の形のレプリケータモデルの方程式を用いる。 dxi = x (X W x ; c) i ij j dt j (53) c は全ポピュレ ーションを一定値にする項である。Wij は非線形の相互作用を表 し 、ポピュレションの比によって変化すると仮定する。その比が１に近いときは 互いに競合し (Wij < 0 Wji < 0 ）、その比がある程度大きくなる (xi =xj > 1:5) と、食う食われるの関係 (xij > 0 xji < 0 ）になり、さらに大きくなると無関係 （ Wij = Wji = 0 ）になるようにする。さらに 、この相互作用の範囲はその平均 ポピュレーション (xi + xj )=2 に比例して大きくなるようにする。各要素を１次元 または２次元空間内に配置して計算機シミュレーションを行った。ど ちらの系で も一様状態は不安定化し 、xi が大きい優者と xi が小さい劣者があらわれ 、その 間に食う食われるの関係がローカルに生じ る。時間とともに 、ローカルな優者間 で競合が生じ 、より大きいスケールで優者間に優劣が生じ 、３層構造になる。こ の過程を繰り返し て、より大きな階層構造も生じ ることがわかった。図 27 は１ 次元系での階層構造形成の例を示している。横軸が i で縦軸が xi 、食う食われる 関係にあるもの間を線で結んでいる。図 27(a) は相互作用の距離のパラメータが 中間的でローカルな階層構造が形成された状態でほぼ定常になった例で、(b) は 遠くまで影響が及ぼす系で 、最上位のものから最下層のものまでの中央主権的な 階層構造ができている。 References 1] 2] 3] 4] 5] R.Lorenz, J. Atom. Sci. 20, 130 (1963). R. M. May, Nature 261, 459 (1976). H. Sakaguchi and K. Tomita, Prog. Theor. Phys. 78, 305 (!987). T. H. Solomon and J. P. Gollub, Phys. Rev. A 38, 6280 (1988) H. Sakaguchi, Phys. Rev. E 65, 067201 (2002). 34 6] 岩波講座 現代物理学の基礎５ 「 統計物理学」岩波書店. 7] H. Sakaguchi, J. Phys. Soc. Jpn. 70, 3247 (2001) 8] E. M. F. Curado and C. Tsallis: J. Phys. A 24, L69 (1991) . 9] R. P. Feynman, R. B. Leighton, M. Sands: Lectures in Physics (AddisonWesley, Reading, 1963) 10] M. O. Magnasco: Phys. Rev. Lett. 71(1993)1477. 11] H. Sakaguchi, J. Phys. Soc. Jpn. 67, 709 (1998). 12] 「 生体分子モーターの仕組み」 共立出版。 13] G. Gallovotti and E. G. Cohen, Phys. Rev. Lett. 74, 2694 (1995). 14] H. Sakaguchi, J. Phys. Soc. Phys. 69, 104 (2000) 15] L. Onsager, Phys. Rev. 65, 117 (1949). 16] H. Sakaguchi, Prog. Theor. Phys. 80, 7 (1988). 17] H. Sakaguchi, Phys. Rev. E 60, 7584 (1998). 18] 岩波講座 現代の物理学 「 高分子物理・相転移ダ イナミクス」岩波書店. 19] H. Sakaguchi, Phys. Rev. E 64, 047101 (2001). 20] I. M. Lifshitz and V. V. Slyozov, J. Phys. Chem. Solid 19, 35 (1961). 21] Y. N. Yang, E. S. Fu and E. D. Williams, Surf. Sci. 356, 101 (1996). 22] H. Sakaguchi and N. Fujimoto, Phys. Rev. E 68, 056103 (2003). 23] Y. Kuramoto, Chemical Oscillations, Waves, and Turbulence SpringerVerlag. 24] H. Sakaguchi and Y. Kuramoto, Prog. Theor. Phys. 76, 305 (1987). 25] K. Wiesenfeld, P. Colet and H. Strogatz, Phys. Rev. E 57, 1563 (1998). 26] H. Sakaguchi and K. Watanabe, J. Phys. Soc.Jpn, 69, 3545 (200). 27] H. Fujisaka and T. Yamada, Prog. Theor. Phys. 69, 32 (1983). 28] M. Rosenblum, A. Pikovsky and J. Kurths, Phys. Rev. Lett. 76, 1804 (1996). 29] H. Sakaguchi, Phys. Rev. E 61,7212 (2000). 30] H. Sakaguchi, S. Shinomoto and Y. Kuramoto, Prog. Theor. Phys. 77,1005 (1987). 31] H. Sakaguchi, Prog. Theor. Phys. 80, 743 (1988). 35 32] P. Bak, C. Tang and K. Wiesenfeld, Phys. Rev. Lett. 59, 381 (1987). 33] 「 地震発生の物理学」大中、松浦著 東京大学出版会. 34] J. M. Carlson and J. S. Langer, Phys. Rev. Lett. 62, 2632 (1989). 35] H. Sakaguchi, J. Phys. Soc. Jpn, 72, 69 (2003). 36] 「 散逸構造」ニコリス、プリゴジン著 岩波書店. 37] M. C. Cross and P. C. Hohenberg, Rev. Mod. Phys. 65,851 (1993). 38] 「 散逸構造とカオス 」森、蔵本著 岩波講座 39] H. Sakaguchi, Prog. Theor. Phys. 95, 823 (1996). 40] H. Sakaguchi, Prog. Theor. Phys. 84, 792 (1990). 41] H. Sakaguchi, Prog. Theor. Phys. 93, 491 (1995). 42] P. G. de Gennes and J. Prost, "The Physics of Liquid Crystals" Clarendon Press. 43] H. Sakaguchi, Physica D 129, 57 (1999). 44] V. Yakhot, Phys. Rev. A 24, 642 (1981). 45] M. Kardar, G. Parisi and Y. Zhang, Phys. Rev. Lett. 56, 889 (1986). 46] K. Sneppen, J. Krug, M. H. Jensen, C. Jayaprakash, and T. Bohr, Phys. Rev. A 46, R7351 (1992). 47] H. Sakaguchi, Phys. Rev. E 62, 8817 (2000). 48] H. Sakaguchi, Prog. Theor. Phys. 62, 8817 (2000). 49] A. Zippelius and E. D. Siggia, Phys. Fluid 26, 2905 (1983). 50] 「パターン形成」蔵本等著 朝倉書店 ４章. 51] "Pattern Formation in Liquid Crystals", A. Buka and L. Kramer SpringerVerlag ６章. 52] M. Kaiser and W. Pesch, Phys. Rev. E 48, 4510 (1993). 53] H. Sakaguchi, Prog. Theor. Phys. 92, 509 (1994). 54] S. Kai, K. Hayashi and Y. Hidaka, J. Phys. Chem. 100, 19007 (1996). 55] H. Sakaguchi, Prog. Theor. Phys. 78, 1046 (1997). 56] A. T. Winfree, "The Geometry of Biological Time" Springer-Verlag. 57] H. Sakaguchi, Prog. Theor. Phys. 82, 7 (1989). 36 58] H. Sakaguchi and T. Fujimoto, Prog. Theor. Phys. 108, 241 (2002). 59] B. R. Aliev and A. V. Panlov, Chaos, Solitons, Fractals, 7, 293 (1996). 60] H. Sakaguchi and T. Fujimoto, Phys. Rev. E 67, 067202 (2003). 61] P. Kolodner, Phys. Rev. Lett. 69, 2519 (1992). 62] P. B. Umbanhowar, F. Melo and H. L. Swinney, Nature 382, 793 (1996). 63] O. Thual and S. Fauve, J. Phys. (Paris) 49, 1829 (1988). 64] H. Sakaguchi and H. Brand, Physica D 117, 95 (1998). 65] H. Sakaguchi and B. Malomed, J.Phys. Soc. Jpn. 72, 1360 (2003). 66] S. Nasuno, M. Sano and Y. Sawada, J. Phys. Soc. Jpn. 58, 1875 (1989). 67] P. Coullet and G. Iooss, Phys. Rev. Lett. 64, 866 (1990). 68] H. Sakaguchi, Prog. Theor. Phys. 87, 1049 (1992). 69] H. Sakaguchi and B. Malomed, Physica D 183, 282 (2003). 70] L. Khaykovich et. al., Science 296, 1290 (2002). 71] H. Sakaguchi and M. Tamura, preprint. 72] B. B. Baizakov, V. V. Konotop and M. Salerno, J. Phys. B 35, 5105 (2002). 73] H. Sakaguchi and B. Malomed, preprint. 74] D. A. Kessler and H. Levine, Phys. Rev. B33, 3352 (1986). 75] T. A. Witten and L. M. Sander, Phys. Rev. Lett. 47, 1400 (1981). 76] H. Sakaguchi, J. Phys. Soc. Jpn. 67, 96 (1998). 77] A. Karma and W. J. Rappel, Phys. Rev. E 53 , R3017 (1996). 78] H. Sakaguchi and S. Tokunaga, Prog. Theor. Phys. 115, 1043 (2003). 79] H. Honjo and S. Ohta, J. Crys. Growth 131, 111 (1993). 80] 「 進化生物学」D.J. フツイマ 蒼樹書房. 81] H. Sakaguchi, Europhys. Lett. 37, 571 (1997). 82] R. Albert and A. Barabasi, Rev. Mod. Phys. 74, 48, 2002. 83] H. Sakaguchi, Phys. Lett. 313, 188 (2003). 37