PCが描く奇妙な画像集(数学的万華鏡と生物形態等の世界)

・インタープリタBASICによるフラクタルとカオスの奇妙な画集。

009. プログラム実行どきのエラー処理

2014-06-13 08:08:33 | 画像作成の説明
プログラムを実行ときにエラーが発生したとき其のエラーの内容及びエラー発生行を表示させるプログラムをサブルーチン化しておくと便利である。

以下はBASIC/98 での其のサブルーチン・プログラムである。

50000 *エラー処理
50010 'file name ER1
50020 IF ERR=6 OR ERR=11 THEN 50080 ELSE 50030
50030 PRINT "「overflow」及び「ゼロで除算」以外のエラー"
50040 PRINT "エラー番号=";ERR
50050 PRINT "エラー発生行=";ERL
50060 PRINT "J=";J;"K=";K
50070 STOP
50080 RESUME NEXT

008. 複素関数:f(Z)^g(Z)の実数部と虚数部の計算

2014-06-13 07:47:33 | 画像作成の説明


f(Z)とg(Z)を任意の複素関数とするとき、f(Z)^g(Z)の実数部と虚数部を求めるプログラムをサブルーチン化しておく。ここで、このサブルーチンの入力値として、f(Z)及びg(Z)の、おのおのの実数部、虚数部は以下のように与えられているとする。

f(Z)の実数部、虚数部を、それぞれ、FZR、FZI とする。
即ち、f(Z)=FZR + i*FZI
g(Z)の実数部、虚数部を、それぞれ、GZR、GZI とする。
即ち、g(Z)=GZR + i*GZI

また、f(Z)は極座標でも表して、f(Z)= =FZR + i*FZI=R*e^iθ とする。
ここで、R=(FZR^2+FZI^2)^0.5 , θ=arctan(FZI/FZR) である。

上記の条件のもとに、f(Z)^g(Z)の実数部と虚数部を、以下のようにして求める。(以下の式において、掛け算の記号(*)は煩雑さを避けるため適宜省略する)

