カオス的アトラクタとリヤプーノフ指数

さて，実際に非周期的状態に遭遇した場合，それがカオスかどうかを判定するにはどうしたらよいの であろうか．以下この判定によく利用されるリヤプーノフ指数について述べる．リヤプーノフ指数はカ オス軌道の接空間での平均的な延び縮みを表す量である．これは軌道に沿った変分の特性指数から求め られる．

簡単のためポアンカレ写像：

T :Rⁿ→Rⁿ; x₀7→x₁=T(x₀) (8.54) で定義される離散力学系の場合を考える．元の方程式は式 (8.19)^{のような周期} L ^{の周期的非自律系と} しよう．自律系の場合は，たとえば軌道を一定時間間隔で離散化すれば同様の議論が可能となる．

いま，任意の初期状態 x₀∈Rⁿ を出る正の半軌道を次式としよう．

γ⁺(x0) =©

x∈Rⁿ|x(k) =T^k(x0), k= 0,1,2, . . .ª

(8.55) 時刻 k^{における状態}x(k) ^の初期値x0 に関する微分は

∂x(k)

∂x0 =DT^k(x0) = ∂

∂x0T^k(x0) = Yk

`=1

DT(x(`)) (8.56)

となる．ここで微分記号はこれまでのように D^あるいは ∂/∂xで表しておいた．いま，式(8.56)^の特 性乗数を考える：

©µ1(k), µ2(k), . . . , µn(k)ª

=©

µi(k)∈C|det(µi(k)In−DT^k(x0))ª

(8.57) そこで，t→ ∞ としたとき次の平均値や極限値があれば，これを軌道(8.55)のリアプーノフ数という．

mi(k) = p^k

|µi(k)|, mi= lim

k→∞mi(k) = lim

k→∞

pk

|µi(k)|, i= 1,2, . . . , n (8.58)

x y

D

Ch⁴₁ Ch⁴₃

Ch⁴₂ Ch⁴₄ I

I²₁

I²₂

(c)

x y

D

Ch²1

Ch²₂ I

(d)

x y

D

Ch

(e)

x y

D

Ch

(f) x

y D

I S²₁

S²₂

(a)

x y

D

S⁴₁

S⁴₄ S⁴₃

S⁴₂ I

I₁²

I²₂

(b)

図8.20 式 (8.53) の相平面図．B0 = −0.075 の場合．(a) B = 0.15，(b) B = 0.185，(c) B= 0.195，(d)B = 0.197，(e)B = 0.199，(f)B = 0.217．

x y

D

I

(g)

図8.21 式(8.53)の相平面図．B0=−0.075の場合．(g)B = 0.23．

また，軌道(8.55) ^{を元の周期}L の周期的非自律系の連続な軌道に返して考えると次の平均指数が定義 できる．

νi(k) = 1

kLlog|µi(k)|, νi= lim

k→∞νi(k) = lim

k→∞

1

kLlog|µi(k)|, i= 1,2, . . . , n (8.59) これをリアプーノフ指数(Lyapounov exponent) という．これらの量は点 x₀ が T の双曲型固定点や 周期点のとき，特性乗数や特性指数に対応している．軌道 (8.55) が不変閉曲線や周期倍分岐の結果得 られた安定なカオスの場合にはリアプーノフ指数は，接空間での解の平均的な伸び縮みを表すと考えら れる．表8.1^参照．

いま，式(8.59)のリアプーノフ指数を大きい順に並べて番号を付けたとしよう：

