第６章

6.1　ジャパニーズ・アトラクタ

【リスト6-1】 ジャパニーズ・アトラクタ
from tkinter import *
import time, math
#■繰返し回数（点の数に一致:多くすると精細な図形になるが遅くなる）
MAX=5000; NSEP=800
def intX(X) : return int(X * 40 + 50.5)
def intY(Y) : return int(150.5-Y * 10)
def errP(str, j,i,B):
    print(str+
　　　　　("の計算でオーバフローが起きました。j=%d, i=%d,  B=%5.2f " %(j,i,B))+
          str+"を0とします。")
def Map(k, B,DT): #■表示用計算
    global MAX, NSEP
    Ctab=["#FF0000","#007700","#0000FF","777700","#770077","#007777",
　　　　　"#770000"]
    T=0; X = 0.2;Y = 0.1; LOOP=MAX+200-1;  NK=1
    for j in range(LOOP): #表示要素数での繰返し
        if j % 200==0:#描画中の表示
            canvas.delete("DT")
            canvas.create_text(100,30,
                               text="描画中です(%2d %%)" % int(j/LOOP*100),
                               font=('Times',10),fill='red', tag="DT")
            canvas.create_text(100,50,text="■"*NK,
                               font=('Times',10),fill='red', tag="DT")
            NK+=1
            if NK>20: NK=1
            root.update()
        for i in range(NSEP):#ジャパニーズ・アトラクタの計算
            XN = X + Y * DT 
            YN = Y + (-k * Y - X * X * X + B * math.cos(T)) * DT
            T += DT ; X = XN ; Y = YN
            if math.isnan(X): errP("X", j,i,B);X=0
            if math.isnan(Y): errP("X", j,i,B);Y=0
        if j>=200:#ある程度(200回)計算した後を表示
            XX=intX(X); YY=intY(Y); C=Ctab[j%7]
            canvas.coords(FID[j-200],XX,YY,XX,YY) #表示
            canvas.itemconfig(FID[j-200], fill=C)
    canvas.delete("DT") #描画中表示削除
    root.update()
#■待ち処理等
def wait():
    tk.update_idletasks();time.sleep(0.0001)
def strdt(DT):
    return "B = %+5.2f" % DT
#■ マウスイベントによる次の描画
def leftMouseDown(event): 
    global B, DB, k, DT; B+=DB
    if B>12.105:B = 11.7
    canvas.itemconfig(strID,text=strdt(B))
    Map(k,B,DT)
#■メイン処理
root=Tk(); root.title("Japanese Attractor") #キャンバスの定義
canvas=Canvas(root, width=250,height=250)
canvas.pack()
FID=[0 for j in range(MAX)]                 #表示用の点の準備
for j in range(MAX):
    FID[j]=canvas.create_rectangle(0,0,0,0,fill="#007700", outline='')
k=0.1; B=11.7; DT=2*math.pi/NSEP; DB=0.02
strID=canvas.create_text(100,10,text=strdt(B), font=('ＭＳ ゴシック',10))
canvas.bind("<Button-1>",leftMouseDown)  #開始
Map(k, B, DT)
root.mainloop()

6.2　ローレンツ・モデル

【リスト6-2】 ローレンツ・モデル
from tkinter import *
import time,math
a = 10; b = 3; dt = 0.004 ; g = 0
#a = 10; b = 3; dt = 0.004 ; g = 10
dg = 1; XS = 30; YS = 20; ZS = 30
def screenX(X, Y) : return 5 * (X - 0.6 * Y) + 200
def screenY(Y, Z) : return 350 - 5 * (Z - 0.6 * Y)
def Line3D(X1, Y1, Z1, X2, Y2, Z2, C, w):
    global canvas
    XS1 = int(screenX(X1, Y1) + 0.5); YS1 = int(screenY(Y1, Z1) + 0.5)
    XS2 = int(screenX(X2, Y2) + 0.5); YS2 = int(screenY(Y2, Z2) + 0.5)
    canvas.create_line(XS1, YS1, XS2, YS2, fill=C, width=w, tag="DT")