f(Z)^g(Z)=(R*e^iθ)^( GZR + i*GZI)
=( R^(GZR+i*GZI) )*( e^(iθ(GZR+iGZI) )
=( R^GZR)* ( e^(-θ*GZI) ) * (R^i*GZI)*(e^iθGZR)
ここで、ZF1= (R^GZR)* ( e^(-θ*GZI) とおけば、
f(Z)^g(Z)=ZF1* (R^i*GZI)*(e^iθGZR)
=ZF1*( e^(logR))^(iGZI) ) * (e^iθGZR) (ここで、logの底はeとする)
    =ZF1*e^i*(logR*GZI+θ*GZR)
ここで、ZF2=logR*GZI+θ*GZR とおけば、
f(Z)^g(Z)=ZF1*e^i*ZF2
=ZFI*(cosZF2+i*sinZF2)
となり、f(Z)^g(Z)の実数部=ZF1*cosZF2
f(Z)^g(Z)の虚数部=ZF1*sinZF2
となる。ここで、ZF1= (R^GZR)* ( e^(-θ*GZI)
        ZF2=logR*GZI+θ*GZR 
である。従って、f(Z)^g(Z)の実数部と虚数部を求めるプログラムをサブルーチン・プロクグラムは以下のようにすればよい。

30000 REM f(Z)^g(Z)の計算:file name FGZ
30010 IF R=0 THEN 30060
30020 ZF1=(R^GZR)*EXP(-TH*GZI)
30030 ZF2=TH*GZR+LOG(R)*GZI
30040 FZGZR=ZF1*COS(ZF2)
30050 FZGZI=ZF1*SIN(ZF2)
30060 RETURN

006. 複素関数:Z^f(Z)の実数部と虚数部の計算

2014-06-13 07:32:08 | 画像作成の説明


Zを複素平面上の点とし、f(Z)を任意の複素関数とするとき、複素関数:Z^f(Z)の実数部と虚数部を求める計算プログラムをサブルーチン化ておく。

そのために、Z^f(Z)の実数部と虚数部を以下のようにして求める。

Z=X+iY=R*e^iθとし,ここで、R=(X^2+Y^2)^0.5 ,θ=arctan(Y/X)から、Rとθは、このサブルーチンの入力値とする。また、任意のf(Z)の実数部(これをREFとする)と虚数部(これをIMFとする)も既に与えられているとする。

即ち、f(Z)=FZR+i*FZI とする。 (ここで、FZR , FZI は既知の実数とする)

(以下の計算式において掛け算の記号(*)は煩雑さを避けるため適宜省略する。)

Z^f(Z)=( R*e^iθ)^( FZR+i*FZI)=( (R*e^iθ)^FZR)*( (R*e^iθ)^iFZI )
=(R^FZR)*(e^(iθ*FZR))*(R^iFZI)*e^(-θ*FZI)
=(R^FZR)*e^(-θ*FZI)*e^(iθ*FZR)*R^(i*FZI)
  ここで、ZF1=(R^FZR)*e^(-θ*FZI) とおくと、

Z^f(Z)=ZF1* ( e^(iθ*FZR)*R^(i*FZI) ) ・・・・・(1)
  ここで、R=e^logR  (但し、logの底はeとする) だから
  R^(i*FZI)=(e^logR)^(i*FZI)=e^(i*logR*FZI)
従って、(1)より、
Z^f(Z)=ZF1*e^i(θ*FZR+logR*FZI)
ここで、ZF2=θ*FZR+logR*FZI とおくと、
Z^f(Z)=ZF1*e^(*ZF2)=ZF1*(cosZF2+i*sinZF2)
=ZF1*cosZF2 + i*(ZF1*sinZF2)
となるから、
Z^f(Z)の実数部= ZF1*cosZF2
Z^f(Z)の虚数部= ZF1*sinZF2
( 但し、ZF1=(R^FZR)*e^(-θ*FZI) , ZF2=θ*FZR+logR*FZI )

となり、Z^f(Z)の実数部と虚数部が求まった。
従って、Z^f(Z)の実数部と虚数部を、それぞれ、ZFZR,ZFZIとすれば、求めるサブルーチン・プログラムは下記のようになる。

20000 REM Z^(f(z)の計算:file name ZFZ
20010 IF R=0 THEN 20060
20020 ZF1=(R^FZR)*EXP(-TH*FZI)
20030 ZF2=TH*FZR+LOG(R)*FZI
20040 ZFZR=ZF1*COS(ZF2)
20050 ZFZI=ZF1*SIN(ZF2)
20060 RETURN

005. 複素関数:Z^Zの実数部と虚数部の計算

2014-06-13 07:12:09 | 画像作成の説明


複素平面上の点Zが与えられていたとき、複素関数:Z^Zの実数部と虚数部の計算をサブルーチン化しておく。 ここで、このサブルーチンを実行するときは、点ZのX,Y,R及びθ(プログラムでの変数名:TH)は、与えられているものとする。また、自然対数eの底(プログラムでの変数名:EE)も与えられているものとする。

X及びYから、R及びθは、R=(X^2+Y^2)^0.5 , θ=arctan(Y/X) より求める。(注:θ=arctan(Y/X) は、あらかじめ、サブルーチン ARCTAN3A で求めておく。)

複素関数:Z^Zの実数部と虚数部の計算方法を以下に示す。
(これも私が数学の公式集を参考にして計算したものであるから誤りがあるかも知れないから御注意!! : 下式で掛け算の記号(*)は適宜付けたり省力したりしている。下式はゴチャゴチャしているが、やっていることは単純なこと。こういうところが数式というものの欠点の一つだと思う。)

Z^Z={R*e^(iθ)}^(X+iY)={R^(X+iY)}*[{e^(iθ)}^{(X+iY)}]
=[(R^X)*{R^(iY)}]*[e^{iθ*(X+iY)}
=[(R^X)*e^(-θY)]*[(R^iY)*(e^iθX)]・・・(1)

ここで、上式を分かりやすくするために、
A=(R^X)*e^(-θY)
とおくと、
Z^Z=A*[(R^iY)*(e^iθX)]・・・(2)

ここで、e^x=R として、x を求める。logの底をeとすると、
logR=log(e^x)=x*log e=x

従って、R=e^x=e^logR だから、(1)式の第2項に、このRを代入すると、
(1)式の第2項=[(R^iY)*(e^iθX)]={(e^logR)^iY}*{(e^iθX)}
={e^i(YlogR)} *{(e^iθX)}
=e^i(YlogR+θX)・・・(3)

ここで、上式を分かりやすくするために、
B=(YlogR+θX)
とおくと、
(1)式の第2項=(3)式=e^iB=cosB+isinB
となるから、(1)式、即ち、Z^Z は、
Z^Z=A*(3)式=A*( cosB+isinB)=AcosB +i(AsinB)
となる。従って、
Z^Zの実数部= AcosB
Z^Zの虚数部= AsinB
となる。ここで、
A=(R^X)*e^(-θY)  , B=(YlogR+θX)

以上より、Z^Zの実数部と虚数部を計算するサブルーチン・プログラムは、下記のようにします。下記のプログラムの行9010,9020,9040,9060及び9070は、使用しているBASICの制約によるもので、本来は不要です。

9000 REM Z^ZのRe.ZとIm.Zの計算
9002 REM ZのX,Y,R,θ(TH)は与えられているものとする。
9004 REM 自然対数の底e:EE=2.718281も与えられているものとする。
9010 IF TH*Y<-709 THEN TH*Y=-709 :GOTO 9030
9020 IF TH*Y>709 THEN TH*Y=709
9030 A1=(R^X)*EE^(-TH*Y)
9040 IF R=0 THEN R=1E^-100
9050 B=Y*LOG(R)+TH*X
9060 IF B>1.0E+9 THEN B=1.0E+9 :GOTO 9080
9070 IF B<-1.0E+9 THEN B=-1.0E+9
9080 ZZR=A1*COS(B)
9090 ZZI=A1*SIN(B)
9100 RETURN

上記のプログラムにおいて、Z^Zの実数部値をZZRとおき、Z^Zの虚数部値をZZIとした、それらをメイン・プログラムへ、RETURNすればよい。