マンデルブロ集合を描画する - Pythonでいろいろやってみる

マンデルブロ集合は漸化式
f:id:T_A_T:20201226205341p:plain
で定義される複素数列で、n→∞の極限において無限大に発散しないという条件を満たす複素数cの集合を指します。

シンプルな定義にも関わらず、マンデルブロ集合を複素平面上に図示することで複雑で美しいフラクタル図形を描画できます。
マンデルブロ集合を求めるプログラムは複素平面において、各点cに対する漸化式を計算し発散するか収束するかで判定します。あらかじめ設定した回数まで漸化式を計算しznがある閾値を越える場合には発散と判定、閾値を越えない場合には収束(マンデルブロ集合)と判定します。

環境

windows10 home
Anaconda 3/ jupyter notebook 5.6.0
Python 3.7.0
OpenCV 4.0.0

コード

画像サイズを600x600とし、幅x・高さyの二重ループ内で対応する複素数が収束するか発散するかを調べています。
yのループはrange(height)つまり0～599の値を取り、複素数の虚部は(y-0.5600)/200として-1.5から1.495まで0.005刻みで変化します。 xのループはrange(width)つまり0～599の値を取り、複素数の虚部は(x-0.7600)/200として-4.2から0.895まで0.005刻みで変化します。このように画像の各画素に複素数cを対応させ、それぞれ場合にznが発散するか収束するかを判定します。
whileループで漸化式をn30まで計算し、計算が完了されれば(ループ終了後のnが30となれば)収束、znが100000000000以上となりn=30となる前にwhileループを抜けた場合は発散と判定します。
収束部の輝度を0、発散部の赤の輝度を255としてマンデルブロ集合を図にします。

import numpy as np
import cv2

height, width = 600, 600
# マンデルブロ集合描画用ndarray
mandelbrot = np.zeros((height, width, 3), np.uint8)  
div = 200

for y in range(height):
    for x in range(width):
        z0 = complex(0, 0)  
        c = complex((x-0.7*width)/div, (y-0.5*height)/div)          
        # z0の絶対値が100000000000を超えるかループが30回になったらwhileループを抜ける
        n = 0
        while abs(z0) < 100000000000 and n < 30:  
            z0 = z0**2 + c  # 漸化式
            n += 1  # ループ回数+1
        if n == 30:  # 漸化式が収束した場合
            mandelbrot[y, x, :] = 0
        else:  # 漸化式が発散した場合
            mandelbrot[y, x, 2] = 255
            
# マンデルブロ集合画像の保存    
cv2.imwrite('mandelbrot.jpg', mandelbrot)

実行結果

マンデルブロ集合(収束部)が黒、それ以外が赤で描画された図形が保存されます。
f:id:T_A_T:20201227071950j:plain

画像を格好よくするために、発散部にグラデーションを付けます。発散時のnを元に、発散部の赤の輝度をnx10により求めます。発散まで回数がかかった場合(nが大きい)は輝度が大きく、すぐに発散した場合(nが小さい)には輝度が小さくなります。

コード

import numpy as np
import cv2

height, width = 600, 600
# マンデルブロ集合描画用ndarray
mandelbrot = np.zeros((height, width, 3), np.uint8)  
div = 200

for y in range(height):
    for x in range(width):
        z0 = complex(0, 0)  
        c = complex((x-0.7*width)/div, (y-0.5*height)/div)          
        # z0の絶対値が100000000000を超えるかループが30回になったらwhileループを抜ける
        n = 0
        while abs(z0) < 100000000000 and n < 30:  
            z0 = z0**2 + c  # 漸化式
            n += 1  # ループ回数+1
        if n == 30:  # 漸化式が収束した場合
            mandelbrot[y, x, :] = 0
        else:  # 漸化式が発散した場合
            mandelbrot[y, x, 2] = n*10
            
# マンデルブロ集合画像の保存  
cv2.imwrite('mandelbrot.jpg', mandelbrot)

実行結果

マンデルブロ集合の周囲が明るく、離れた位置は暗い画像が保存されます。マンデルブロ集合(収束部)近くは複素数の絶対値が比較的小さく漸化式の発散までのnが大きくなりますが、マンデルブロ集合から離れた個所は複素数cが大きいためすぐに発散しnは小さいので輝度も小さくなります。
f:id:T_A_T:20201227082010j:plain 　

計算する複素数の範囲を変更します。
複素数の虚部は0.15から0.524375まで変化させます。複素数の実部は-1.425から-1.050625まで変化させます。

コード

import numpy as np
import cv2

height, width = 600, 600
# マンデルブロ集合描画用ndarray
mandelbrot = np.zeros((height, width, 3), np.uint8)  
div = 1600

for y in range(height):
    for x in range(width):
        z0 = complex(0, 0)  
        c = complex((x-3.8*width)/div, (y+0.4*height)/div)          
        # z0の絶対値が100000000000を超えるかループが30回になったらwhileループを抜ける
        n = 0
        while abs(z0) < 100000000000 and n < 30:  
            z0 = z0**2 + c  # 漸化式
            n += 1  # ループ回数+1
        if n == 30:  # 漸化式が収束した場合
            mandelbrot[y, x, :] = 0
        else:  # 漸化式が発散した場合
            mandelbrot[y, x, 2] = n*10
            
# マンデルブロ集合画像の保存      
cv2.imwrite('mandelbrot.jpg', mandelbrot)