def Coord3D():
    Line3D(-30,   0,   0, 50,  0,  0, "#FF0000",2)
    Line3D(  0, -30,   0,  0, 50,  0, "#007700",2)
    Line3D(  0,   0, -30,  0,  0, 50, "#000077",2)
def dspValue():
    global g
    canvas.create_text(250, 20 ,text="g = %7.4f" % g,
                       font=('ＭＳ ゴシック',10),fill="#FFFF00", tag="DT")
    canvas.create_text(250, 40 ,
                       text="黄緑：XY平面への投影,青色:XZ平面への投影, 黄色:3D座標の位置",
                       font=('ＭＳ ゴシック',10),fill="#FFFF00", tag="DT")
#■表示用計算
def cmpIter():
    global XS, YS, ZS, dt, a, b,g; canvas.delete("DT")
    dspValue(); X= XS; Y = YS; Z = ZS; Coord3D()
    for i in range(6000):
        XN = X + dt * a * (Y - X)
        YN = Y + dt * (-X * Z + g * X - Y)
        ZN = Z + dt * (X * Y - b * Z)
        Line3D(X, Y, 0, XN, YN,  0, "#00FF00", 1)
        Line3D(X, 0, Z, XN,  0, ZN, "#0000FF", 1)
        Line3D(X, Y, Z, XN, YN, ZN, "#FFFF00", 1)
        X = XN; Y = YN; Z = ZN
    root.update()
#■待ち処理等
def wait(): tk.update_idletasks(); time.sleep(0.0001)
def strdt(DT): return "B = %+5.2f" % DT
#■ マウスイベントによる次の描画
def leftMouseDown(event): 
    global g, dg
    if abs(g-30)<0.00001: g = 0
    g = g + dg ; cmpIter()
#■メイン処理
root=Tk()
root.title("Lorenz Model")
canvas=Canvas(root, width=500,height=500,bg="#000000")
canvas.pack()
#画面キャプチャするときは以下の２行のコメントを外す
#canvas.bind("<Button-1>",leftMouseDown)
#root.mainloop()
while True:
    cmpIter()
    root.update_idletasks()
    root.update(); time.sleep(0.01)
    if abs(g-30)<0.00001: break
    g = g + dg ; cmpIter()

6.3　ロジスティック写像

【リスト6-3】 ロジスティック写像
from tkinter import *
BITX= 500 # X方向ビット数
BITY= 500 # Y方向ビット数
def array(N1,N2):
    return [[0 for i in range(N2)] for j in range(N1)]
SCX=[ 300,  600, 1200, 1800, 2400]                      # X方向スケール
SCY=[ 200,  400,  800, 1200, 1600]                      # Y方向スケール
X0 =[ 800, 1900, 4200, 6200, 8400]                      # X始点
Y0 =[ 300,  450,  800, 1100, 1500]                      # Y始点
ast =[ 2.9,  3.15,  3.5, 3.4, 3.45]                     # 開始 a の値
da =[0.00001, 0.00001, 0.000001, 0.000001, 0.000001]    # a の値の刻み幅
XS=0.1                                                  # Xの開始位置
def createBM(): # ビットマップ生成(1×1の正方形をビットに代替させる)
    global canvas, BITX, BITY, BM
    BM=array(BITX,BITY)
    for i in range(BITX):
        for j in range(BITY):
            BM[i][j]=canvas.create_rectangle(i,j,i,j,
                                             fill="#000077", outline='')
def clearBM(): # ビットマップクリア
    global canvas, BITX, BITY,BM
    for i in range(BITX):
        for j in range(BITY):
            if canvas.itemcget(BM[i][j],option='fill')!="#000077":
                canvas.itemconfig(BM[i][j],fill="#000077")
