2009年2月28日土曜日

カオス #6 畑政義写像

畑政義写像再び.有名な図形を畑政義写像でも描ける例が見つかったので,今回は,それについて.

C曲線(Levy C curve)
ググれば一般的なC曲線が分かると思いますが,'C'の形をしています.この曲線は畑政義写像だと以下のパラメータで再現できます.

a=0
b=-0.5-0.5i
c=0
d=-0.5+0.5

このパラメータでプロットしていくと,Fig.1のようになります.


Fig.1 Levy C curve

コッホ曲線(Koch curve)
あの,雪印の一部みたいな曲線です.下のパラメータで実現できました(Fig.2).

a=0
b=0.5+0.25i
c=0
d=0.5-0.25i


Fig.2 Koch curve

2009年2月25日水曜日

カオス #5 バーニングシップ・フラクタル

今回は,名前がかっこいいバーニングシップ・フラクタルです.

バーニングシップ・フラクタル(Burning Ship fractal)

バーニングシップ・フラクタルは次の漸化式
zn+1=(|Re{zn}|+i|Im{zn}|)2+C
z0=0
において.lim(n→∞)znが発散しなかった時の複素数C全体の集合.

ほとんどマンデルブロ集合と同じです.実部と虚部のそれぞれで絶対値を取っている処だけが異なっています.

発散しない処を黒,発散する処をその速度に応じて色付けしたものが,次の図です.


Fig.1 Burning Ship fractal

マンデルブロ集合とは全く違います.何と無く燃えている船のように見えます.
中央左の小さなわだかまりに注目してください.ここを,ズームすると,下の図になります.


Fig.2 Burning Ship fractal(zoom)

ここでも,燃えている船が見えます.

この付近の発散しているところに注目すると,面白い幾何学模様が見えます.


Fig.3 Burning Ship fractal(zoom)

異なる場所を拡大するとFig.4のようになります.


Fig.4 Burning Ship fractal(zoom)

建築物のようにも見えます.歪んだタワー,あるいは鉄橋とも言えます.仄かに,建物の影さえ見えるような気がします.単純な式なのに,これ程異様で人工的なものを思わせる図形が出来るのが非常に不思議です.

参考文献

Wikipediaの執筆者たち.“バーニングシップ・フラクタル”.<http://ja.wikipedia.org/wiki/%E3%83%90%E3%83%BC%E3%83%8B%E3%83%B3%E3%82%B0%E3%82%B7%E3%83%83%E3%83%97%E3%83%BB%E3%83%95%E3%83%A9%E3%82%AF%E3%82%BF%E3%83%AB>.Wikipedia.(参照2009年2月22日)

2009年2月22日日曜日

カオス #4 マンデルブロ集合

きっと,誰もが見たことある,そう信じているのが下図です.


Fig.1 Mandelbrot set

この図はマンデルブロ集合を可視化したものです.

マンデルブロ集合とは,ブノワ・マンデルブロ(Benoît B. Mandelbrot)が発見したものです.彼は,フラクタルの父として著名です.

定義:
マンデルブロ集合とは,次の漸化式

zn+1=zn2+C
z0=0

において,lim(n→∞)znが発散しなかった時の複素数C全体の集合.

というわけで,たった数行で説明できるものです.複素数を知っている人なら,理解できます.Fig.1で言うと,黒い部分がマンデルブロ集合です.カラフルな所はマンデルブロ集合ではなく,znの発散速度に応じて,色づけしたものです.

ちなみに,実際の計算時には
zn=xn+yni
C=a+bi
=>
xn+1=xn2-yn2+a
yn+1=2xnyn+b
としています.

z2ときたら,z3がどうなるか知りたくなります.というわけで,指数を3,4とした時のものが下になります.


Fig.2 zn+1=zn3+C


Fig.3 zn+1=zn4+C

指数が2の時には一つだったダルマ状のモノが段々増えていきます.nを指数とすると,n-1個のダルマが出来るようです.

ついでに指数が1の時も…


Fig.4 zn+1=zn+C