ν₁> ν₂>· · ·> ν_n (8.60) そこでまず，最大リアプーノフ指数ν₁ の計算をしてみよう．長さ 1の勝手なベクトルe(0)∈Rⁿ を1 つ選んで，このベクトルに線形写像 DT(x(`)), `= 1,2, . . . を作用させるとどのように伸びるかをみて

表8.1 ポアンカレ写像の不変集合とリアプーノフ指数．

ポアンカレ写像T ^{の不変集合 リアプーノフ指数 周期} L ^{の周期的非自律系} タイプ _kD, _kI の固定点 k 個の正の指数 k 次元不安定周期解

または周期解 νi>0, i= 1, . . . , k νi<0, i=k+ 1, . . . , n

安定な不変閉曲線(ICC) ν1= 0, 安定な準周期（２重周期）解 ν_i<0, i= 2, . . . , n

安定なカオス ν₁>0, 1次元的に伸び，n−1 次元 ν_i<0, i= 2, . . . , n 的に縮む性質のカオス

みよう．

f(0) = e(0), ke(0)k= 1

f(1) = DT(x₀)e(0), ke(1)k= ke(0)k kf(1)kf(1)

· · ·

f(k) = DT(T^k−1x0))e(k−1), ke(k)k= ke(0)k kf(k)kf(k)

(8.61)

とおこう．このとき e(1) = ke(0)k

kf(1)kf(1) = ke(0)k

kf(1)kDT(x0)e(0) f(2) = DT(T(x₀))e(1) = ke(0)k

kf(1)kDT(T(x₀))◦DT(x₀)e(0) = ke(0)k

kf(1)kDT²(x₀)e(0) e(2) = ke(0)k

kf(2)kf(2) = ke(0)k²

kf(1)kkf(2)kDT²(x₀)e(0)

(8.62)

となるから，一般に次の関係を得る．

f(k) = ke(0)k^k−1 Q_k−1

i=1 kf(i)kDT^k(x₀)e(0) e(k) = ke(0)k

kf(k)kf(k) = ke(0)k^k Q_k

i=1kf(i)kDT^k(x₀)e(0)

(8.63)

そこで

1 = ke(k)k

kf(k−1)k = ke(0)k kf(k)k

kDT^k(x0)e(0)k kDT^k−1(x₀)e(0)k

より

kf(k)k

ke(0)k = kDT^k(x₀)e(0)k

kDT^k−1(x0)e(0)k (8.64)

の関係を得る．これらの関係式より Xk

i=1

logkf(k)k ke(0)k =

Xk

i=1

kDTⁱ(x₀)e(0)k

kDTⁱ⁻¹(x0)e(0)k = log Yk

i=1

kDTⁱ(x₀)e(0)k

kDTⁱ⁻¹(x0)e(0)k = logkDT^k(x₀)e(0)k ke(0)k

となる．一方e_k は kが大きくなるにつれて次第に伸びが一番大きい方向に向くことから ν₁= lim

k→∞

1 kL

Xk

i=1

logkf(i)k

ke(0)k = lim

k→∞

1 kL

Xk

i=1

logkDT^k(x₀)e(0)k

ke(0)k (8.65)

を得る．以上のことから最大リアプーノフ指数は，式 (8.61)でベクトルの長さを正規化しながら求め，

式 (8.65)を用いて計算すればよいことが分かる．

次に，n個のリアプノフ指数を求めることを考えよう．このためにはn個の一次独立なベクトルに線 形写像 DT(x(`)), ` = 1,2, . . . を作用させ，適当に長さを正規化しながら n 次元体積要素の変化を見 ればよい．

そこで，まず長さが 1 の n 個の独立な初期ベクトルを e₁(0), e₂(0), . . . , e_n(0) とし，これらの DT(x₀) による像を

f₁(1) =DT(x₀)e₁(0), . . . , f_n(1) =DT(x₀)e_n(0) (8.66) とする．f<sub>`</sub>(1), `= 1, . . . , n が張る体積要素の体積を保ち，かつ各伸びの部分空間に次の繰り返しの初 期値となるベクトル e1(1), e2(1), . . . , en(1) を取りたい．これにはグラム・シュミットの正規直交化 の方法を用いるとよい．すなわちベクトル

g₁(1) = f₁(1), e₁(1) = g₁(1) kg1(1)k

g2(1) = f2(1)−(f2(1), e1(1))e1(1), e2(1) = g2(1) kg2(1)k

· · ·

g_n(1) = f_n(1)−

n−1X

i=1

(f_n(1), e_i(1))e₁(1), e_n(1) = g_n(1) kgn(1)k

(8.67)

を構成する．そこで式 (8.66)^と式 (8.67)^によって

ei(k), fi(k), gi(k); i= 1,2, . . . , n, k= 1,2, . . . を作る．この作り方から

¯¯

¯det[f₁(k)f₂(k) · · · f_n(k)]

¯¯

¯ =

¯¯

¯det[g₁(k)g₂(k) · · · g_n(k)]

¯¯

¯

= kg₁(k)kkg₂(k)k · · · kg_n(k)k

(8.68)

-0.02 0.02

-0.04 0

0.12 0.14 0.16 0.18 0.2 0.22

0.10

B

ν

1

P

P

⁸

P

⁴

P

²

E

カオス・アトラクタ

0

D

₀

D

²

0D⁴

0

D P

^∞

-0.06 0.04 0.06 0.08

図8.22 式(8.53)の基本調波周期解に関する固定点の分岐図，k= 0.1の場合．

となっていることが容易にわかる．これはリアプノフ数 (8.59)^{の積となっている．}

¯¯

¯det[f1(k)f2(k) · · · fn(k)]

¯¯

¯ = m1(k)·m2(k)· · · · ·mn(k)

= e^ν1(k)kLe^ν2(k)kL· · ·e^νn(k)kL =e^{(ν1(k)+ν2(k)+···+νn(k))kL} (8.69) したがって，式(8.68) ^と(8.69) ^{の最後の項より}

ν1(k) +ν2(k) +· · ·+νn(k) = 1 kL

Xn

i=1

logkgi(k)k (8.70)

を得る．ところで k^{が十分大きくなると} gi(k) は平均的な伸びの部分空間を張っているので νi(k) = 1

kLlogkgi(k)k, i= 1,2, . . . , n (8.71) を得る．これは

νi= lim

k→∞νi(k) = lim

k→∞

1

kLlogkgi(k)k, i= 1,2, . . . , n (8.72) の近似値となっている．式(8.71)は具体的にリアプーノフ指数を計算する手順を与えている．

【例 8.10】最大リアプーノフ数の計算例

例 8.9でみた周期倍分岐からカオス状態に至るパラメータの変化に対するアトラクタの最大リアプー ノフ指数を計算しておこう．図 8.20 ^は，式 (8.53) ^において b0 = −0.075 ^{と固定し，}B ^を 0.1 ^から 0.22 付近まで変えたとき，各B の値で得られた最大リアプーノフ指数を描いたものである．

まず，B =P ^{において安定固定点} 0D が周期倍分岐により逆不安定固定点に変化し，P < B < P² において安定2 周期点 ₀D² が存在する．B=P²でこの安定な 2周期点が周期倍分岐する．同様な周 期倍分岐列が続き，B > P^∞ で最大リアプーノフ数が正となるアトラクタすなわちカオス状態が観察 されることとなる．B =E ではカオス集合がサドルのω 枝に接触し，アトラクタの性格を失う．その 結果，系の状態は過渡状態を経て共振状態に対応する安定固定点に落ち着く．

なお，図中安定状態の区間において最大リアプーノフ指数が一定となっている区間がみられる．この 区間では固定点，周期点の特性乗数が複素数となっている．この定数値は，リュービルの公式から

µ1µ2=|µ|²e^−2πk = (e^2πν)² したがって

ν1=ν2= −2πk 4π =−k

2 =−0.05 となっている． ■

ドキュメント内 Phase space Attractor mutual interaction selforganization jump hysteresis stabilize or destabilize Synchronization Self-excited Oscillation Nonlinear (ページ 190-196)