def intX(X): global ID, SCX, X0; return int(X * SCX[ID] - X0[ID]) # 画面X座標
def intY(X): global ID, SCX, Y0; return int(Y0[ID] - X * SCY[ID]) # 画面Y座標
def cmpIter(): # 写像計算と描画
    global XS, ID,BM
    a=ast[ID]; X= XS ; N=int((4-a)/da[ID]); checkN=int(N/20); k=0; NK=1
    while a<4: #繰返し計算
        k+=1
        if (k % checkN) ==0: # 描画中の表示
            canvas.delete("DT")
            canvas.create_text(200,20, text="描画中です(%2d %% a=%5.2f)"
　　　　　　　　　　　　　　　 % (int(k/N*100),a),
                               font=('Times',10),fill='yellow', tag="DT")
            canvas.create_text(200,40,text="■"*NK,
                               font=('Times',10),fill='yellow', tag="DT")
            NK+=1
            if NK>20: NK=1
            root.update()
        i=intX(a); j=intY(X) # 描画 
        if i>=0 and i<BITX and j>=0 and j<BITY: 
	    canvas.itemconfig(BM[i][j],fill="#FFFFFF")
        X = a*(1-X)*X; a+=da[ID]
    canvas.delete("DT")
def leftMouseDown(event): # マウスイベントによる次の描画
    global ID; ID+=1
    if ID>4: ID=0
    clearBM(); cmpIter()
root=Tk(); root.title("ロジスティック写像")
canvas=Canvas(root,width=BITX,height=BITY); canvas.pack()
ID=0
createBM(); cmpIter(); root.update() # 初期描画
canvas.bind("<Button-1>",leftMouseDown)
root.mainloop()

6.4 エノン写像

【リスト6-4】 エノン写像
from tkinter import *
import time,math
BITX=200; BITY=200; a=1.5; b=0.25;XS=1; YS=1; ID=0
def array(N1,N2): return [[0 for i in range(N1)] for j in range(N2)]
BM=array(BITX,BITY)
#■ビットマップ生成
def initBM():#(ビットマップ1ピクセルを1×1のサイズの正方形で代替する)
    global BM, BITX,BITY, canvas
    for i in range(BITX):
        for j in range(BITY):
            BM[i][j]=canvas.create_rectangle(i,j,i,j, fill="#000077", outline='')
#■ビットマップクリア
def clearBitmap():#1×1のサイズの正方形の色を変更することでクリアとする
    global BM, BITX,BITY, canvas
    for i in range(BITX):
        for j in range(BITY):#以下、同じ色なら変更しないことで高速化
            if canvas.itemcget(BM[i][j],option='fill')!="#000077":
                canvas.itemconfig(BM[i][j],fill="#000077")
SCX = [  70, 200,  400,  800,  1600] # X方向スケール
SCY = [ 140, 400,  800, 1600,  3200] # Y方向スケール
X0  = [ 100, -70, -300, -800, -1750] # 画面上の原点 X
Y0  = [ 100, 100,  120,  145,   170] # 画面上の原点 Y
#■画面座標系に変換
def intX(X): global SCX, X0,ID; return int(X * SCX[ID] + X0[ID])
def intY(Y): global SCY, Y0,ID; return int(Y0[ID] - Y * SCY[ID])
#■表示用計算
def cmpIter():
    global XS,YN,canvas, BM,a,b
    canvas.create_text(100,20,text="画面をクリアしています" ,
                       font=('Times',10),fill='yellow', tag="DT")
    root.update()
    clearBitmap()
    canvas.delete("DT")
    canvas.create_text(100,20,text="描画中です" ,
                       font=('Times',10),fill='yellow', tag="DT")
    root.update()
    #表示処理
    X= XS;Y=YS;LOOP=100000
    for k in range(LOOP):
        XN=1-a*X*X+Y;YN=b*X
        i=intX(XN); j=intY(YN)
        if i>=0 and i<BITX and j>=0 and j<BITY:#同じ色なら変更する必要なし
            if canvas.itemcget(BM[i][j],option='fill') !="#FFFFFF":
                canvas.itemconfig(BM[i][j],fill="#FFFFFF")
        X=XN; Y=YN;
    #表示終了処理
    canvas.delete("DT")
    root.update()
#■ マウスイベントによる次の描画
def leftMouseDown(event): 
    global ID; ID+=1
    if ID>4:ID=0
    cmpIter()