これらの図形は自己相似っぽいフラクタルです.厳密には,拡大すると少しずつ違った図形となるらしいです.

ちょっと配色を変えたり,拡大したものを下にチョコチョコと載せておきます.

ギャラリー


Fig.1の一部を何倍かに拡大したもの.同じような形が出てきます.絵としてはこちらのほうが綺麗です.


Fig.1の一部を何倍かに拡大したもの.


割と綺麗に描画出来てると思います.


配色を変更.


こんな風に見える場所もあります.

次は,さらに指数が小数だとどうなるかを紹介しようかと思います.

参考文献

Wikipediaの執筆者たち.“マンデルブロ集合”.<http://ja.wikipedia.org/wiki/%E3%83%9E%E3%83%B3%E3%83%87%E3%83%AB%E3%83%96%E3%83%AD%E9%9B%86%E5%90%88>.Wikipedia.(参照2009年2月22日)

森川浩.“マンデルブロー集合”.<http://www2.neweb.ne.jp/wc/morikawa/man.html>.数理科学美術館.(参照2009年2月22日)

2009年2月19日木曜日

カオス #3 畑政義写像

今回は,畑政義写像です.最近知りました.

畑政義写像
この写像は,以下の写像で定義されています.
f1(z)=az+bz~
f2(z)=c(z-1)+d(z~-1)+1
(z=x+yi => z~=x-yi)

この写像に対して適当な初期値z0をとり,
zn+1=f1(zn), f2(zn)
を計算していきます.

f1,f2はぱっと見は割と簡単な写像ですが,組み合わせることにより,非常に複雑になるようです.

実際に計算し,z=x+yiをプロットしたものを以下に示します.単にプロットしただけでは,あまり面白みがないので,密度による色付けをしました.


Fig.1 a=0.7+0.2i,b=0.0,c=0.0,d=2/3,z0=1.0

葉っぱみたいです.下の参考文献を参考にしたのですが,予想以上に綺麗なものが出来ました.


Fig.2 a=0.0,b=0.45+0.5i,c=0.0,d=0.45-0.5i,z0=1.0

面白い幾何学模様です.この図形の次元は約1.736だそうです(参考文献参照).

パラメータが非常に多いので,色々模索する価値がありそうです.

参考文献
芹沢 浩.“2.畑政義の写像によるフラクタル図形”.<http://www001.upp.so-net.ne.jp/seri-cf/gallery/f09.html内>.カオス&フラクタル紀行.(参照2009年2月19日)

2009年2月6日金曜日

カオス #2 ロジスティック写像

こんばんは,今回もまた,カオス.あの有名な写像について.

ロジスティック写像(logistic map)
ロジスティック写像は,

Xn+1=a*Xn*(1-Xn) 0<=a<=4, 0<=X0<=1

で表される写像です.
これまた簡単な2次式です.

もともと,ロジスティック方程式という微分方程式がありました.

dN/dt=r(K-N)N
ロジスティック方程式

N - 個体数
r - 相対内的増加率
K - 環境収容力

この方程式は,個体群生態学において,個体群成長のモデルとして考案された微分方程式で,次のような性質を持ちます.

・個体数N=0では,増加率は0.
・個体数が増加するにつれて,増加率は減少する.
・環境の収納可能個体数には限界があるので,N=Kの時,増加率は0.

で,この微分方程式を解くと,

N=K/(1+exp(rK(t0-t)))

となります(t0は初期値).ちなみの上の解は,N→0(t→-∞),N→K(t→∞)となります.

というわけで,微分方程式の時はわけもなく解くことができました.

しかし,これを差分方程式にすると,途端に難しくなります.
ロジスティック写像の振る舞いは,aによって,次のようになります.

0<=a<=1
Xnは0へ単調収束.
1<a<=2
Xnは(a-1)/aへ単調収束.
2<a<=3
Xnは(a-1)/aへ振動しながら収束.
3<a<=3.5699456
Xnは2のべき乗個の周期点を振動.
3.5699456<a<=4
Xnは不規則で,特定の周期をもたない(カオス領域).
ただし,Xnが周期的になる領域が存在(周期性の窓).

