分形造型(一)
一、引言
欧氏几何方法
以规则几何图形为研究对象,物体形状由方程来描述
欧氏几何主要描述工具
点,直线,平滑的曲线、平面
欧氏几何不能模拟自然景物
山脉、云、树、草、火、雾、波浪等
为解决复杂图形生成问题,分形造型应运而生
二、分形的历史
分形(Fractal)一词,是曼德布罗特创造出来的,其原意具有不规则、支离破碎等意义
分形几何学是一门以非规则几何形态为研究对象的几何学
分形理论创始人是美籍法国数学家B.B.曼德布罗特。他1924年11月20日生于波兰华沙。1985年获巴纳德(Barnard)奖章。1986年获富兰克林(Franklin)奖章。1967年发表于美国《科学》杂志上的“英国的海岸线有多长”的划时代论文,是他的分形思想萌芽的重要标志。1975年冬他创造了fractal一词。著作《分形:形状、机遇和维数》法文版于同年出版。
1967年B.B.Mandelbrot提出了“英国的海岸线有多长?”的问题,以1km为单位测量海岸线,得到的近似长度将短于1km的迂回曲折都忽略掉了,若以1m为单位测量,则能测出被忽略掉的迂回曲折,长度将变大,测量单位进一步变小,测得的长度将愈来愈大,这些愈来愈大的长度将趋近于一个确定值,这个极限值就是海岸线的长度。海岸线的长度是不确定的,长度依赖于测量单位。
三、分形几何与传统几何相比
-
从整体上看,分形几何图形是处处不规则的。例如,海岸线和山川形状,从远距离观察,其形状是极不规则的。
-
在不同尺度上,图形的规则性又是相同的。上述的海岸线和山川形状,从近距离观察,其局部形状又和整体形态相似,它们从整体到局部,都是自相似的。
四、分数维D定义
N相当于把原图复制出若干份
如果某图形是由把原图缩小为\(\frac 1 S\)的相似的N个图形所组成,有:
\]
的关系成立,则指数D称为相似性维数,D可以是整数,也可以是分数。
五、经典分形
Koch曲线
import turtle as t
def koch(size, n):
if n == 0:
t.fd(size)
else:
for angle in [0, 60, -120, 60]:
t.left(angle)
koch(size / 3, n - 1)
if __name__ == \'__main__\':
t.setup(1000, 1000)
t.pen(speed=0, pendown=False, pencolor=\'blue\')
a, n = 400, 4
t.goto(-a / 2, a / 2 / pow(3, 0.5))
t.pd()
for i in range(3):
koch(a, n)
t.right(120)
t.ht()
t.done()
原来的一条边变成4条边,所以N=4
现在的一个小三角的边长是原来的\(\frac 1 3\),所以S=\(\frac 1 3\)
所以\(D=\frac {ln4}{ln3}=1.2618\),D介于1和2之间
Cantor集
Cantor集是由不断去掉线段中间的三分之一开区间而得出。
#!/apps/bin/python
from tkinter import *
class Cantor(Frame):
limit = 1
def __init__(self, master=None):
Frame.__init__(self, master)
self.draw = Canvas(self, width=800, height=500)
self.grid()
self.createWidgets()
def createWidgets(self):
# self.quitButton = Button(self, text="QUIT", command=self.quit)
# self.quitButton.pack(side=BOTTOM, fill=BOTH)
self.draw.pack(side=LEFT)
self.drawCanvas(100, 100, 700, 100)
def drawCanvas(self, ax, ay, bx, by):
self.draw.create_line(ax, ay, bx, by)
if (bx - ax) > self.limit:
cx = ax + (bx - ax) / 3
cy = ay + 50
dx = bx - (bx - ax) / 3
dy = by + 50
ay = ay + 50
by = by + 50
self.drawCanvas(ax, ay, cx, cy)
self.drawCanvas(dx, dy, bx, by)
app = Cantor()
app.master.title("Cantor (recursive)")
app.mainloop()
一条线变成2条线,所以D=2,长度是原来的\(\frac 1 3\),所以S=\(\frac 1 3\)
所以
\]
Sierpinski垫片
将一个等边三角形四等分,得到四个小等边三角形,去掉中间一个,保留它的三条边,将剩下的三个小等边三角形在分别进行四等分,并分别去掉中间一个,保留它们的边。
"""
功能:绘制谢尔宾斯基三角形
环境:python3.7
日期:2019/1/14 21:49
作者:指尖魔法师
版本:1.0
"""
import turtle as t
def sanjiaoxing(san):
"""
传入三个点坐标,绘制三角形
"""
t.penup()
t.goto(san[0])
t.pendown()
t.goto(san[1])
t.goto(san[2])
t.goto(san[0])
def get_mid(a, b):
"""
计算返回2个点的中间点坐标
"""
x = (a[0] + b[0]) / 2
y = (a[1] + b[1]) / 2
return [x, y]
def draw_san(size, i):
"""
绘制谢尔宾斯基三角形函数
:param size: 三个点坐标列表
:param i: 递归次数
"""
# 绘制三角形
sanjiaoxing(size)
if i > 0:
# 绘制左边小三角形
size2 = [size[0], get_mid(size[0], size[1]), get_mid(size[0], size[2])]
draw_san(size2, i - 1)
# 绘制上边的小三角形
size3 = [get_mid(size[0], size[2]), get_mid(size[1], size[2]), size[2]]
draw_san(size3, i - 1)
# 绘制右边的小三角形
size4 = [get_mid(size[0], size[1]), size[1], get_mid(size[1], size[2])]
draw_san(size4, i - 1)
def main():
"""
主函数
"""
# 打印图形标题
t.penup()
t.left(90)
t.forward(350)
t.pendown()
t.write("谢尔宾斯基三角形", False, align="center", font=("宋体", 20, "normal"))
t.speed(5)
# 初始三角形坐标
points = [[-200, 0], [200, 0], [0, 300]]
# 递归5次
count = 5
# 调用绘制谢尔宾斯基三角形函数
draw_san(points, count)
t.exitonclick()
if __name__ == \'__main__\':
main()
\]
Sierpinski海绵
import javax.swing.*;
import java.awt.*;
import java.util.stream.IntStream;
public class Sponge extends JFrame implements Runnable {
public static void main(String[] args) {
Sponge mf = new Sponge();
mf.initUI();
}
Point[] pn = new Point[20];
public void paint(Graphics g) {
super.paint(g);
//正方形的顶点坐标
int x = 100, y = 500;
int d = 120, dx = 81, dy = 60;
Point p0 = new Point(x, y);
drawBottomLevel(g, p0, d, dx, dy);
Point p1 = new Point(x, y - d);
drawMidLevel(g, p1, d, dx, dy);
Point p2 = new Point(x, y - 2 * d);
drawTopLevel(g, p2, d, dx, dy);
IntStream.range(0, pn.length).forEach(i -> {
Point ptem = pn[i];
if (ptem != null) {
Point ptem0 = new Point(pn[i].x, pn[i].y + 2 * d / 3);
drawBottomLevel(g, ptem0, d / 3, dx / 3, dy / 3);
Point ptem1 = new Point(pn[i].x, pn[i].y + d / 3);
drawMidLevel(g, ptem1, d / 3, dx / 3, dy / 3);
Point ptem2 = new Point(pn[i].x, pn[i].y);
drawTopLevel(g, ptem2, d / 3, dx / 3, dy / 3);
}
});
}
//根据起点p0得到另外6个点
public Point[] getPointByP0(Point p0, int d, int dx, int dy) {
Point p1 = new Point(p0.x + dx, p0.y - dy);
Point p2 = new Point(p0.x + dx + d, p0.y - dy);
Point p3 = new Point(p0.x + d, p0.y);
Point p4 = new Point(p0.x + d, p0.y + d);
Point p5 = new Point(p0.x, p0.y + d);
Point p6 = new Point(p0.x + dx + d, p0.y + d - dy);
Point[] ps = new Point[7];
ps[0] = p0;
ps[1] = p1;
ps[2] = p2;
ps[3] = p3;
ps[4] = p4;
ps[5] = p5;
ps[6] = p6;
return ps;
}
int count = 0;//递归次数
public void draw(Graphics g, Point p0, Point p1, Point p2, Point p3, Point p4, Point p5, Point p6) {
if (count < 20) {
pn[count] = p0;
}
count++;
Graphics2D gD = (Graphics2D) g;
//填充面的颜色
Polygon pon1 = new Polygon();
//填充第一个可见面
pon1.addPoint(p0.x, p0.y);
pon1.addPoint(p1.x, p1.y);
pon1.addPoint(p2.x, p2.y);
pon1.addPoint(p3.x, p3.y);
gD.setColor(new Color(0xFFC0CB));
gD.fillPolygon(pon1);
//填充第二个可见面
Polygon pon2 = new Polygon();
pon2.addPoint(p0.x, p0.y);
pon2.addPoint(p3.x, p3.y);
pon2.addPoint(p4.x, p4.y);
pon2.addPoint(p5.x, p5.y);
gD.setColor(new Color(0x39c5bb));
gD.fillPolygon(pon2);
//填充第三个可见面
Polygon pon3 = new Polygon();
pon3.addPoint(p2.x, p2.y);
pon3.addPoint(p3.x, p3.y);
pon3.addPoint(p4.x, p4.y);
pon3.addPoint(p6.x, p6.y);
gD.setColor(new Color(0xFFA500));
gD.fillPolygon(pon3);
}
public void drawBottomLevel(Graphics g, Point p0, int d, int dx, int dy) {
Point[] ps1 = getPointByP0(p0, d, dx, dy);
Point[] ps2 = getPointByP0(ps1[3], d, dx, dy);
Point[] ps3 = getPointByP0(ps2[3], d, dx, dy);
Point[] ps4 = getPointByP0(ps3[1], d, dx, dy);
Point[] ps5 = getPointByP0(ps4[1], d, dx, dy);
Point[] ps6 = getPointByP0(ps1[1], d, dx, dy);
Point[] ps7 = getPointByP0(ps6[1], d, dx, dy);
Point[] ps8 = getPointByP0(ps7[3], d, dx, dy);
draw(g, ps7[0], ps7[1], ps7[2], ps7[3], ps7[4], ps7[5], ps7[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps8[0], ps8[1], ps8[2], ps8[3], ps8[4], ps8[5], ps8[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps5[0], ps5[1], ps5[2], ps5[3], ps5[4], ps5[5], ps5[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps4[0], ps4[1], ps4[2], ps4[3], ps4[4], ps4[5], ps4[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps6[0], ps6[1], ps6[2], ps6[3], ps6[4], ps6[5], ps6[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps1[0], ps1[1], ps1[2], ps1[3], ps1[4], ps1[5], ps1[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps2[0], ps2[1], ps2[2], ps2[3], ps2[4], ps2[5], ps2[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
draw(g, ps3[0], ps3[1], ps3[2], ps3[3], ps3[4], ps3[5], ps3[6]);
try {
Thread.sleep(100);
} catch (InterruptedException e) {
e.printStackTrace();
}
}
public void drawMidLevel(Graphics g, Point p0, int d, int dx, int dy) {
Point[] ps1 = getPointByP0(p0, d, dx, dy);
Point[] ps2 = getPointByP0(ps1[3], d, dx, dy);
Point[] ps3 = getPointByP0(ps2[3], d, dx, dy);
Point[] ps4 = getPointByP0(ps3[1], d, dx, dy);
Point[] ps5 = getPointByP0(ps4[1], d, dx, dy);
Point[] ps6 = getPointByP0(ps1[1], d, dx, dy);
Point[] ps7 = getPointByP0(ps6[1], d, dx, dy);
Point[] ps8 = getPointByP0(ps7[3], d, dx, dy);
draw(g, ps7[0], ps7[1], ps7[2], ps7[3], ps7[4], ps7[5], ps7[6]);
draw(g, ps5[0], ps5[1], ps5[2], ps5[3], ps5[4], ps5[5], ps5[6]);
draw(g, ps1[0], ps1[1], ps1[2], ps1[3], ps1[4], ps1[5], ps1[6]);
draw(g, ps3[0], ps3[1], ps3[2], ps3[3], ps3[4], ps3[5], ps3[6]);
}
public void drawTopLevel(Graphics g, Point p0, int d, int dx, int dy) {
Point[] ps1 = getPointByP0(p0, d, dx, dy);
Point[] ps2 = getPointByP0(ps1[3], d, dx, dy);
Point[] ps3 = getPointByP0(ps2[3], d, dx, dy);
Point[] ps4 = getPointByP0(ps3[1], d, dx, dy);
Point[] ps5 = getPointByP0(ps4[1], d, dx, dy);
Point[] ps6 = getPointByP0(ps1[1], d, dx, dy);
Point[] ps7 = getPointByP0(ps6[1], d, dx, dy);
Point[] ps8 = getPointByP0(ps7[3], d, dx, dy);
draw(g, ps7[0], ps7[1], ps7[2], ps7[3], ps7[4], ps7[5], ps7[6]);
draw(g, ps8[0], ps8[1], ps8[2], ps8[3], ps8[4], ps8[5], ps8[6]);
draw(g, ps5[0], ps5[1], ps5[2], ps5[3], ps5[4], ps5[5], ps5[6]);
draw(g, ps4[0], ps4[1], ps4[2], ps4[3], ps4[4], ps4[5], ps4[6]);
draw(g, ps6[0], ps6[1], ps6[2], ps6[3], ps6[4], ps6[5], ps6[6]);
draw(g, ps1[0], ps1[1], ps1[2], ps1[3], ps1[4], ps1[5], ps1[6]);
draw(g, ps2[0], ps2[1], ps2[2], ps2[3], ps2[4], ps2[5], ps2[6]);
draw(g, ps3[0], ps3[1], ps3[2], ps3[3], ps3[4], ps3[5], ps3[6]);
}
public void initUI() {
this.setTitle("门格海绵");
this.setSize(1000, 800);
this.setDefaultCloseOperation(3);
this.setResizable(false);
this.setLocationRelativeTo(null);
this.setVisible(true);
}
@Override
public void run() {
}
}
Julia集
Python实现
import matplotlib.pyplot as plt
import numpy as np
from matplotlib import cm
import matplotlib.backend_bases as mbb
import time
class MandelBrot():
def __init__(self, x0, x1, y0, y1, n):
self.oriAxis = np.array([x0, x1, y0, y1]) # 初始坐标
self.axis = self.oriAxis
self.nx, self.ny, self.nMax = n, n, n # x,y方向的网格划分个数
self.nIter = 100 # 迭代次数
self.n0 = 0 # 预迭代次数
self.z = genZ(self.oriAxis, self.nx, self.ny)
self.DrawMandelbrot()
def DrawMandelbrot(self):
mBrot = getJulia(self.z, self.z, self.nIter)
self.fig, ax = plt.subplots()
plt.imshow(mBrot, cmap=cm.jet, extent=self.axis)
plt.gca().set_axis_off()
self.fig.canvas.mpl_disconnect(self.fig.canvas.manager.key_press_handler_id)
self.fig.canvas.mpl_connect(\'button_press_event\', self.OnMouse)
self.fig.canvas.mpl_connect(\'button_release_event\', self.OnRelease)
self.fig.canvas.mpl_connect(\'scroll_event\', self.OnScroll)
plt.show()
def DrawJulia(self, c0):
z = genZ([-2, 2, -2, 2], 800, 800)
julia = getJulia(z, c0, self.nIter)
jFig, jAx = plt.subplots()
plt.cla()
plt.imshow(julia, cmap=cm.jet, extent=self.axis)
plt.gca().set_axis_off()
plt.show()
jFig.canvas.draw_idle()
# 滚轮缩放
def OnScroll(self, evt):
x0, y0 = evt.xdata, evt.ydata
if evt.button == "up":
self.axis = (self.axis + [x0, x0, y0, y0]) / 2
elif evt.button == \'down\':
self.axis = 2 * self.axis - [x0, x0, y0, y0]
z = genZ(self.axis, self.nx, self.ny)
mBrot = getJulia(z, z, self.nIter)
plt.cla()
plt.imshow(mBrot, cmap=cm.jet, extent=self.axis)
plt.gca().set_axis_off()
mBrot[mBrot < 1] == self.n0 + self.nIter
self.n0 = int(np.min(mBrot))
self.fig.canvas.draw_idle()
pass
def OnMouse(self, evt):
self.xStart = evt.xdata
self.yStart = evt.ydata
self.fig.canvas.draw_idle()
def OnRelease(self, evt):
x0, y0, x1, y1 = self.xStart, self.yStart, evt.xdata, evt.ydata
if evt.button == mbb.MouseButton.LEFT:
self.DrawJulia(x1 + y1 * 1j) # 如果释放的是左键,那么就绘制Julia集并返回
return
# 右键拖动,可以对Mandelbrot集进行真实的放大
self.axis = np.array([min(x0, x1), max(x0, x1),
min(y0, y1), max(y0, y1)])
nxny = self.axis[[1, 3]] - self.axis[[0, 2]]
self.nx, self.ny = (nxny / max(nxny) * self.nMax).astype(int)
z = genZ(self.axis, self.nx, self.ny)
n = 100 # n为迭代次数
mBrot = getJulia(z, z, n)
plt.cla()
plt.imshow(mBrot, cmap=cm.jet, extent=self.axis)
plt.gca().set_axis_off()
mBrot[mBrot < 1] == self.n0 + n
self.n0 = int(np.min(mBrot))
self.fig.canvas.draw_idle()
def genZ(axis, nx, ny):
x0, x1, y0, y1 = axis
x = np.linspace(x0, x1, nx)
y = np.linspace(y0, y1, ny)
real, img = np.meshgrid(x, y)
z = real + img * 1j
return z
def getJulia(z, c, n, n0=0, m=2):
t = time.time()
c = np.zeros_like(z) + c
out = abs(z)
for _ in range(n0):
z = z * z + c
for i in range(n0, n0 + n):
absz = abs(z)
z[absz > m] = 0
c[absz > m] = 0
out[absz > m] = i
z = z * z + c
print("time:", time.time() - t)
return out
if __name__ == "__main__":
x, y = 0, 0
brot = MandelBrot(-2, 1, -1.5, 1.5, 1000)
介绍
在复平面上,水平的轴线代表实数,垂直的轴线代表虚数。每个Julia集合(有无限多个点)都决定一个常数C,它是一个复数。在复平面上任意取一个点,其值是复数Z。将其代入下面方程中进行反复迭代运算:
\]
Julia集的生成过程
假设显示器的分辨率是\(a*b\)点,可显示的颜色为\(K+1\)种,分别以数字\(0,1,2,\cdots,K\)表示,且以0表示黑色。
从\(Z_i\)到\(Z_{i+1}\)的迭代过程就是
\]
步骤0:
选定参数\(x_{min}=y_{min}=-1.5, x_{max}=y_{max}=1.5\),又取M=100
令
\]
对所有的点\((n_x,n_y),n_x=0,1,\cdots,a-1\)及\(n_y=0,1,\cdots,b-1\)完成下面的循环
步骤1:
令
\]
步骤2:
由迭代过程从\((x_i,y_i)\)算出\((x_{i+1},y_{i+1})\),并计算\(k:=k+1\)
步骤3:
计算\(r=x_i^2+y_i^2\)。
如果\(r>M\),则选择颜色\(k\),转至步骤4;
如果\(k=K\),则选择颜色\(0\),转至步骤4;
如果\(r≤M\),且\(k<K\),则转至步骤2。
步骤4:
对点\((n_x,n_y)\)显示颜色\(k\)并转至下一点,再从头作步骤1