#■メイン処理
root=Tk();root.title("Henon Map")
canvas=Canvas(root, width=BITX,height=BITY); canvas.pack()
initBM(); root.update()
cmpIter()
canvas.bind("<Button-1>",leftMouseDown)
root.mainloop()

6.5 グモウスキーとミラの写像

【リスト6-5】 グモウスキーとミラの写像
from tkinter import *
import time
#■繰返し回数（点の数に一致:多くすると精細な図形になるが遅くなる）
MAX=5000
#■表示用計算
def Map(mu, a, b):
    X=0.1; XX=X*X; Y=0.0; mm2=2*(1.0-mu)
    G = mu * X + mm2 * XX / (1.0 + XX)
    for i in range(MAX):
        XP = int(500.0 - X * 10); YP = int(200.0 - Y * 10)
        canvas.coords(FID[i],XP,YP,XP,YP) #表示
        XN = X; X = Y + a * (1.0 - b * Y * Y) * Y + G
        XX=X*X
        G = mu * X + mm2 * XX / (1.0 + XX); Y = -XN + G
    tk.update()
#■待ち処理等
def wait():
    tk.update_idletasks(); time.sleep(0.00001)
def strdt(mu): return "μ = {:+f}".format(float(mu))
#■メイン処理
#キャンバスの定義
tk=Tk(); tk.title("Gumowski-Mira Map")
canvas=Canvas(tk, width=1000,height=400,bg='white'); canvas.pack()
#表示用の点の準備
FID=[0 for j in range(MAX)]
for j in range(MAX):
    FID[j]=canvas.create_rectangle(0,0,0,0,fill="#007700", outline='')
#アニメーション
mu=-1; a=0.008; b=0.05; dm=0.0001
strID=canvas.create_text(100,10,text=strdt(mu), font=('ＭＳ ゴシック',10))
while True:
    mu+=dm
    if mu>1.0: mu=-1.0
    canvas.itemconfig(strID,text=strdt(mu))
    Map(mu,a,b); wait()
        
【リスト6-6】　式(6)によるグモウスキーとミラの写像
from tkinter import *
import time
#■繰返し回数（点の数に一致:多くすると精細な図形になるが遅くなる）
MAX=10000
lamb=2.4;  a=0.4;
#■表示用計算
def fun(x):
    global lamb, a; xx=abs(x); iflag=1
    if x<0: iflag=-1; 
    if xx<a: return (0.2*lamb*xx)*iflag
    b=0.2*lamb*a;
    if xx<(1-a):
        return iflag*(0.45*lamb*(xx*xx-a*a-xx+a)/
                      (0.5-a)+1.1*lamb*(xx-a)+b)
    c=1.1*lamb*(1-2*a)+b;
    return iflag*(2*lamb*(xx+a-1)+c)
def Map():
    global lamb; XX=0.2; YY=0.2;
    for i in range(10000):
        X=XX; Y=YY;	XX=(1-lamb)*X+Y; YY=Y-fun(X);
        XP=int(XX*80+200); YP=int(200-YY*80)
        canvas.coords(FID[i],XP,YP,XP,YP)
    tk.update()
#■待ち処理等
def wait():
    tk.update_idletasks(); time.sleep(0.00001)
def strdt(mu): return "λ = %8.4f" % mu
#■メイン処理
tk=Tk(); tk.title("Gumowski-Mira Other Map No1")
canvas=Canvas(tk, width=400,height=400,bg='white');canvas.pack()
FID=[0 for j in range(MAX)]
for j in range(MAX):
    FID[j]=canvas.create_rectangle(0,0,0,0,fill="#007700", outline='')
dl=0.0001
strID=canvas.create_text(100,10,text=strdt(lamb), font=('ＭＳ ゴシック',10))
while True:
    lamb+=dl
    if lamb>2.6: lamb=2.4
    canvas.itemconfig(strID,text=strdt(lamb))
    Map(); wait()

