python代码画一个图形_如何用Python画各种著名数学图案

大数据文摘作品，转载具体要求见文末

编译团队：Aileen，徐凌霄

用Python绘制著名的数学图片或动画，展示数学中的算法魅力。

本项⽬目将持续更更新，数学中有着太多有趣的事物可以⽤用代码去展示。欢迎提出建议和参与建设！

后台回复“数学”查看完整代码集哦

Mandelbrot 集

代码：46 lines (34 sloc) 1.01 KB'''A fast Mandelbrot set wallpaper rendererreddit discussion: https://www.reddit.com/r/math/comments/2abwyt/smooth_colour_mandelbrot/'''import numpy as npfrom PIL import Imagefrom numba import jitMAXITERS = 200RADIUS = 100@jitdef color(z, i):v = np.log2(i + 1 - np.log2(np.log2(abs(z)))) / 5if v < 1.0:return v**4, v**2.5, velse:v = max(0, 2-v)return v, v**1.5, v**3@jitdef iterate(c):z = 0jfor i in range(MAXITERS):if z.real*z.real + z.imag*z.imag > RADIUS:return color(z, i)z = z*z + creturn 0, 0 ,0def main(xmin, xmax, ymin, ymax, width, height):x = np.linspace(xmin, xmax, width)y = np.linspace(ymax, ymin, height)z = x[None, :] + y[:, None]*1jred, green, blue = np.asarray(np.frompyfunc(iterate, 1, 3)(z)).astype(np.float)img = np.dstack((red, green, blue))Image.fromarray(np.uint8(img*255)).save('mandelbrot.png')if __name__ == '__main__':main(-2.1, 0.8, -1.16, 1.16, 1200, 960)

多米诺洗牌算法

代码链接：https://github.com/neozhaoliang/pywonderland/tree/master/src/domino

正二十面体万花筒

代码：53 lines (40 sloc) 1.24 KB'''A kaleidoscope pattern with icosahedral symmetry.'''import numpy as npfrom PIL import Imagefrom matplotlib.colors import hsv_to_rgbdef Klein(z):'''Klein's j-function'''return 1728 * (z * (z**10 + 11 * z**5 - 1))**5 / \(-(z**20 + 1) + 228 * (z**15 - z**5) - 494 * z**10)**3def RiemannSphere(z):''' map the complex plane to Riemann's sphere via stereographic projection '''t = 1 + z.real*z.real + z.imag*z.imagreturn 2*z.real/t, 2*z.imag/t, 2/t-1def Mobius(z):''' distort the result image by a mobius transformation '''return (z - 20)/(3*z + 1j)def main(imgsize):x = np.linspace(-6, 6, imgsize)y = np.linspace(6, -6, imgsize)z = x[None, :] + y[:, None]*1jz = RiemannSphere(Klein(Mobius(Klein(z))))# define colors in hsv spaceH = np.sin(z[0]*np.pi)**2S = np.cos(z[1]*np.pi)**2V = abs(np.sin(z[2]*np.pi) * np.cos(z[2]*np.pi))**0.2HSV = np.dstack((H, S, V))# transform to rgb spaceimg = hsv_to_rgb(HSV)Image.fromarray(np.uint8(img*255)).save('kaleidoscope.png')if __name__ == '__main__':import timestart = time.time()main(imgsize=800)end = time.time()print('runtime: {:3f} seconds'.format(end - start))

Newton 迭代分形

代码：46 lines (35 sloc) 1.05 KBimport numpy as npimport matplotlib.pyplot as pltfrom numba import jit# define functions manually, do not use numpy's poly1d funciton!@jit('complex64(complex64)', nopython=True)def f(z):# z*z*z is faster than z**3return z*z*z - 1@jit('complex64(complex64)', nopython=True)def df(z):return 3*z*z@jit('float64(complex64)', nopython=True)def iterate(z):num = 0while abs(f(z)) > 1e-4:w = z - f(z)/df(z)num += np.exp(-1/abs(w-z))z = wreturn numdef render(imgsize):x = np.linspace(-1, 1, imgsize)y = np.linspace(1, -1, imgsize)z = x[None, :] + y[:, None] * 1jimg = np.frompyfunc(iterate, 1, 1)(z).astype(np.float)fig = plt.figure(figsize=(imgsize/100.0, imgsize/100.0), dpi=100)ax = fig.add_axes([0, 0, 1, 1], aspect=1)ax.axis('off')ax.imshow(img, cmap='hot')fig.savefig('newton.png')if __name__ == '__main__':import timestart = time.time()render(imgsize=400)end = time.time()print('runtime: {:03f} seconds'.format(end - start))

李代数E8 的根系

代码链接：https://github.com/neozhaoliang/pywonderland/blob/master/src/misc/e8.py