今回はgnuplotで再帰定義の方法とその応用として漸化式やフラクタルのプロットについて紹介する。
1. 再帰定義
再帰定義を行うには、関数を定義する際にその関数自体を含むようにする。例えば階乗N!はf(n)=n*f(n-1)と書けるので以下のように再帰定義することができる。 ただし、n=0のときに値が必要になるので、前回紹介した3項演算子を用いてn=0に対して場合分けを行う。
fac(n)=(n==0)?1:n*fac(n-1) print fac(1) >1 print fac(4) >24
プロットする場合は、xが整数以外も含む可能性があるため、int(x)として引数が必ず整数になるようにする。 また再帰定義の計算は時間がかかるため、サンプル点を減らす工夫をする。
set xrange[0:10] set samples 11 set logscale y f(x)=(int(x)==0)?1:int(x)*fac(int(x)-1) plot f(x)
2. 漸化式への応用
フィボナッチ数列は以下の漸化式で与えられる。
ただしF0=0, F1=1。 2つの項を設定しないといけないため、3項演算子を2回用いる。
fibo(x)=(int(x)==0)?0:((int(x)==1)?1:fibo(int(x)-1)+fibo(int(x)-2)) unset logscale set xrange[0:11] set samples 12 plot fibo(x) with points
x=11までのプロットだが、再帰定義関数を用いると計算に若干時間がかかる。
3. フラクタルへの応用
マンデルブロ集合のときは初期値z0=0.0と固定していたが、ジュリア集合では初期値はその座標値で与えれる。
関数mandelはx,yはそのままにしてzの引数を再帰定義している。計算を打ち切る条件としてzの絶対値が2を超えるか、代入回数が200以上になるという条件を設定している。
mandel(x,y,z,n)=(abs(z)>2.0 || n>=200) ? n: mandel(x,y,z*z+complex(x,y),n+1) set xrange [-0.5:0.5] set yrange [-0.5:0.5] set logscale z set isosample 50 set hidden3d set contour a= -0.37 b= -0.612 splot mandel(a,b,complex(x,y),0) notitle
A=[-0.37,-0.612]の場合のプロットを下に示す。
bをわずかに変化させたA=[-0.37,-0.6]の場合のプロットを以下に示すが、bがわずかに異なるだけでプロットが大きく変化していることがわかる。
4. まとめ
今回はgnuplotで再帰定義の方法とその応用として漸化式やフラクタルのプロットについて紹介した。再帰定義は漸化式などそのまま扱えるので便利な反面、繰り返し計算回数が増えると時間がかかるというデメリットがある。その場合はCなど別のプログラムで作成して、プロットのみを行ったほうがよい。