【リスト6-7】　式(7)によるグモウスキーとミラの写像
from tkinter import *
import math
#■繰返し回数（点の数に一致:多くすると精細な図形になるが遅くなる）
MAX=10000
#fai=1.2566; sigma=1.15      #例題用の値
#fai=2.0944; sigma=1.0730
#fai=1.5708; sigma=1.1500
fai=1.8850; sigma=1.1300
#■表示用計算
def Map():
    global fai, sigma; XX=0.2; YY=0.2
    for i in range(10000):
        X=XX; Y=YY
        XX=sigma*(X*math.cos(fai)-Y*math.sin(fai))+X*Y
        YY=sigma*(X*math.sin(fai)+Y*math.cos(fai))+X*X+Y*Y*Y
        XP=int(XX*100+200); YP=int(200-YY*100)
        canvas.coords(FID[i],XP,YP,XP,YP)
    tk.update()
#■メイン処理
#キャンバスの定義
tk=Tk(); tk.title("Gumowski-Mira Other Map No2")
canvas=Canvas(tk, width=400,height=400,bg='white');canvas.pack()
#表示用の点の準備
FID=[0 for j in range(MAX)]
for j in range(MAX):
    FID[j]=canvas.create_rectangle(0,0,0,0,fill="#007700", outline='')
Map(); tk.mainloop()

6.6 リアプノフ・フラクタル

【リスト6-6】 リアプノフ・フラクタル
# ■ Lyapunov fractalの描画
from tkinter import *
import math
#MD=5
#STR=[0,0,1,0,1] # MDの個数　(0がA, 1がBを示す。この場合AABABを示す。
MD=7
STR=[0,1,1,0,1,0,1] # MDの個数　(0がA, 1がBを示す。この場合ABBABABを示す。
#MD=2
#STR=[0,1] # MDの個数　(0がA, 1がBを示す。この場合ABを示す。

# ■ Lyapunov指数の計算
def lamda(a, b):
    global MD,STR
    try:
        x0=0.5; N=1000; lam=0
        for n in range(1,N):
            if STR[n % MD]==0: rn = a
            else :rn = b
            xn = rn * x0 * (1-x0)；rr=abs(rn*(1-2*xn))
            if rr !=0 : lam+= math.log(rr)
            x0=xn
        return lam/N
    except :
        print("error", a,b)； return 0

MAX=200; DV=int(MAX/4)
root=Tk(); root.title("Lyapunov")
canvas=Canvas(root,width=MAX,height=MAX); canvas.pack()
k=0; NK=1; checkN=int(MAX/20)
for i in range(MAX):
    k+=1
    if (k % checkN) ==0: # 描画中の表示
        canvas.delete("DT")
        canvas.create_text(100,20,text="描画中です(%2d %%)" % int(k/MAX*100),
                           font=('Times',10),fill='#FF0000', tag="DT")
        canvas.create_text(100,40,text="■"*NK,
                           font=('Times',10),fill='#FF0000', tag="DT")
        NK+=1
        if NK>20: NK=1
        root.update()
    a=i/DV
    for j in range(MAX):
        b=j/DV
        if a==0 or b==0 :lm=0 
        else : lm=lamda(a,b)
        if lm>0:
            CC=int(255*lm)
            if CC>255: CC=255
            C="#FF%02X%02X" % (CC,CC)
        else:
            CC=int(255*math.exp(lm))
            if CC>255: CC=255
            C="#%02X%02XFF" % (CC,CC)
        canvas.create_rectangle(i,j,i,j,fill=C,outline="")
canvas.delete("DT")
root.update()
root.mainloop()

6.7 燃える船フラクタル

【リスト6-6】 リアプノフ・フラクタル
from tkinter import *
#import time

#■配列宣言
#
#  N1,N2  : 配列の大きさ
#  戻り値 : ゼロ設定された配列(A[N1][N2]を宣言したこととなる）
#           実装的には[N2個のリスト，N2個のリスト，・・・]となる。
#
def array2(N1,N2):
    return [[0 for j in range(N2)]for k in range(N1)]