文章や式だけでは,さっぱり分からんので,可視化しましょう.

下に示す分岐図は,縦軸にXn,横軸にaを取った時,(a, Xn)をプロットしたものです.


Fig.1 ロジスティック写像の分岐図(0<=a<=4,0<=Xn<=1)

aが3.5699456(ファイゲンバウム点といいます)を超えたあたりから,カオスな事になっていきます.
このカオス領域では,初期値X0がほんのわずか異なっただけで,将来の値Xnは大きく変わってしまいます.この写像にこれ以上深入りするのは止めにしておきましょう.

ところで,Fig.1はあまり綺麗ではありません.画像が全体的に滑らかでないです.ちょこっと調べた結果,分布の密度を画像のピクセル値に反映すればよいようです.

実際にやってみたところ,下のような図が出来ました.


Fig.2 ロジスティック写像の分岐図(分布密度版,0<=a<=4,0<=Xn<=1)

さらに右端(3.5<=a<=4)を拡大してみるとFig.3,Fig.4のようになります.


Fig.3 ロジスティック写像の分岐図(分布密度版,3.5<=a<=4,0<=Xn<=1)


Fig.4 ロジスティック写像の分岐図(分布密度版,3.5<=a<=4,0<=Xn<=1)

Fig.3,Fig.4は配色を変えてあります.
というわけで,単純な2次式から,こんな美しい画像が生成できました.

参考文献

Wikipediaの執筆者たち.“カオス理論”.Wikipedia.<http://ja.wikipedia.org/wiki/%E3%82%AB%E3%82%AA%E3%82%B9%E7%90%86%E8%AB%96>.(参照2009年2月6日)

森川浩.“ロジスティック写像”.<http://www2.neweb.ne.jp/wc/morikawa/log.html>.数理科学美術館.(参照2009年2月6日)

2009年1月23日金曜日

浮動小数点

ちょっと,浮動小数点の演算でハマっていた所があったのでメモがてら.

とある事情で,文字認識実験をしていました.その実験は,平仮名46文字を認識するというもので,認識手法の中に重み付きユークリッド距離というものがありました.重み付きユークリッド距離とは,

D(x, u) = √(∑(x[i]-u[i])^2/v[i])

x : 未知の文字の特徴ベクトル
u : ある字種の平均特徴ベクトル
v : uの分散ベクトル

というもので,適当にJavaで組んであるものを動かしていました.ところが,ある文字に対しては,非常に良い認識率が得られるが,ある文字は全く認識できない(認識率0%),という事態になってしまいました.
ファイルの読み込みや,計算式を入念にチェックしましたが,原因が見つかりませんでした.

よくよく考え,D(x, u)を出力してみました.

結果
NaN

ここでやっと原因が分かりました.単に0.0f/0.0fという演算が発生していたというだけのことです(上の式の場合x[i]-u[i]=0.0fかつ,v[i]=0.0fとなっていた).

そんなわけで,この問題は解決しましたが,浮動小数点における特殊な数について少しばかり調べてみました.上の実験では,Java,float型を使用してましたので,それに限りますが.
ちなみに,JDKのバージョンは1.6.0_06です.

Float.POSITIVE_INFINITY
正の無限大です.Floatクラスでの定義は,
public static final float POSITIVE_INFINITY = 1.0f / 0.0f;
でした.

Float.NEGATIVE_INFINITY
正の無限大です.Floatクラスでの定義は,
public static final float NEGATIVE_INFINITY = -1.0f / 0.0f;
でした.

Float.NaN
Not a Number(=非数)です.Floatクラスでの定義は,
public static final float NaN = 0.0f / 0.0f;

NaNの発生条件には以下のようなものがあります.

異符号の無限大同士の加算
同符号の無限大同士の減算
0と無限大の乗算
0同士の除算
無限大同士の除算
被除数が無限大の剰余
除数が0の剰余
負数の平方根
値が未定義となる逆三角関数,対数関数

