分形与混沌之-Mandelbrot集合

自然界的很多事物，例如树木、云彩、山脉、闪电、雪花以及海岸线等等都呈现出传统的几何学不能描述的形状。这些形状都有如下的特性：
- • 有着十分精细的不规则的结构
- • 整体与局部相似，例如一根树杈的形状和一棵树很像
分形几何学就是用来研究这样一类的几何形状的科学，借助计算机的高速计算和图像显示，使得我们可以更加深入地直观地观察分形几何。

Mandelbrot(曼德布洛特)集合是在复平面上组成分形的点的集合。
Mandelbrot集合的定义(摘自维基百科) Mandelbrot集合可以用下面的复二次多项式定义：

f_c(z) = z^2 + c

其中c是一个复参数。对于每一个c，从z=0开始对函数 fc(z) 进行迭代。序列 (0,fc(0),fc(fc(0)),fc(fc(fc(0))),…) 的值或者延伸到无限大，或者只停留在有限半径的圆盘内。 Mandelbrot集合就是使以上序列不发散的所有c点的集合。

从数学上来讲，Mandelbrot集合是一个复数的集合。一个给定的复数c或者属于Mandelbrot集合，或者不是。
用程序绘制Mandelbrot集合时不能进行无限次迭代，最简单的方法是使用逃逸时间(迭代次数)进行绘制，具体算法如下：

• 判断每次调用函数 fc(z) 得到的结果是否在半径R之内，即复数的模小于R
• 记录下模大于R时的迭代次数
• 迭代最多进行N次
• 不同的迭代次数的点使用不同的颜色绘制

输出图如下:

代码如下:

# -*- coding: utf-8 -*-
# by whyx 2016/6
 
import numpy as np
import pylab as pl
from math import log
import time
from matplotlib import cm
 
def iter_point(c):
     z = c
     for i in xrange(1, 100): # 最多迭代100次
        if abs(z)>2.0: break # 半径大于2则认为逃逸
        z = z*z+c
     return i # 返回迭代次数
 
escape_radius = 2.0
iter_num = 100
 
def smooth_iter_point(c):
    z = c
    for i in xrange(1, iter_num):
        if abs(z)>escape_radius: break
        z = z*z+c
 
    absz = abs(z)
    if absz > 2.0:
        mu = i - log(log(abs(z),2),2)
    else:
        mu = i
    return mu # 返回正规化的迭代次数
 
 
def draw_mandelbrot(cx, cy, d):
 
     # 绘制点(cx, cy)附近正负d的范围的Mandelbrot
 
    (x0, x1, y0, y1) = (cx-d, cx+d, cy-d, cy+d)
 
    y, x = np.ogrid[y0:y1:200j, x0:x1:200j]
 
    c = x + y*1j
    start = time.clock()
    mandelbrot = np.frompyfunc(smooth_iter_point,1,1)(c).astype(np.float)
    print "time=",time.clock() - start  #时间单位 秒
    pl.imshow(mandelbrot, cmap=cm.Blues_r, extent=[x0,x1,y0,y1])
    pl.gca().set_axis_off()
 
x,y = 0.27322626, 0.595153338
 
pl.subplot(231)
draw_mandelbrot(-0.5, 0, 1.5)
for i in range(2,7):
     pl.subplot(230+i)
     draw_mandelbrot(x, y, 0.2**(i-1))
 
pl.subplots_adjust(0.02, 0, 0.98, 1, 0.02, 0)
pl.show()