#■画素初期設定(画素をピクセル1の正方形で代表させる)
#
#  canvas : Canvasクラスオブジェクト
#  NX     : X方向最大数+1を指定すること
#  NY     : Y方向最大数+1を指定すること
#  戻り値 : 画素IDの配列 
#
def initDt(canvas, NX,NY):
    A=array2(NX,NY)
    for i in range(NX):
        for j in range(NY):
            A[i][j]=canvas.create_rectangle(i,j,i,j,fill='blue', outline='')   
    return A

#■Tk初期設定
#
#  Title  : 最初に表示するウィンドウのタイトル
#  戻り値 : Tkクラスオブジェクト 
#
def initTk(Title):
    tk=Tk()
    tk.title("Mandelbrot")
    tk.resizable(0,0)
    return tk

#■キャンバス初期設定
#
#  tk : Tkクラスオブジェクト
#  W  : キャンバスの幅
#  H  : キャンバスの高さ
#  戻り値 : Canvasクラスオブジェクト
#
def initCanvas(tk, W,H):
    canvas=Canvas(tk,width=W,height=H,highlightthickness=0)
    canvas.pack()
    return canvas
 
#■繰返し演算
#
#  C    　:  漸化式　Zn = Zn^2+C のC
#  LOOP 　:  収れん判定のための繰返し回数
#  戻り値 :  収れん回数/発散検出回数
#　　　　　　0 のとき，LOOP回数では収れん/発散を判定できなかったことを示す。
#
def Iter(C, LOOP):
    Z=complex(0)
    for i in range(1,LOOP+1):
        # 以下を ZZ=Z*Z*Z+C としても興味深い画像になる
        Z=complex(abs(Z.real), abs(Z.imag))
        ZZ=Z*Z+C
        if(abs(ZZ)>=4):
            return -i
        if(abs(ZZ-Z)<0.001):
            return i
        Z=ZZ
    return 0

def drawFrac(canvas, CID,DX,X0,Y0, XE,YE, LOOP,MD,Ctab1,Ctab2):
    X=X0
    XMAX=len(CID); YMAX=len(CID[0])
    while X<=XE:
        #実行途上表示
        if int(X/DX)%5==0:
            tk.update()
        Y=Y0
        while Y<=YE:
            #繰返し回数による色の選定
            k=Iter(complex(X,Y),LOOP)
            if k>=0:
                CL=Ctab1[k % MD]
            else:
                #-----------------------------------------
                #以下3行、発散回数で色分け
                #CC=-k
                #if CC>255:CC=255
                #CL="#%02X%02X%02X" % (CC,0,int(128-CC/2))
                #-----------------------------------------
                #以下1行、発散回数の剰余で色分け
                CL=Ctab2[k % MD]
            #表示用座標設定と色変更
            XX=int((X-X0)/DX+0.5)
            YY=int((Y-Y0)/DX+0.5)
            if XX>=0 and XX<XMAX and YY>=0 and YY<YMAX:
                canvas.itemconfig(CID[XX][YY],fill=CL)
            #Y増分
            Y+=DX
        #X増分
        X+=DX

#■実行メイン
Ctab1=['black','#000077','#0000FF','#0077FF','#00FFFF']
Ctab2=['#773300','#770000','#FF0000','#FF7700','#FFFF00']
#Ctab2 =['#FF0000','#FF0000','#FF0000','#FF0000','#FF0000']
tk=initTk("Burning Ship")
canvas=initCanvas(tk,400,400)
CID=initDt(canvas,401,401)
tk.title("Burning Shipの描画中")
#drawFrac(canvas,CID, DX , X0 , Y0   , XE ,  YE, LOOP,MD, Ctab1,Ctab2):
drawFrac(canvas,CID,0.01,-2.5,-2.0  , 1.5,  2 , 200 , 5, Ctab1,Ctab2)
#drawFrac(canvas, CID,0.0004,-0.84, -1.00, -0.68, -0.84, 500 , 5, Ctab1,Ctab2)
#drawFrac(canvas,CID,0.00025,-1.8,-0.08  , -1.7,  0.02 , 500 , 5, Ctab1,Ctab2)
tk.update()
tk.title("Burning Shipの描画終了")
tk.mainloop()