また,この定数は特殊であり,
NaN op something (op -> "==", "<", ">", "<=", ">=")
は必ずfalseであり,
NaN != something
は必ずtrueとなります.

実際,あるfloat値がNaNかどうかを返すメソッドboolean isNaN(float)は以下のように定義されていました.

static public boolean isNaN(float v) {
return (v != v);
}

vがNaNの時にだけ,v != vはtrueとなります.

今日はこんなところで.

2009年1月17日土曜日

カオス #1 エノン写像

近頃,カオス理論に興味を持っています.理論は非常に難しいですが,それらを可視化した画像が非常に綺麗ですので,それらをちょこちょこっと紹介しようと思います.

エノン写像(Hénon map)
エノン写像とは,フランスの天文学者Hénonが考案したモデルで,次の式で定義される写像です.

x(t+1)=1-ax(t)^2+by(t)
y(t+1)=x(t)

手計算でも,コンピュータでも計算できる写像です.
で,下の図がa=1.4,b=0.3とした時のアトラクターです.


Fig.1 Hénon attractor(a=1.4, b=0.3)

アトラクター(attractor)とは,「力学系が十分長い時間を経た後に発展する集合」です.雑に言うと,ある初期値x(0),y(0)から,x(t),y(t)をどんどん計算していき,その時の(x(t), y(t))をプロットしたものです(実際には初期値から1000回程度計算させ,それ以降の(x, y)をプロットしています).写像を適用していった時の軌道みたいなものです.

このエノン写像.特徴は以下のとおりです(勿論,ここに挙げたものはほんの一部です).

・初期値x(0),y(0)を変化させても,アトラクターはほとんど変化しない.
・パラメータa,bにより,アトラクターは大きく変化する.
・自己相似性(フラクタル)をもつ.

上の図では,太い曲線があるように見えますが,実際は無数の曲線で構成されています.理論的にはどれだけ拡大しても,無限に曲線が見えますが,計算時間が有限ですので,途中で消えてしまいます.

最後に,少しずつパラメータを変えていくことにより,生成されるアトラクタのアニメーションを.



今回はこんな処で.

2009年1月13日火曜日

電子工作 #3 矩形波発生器

こんばんは,ここ最近色々弄っていたものについて.

まず,この頃,適当な部品を組み合わせ,電子楽器もどきを作ろうと思っていました.
音を鳴らすためには,適当な発振器で信号を生成して,それをスピーカに通せばそれだけで音が鳴ります.
で,音階も制御できないとマズイです,面白くないです.音階は発振器の周波数に対応するので,てっとり早く発振器の周波数を変えられる必要があります.
という事で,前回用いた発振回路と同じものを用いれば良いという結論に至りました(RまたはCを変えるだけで周波数を変えられます).
問題は,音の数ですが,3オクターブ位出せればいいかと思っていました.が,仮に一つの音に一つの発振回路を割り当てた場合,かなり費用がかかります.また,実際の回路制作にも負担がかかります.散々色々考えた挙句,下のような回路にしました.



可変抵抗+スイッチの部分は,実際にはもっと沢山(12個位)あります.また,トランジスタ付近は,増幅回路ですので無くても動きます.
ポイントは発振回路を一つにして,回路中の抵抗値と,静電容量(コンデンサ)をスイッチで変更できるようにしたことです.T=3RCから,f=1/3RC,つまり,静電容量が1/2倍になると周波数が2倍になり,1オクターブ上昇するということです.ですので,コンデンサ部分でオクターブを制御して,抵抗部分でドレミファ…を制御するということになります.

制作ですが今回は,ブレッドボードに組んだだけです(実験的に設計しただけですので).

この回路を作成した後,大変だったのは,音調整です.可変抵抗で地道に音合わせをする必要があるので,中々面倒でした.
実際に音を鳴らしてみましたが,矩形波による「ビー」という音が割と明瞭に聞こえます.増幅回路を通すとかなり大きな音が出ます.

又機会があれば,実際に基盤に組むか,その他の電子楽器もどきでも作ってみようかと思います.