Python分形几何可视化—— 复数迭代、L系统与生物分形模拟
本节将深入探索分形几何的奇妙世界，实现Mandelbrot集生成器和L系统分形树工具，并通过肺部血管分形案例展示分形在医学领域的应用。我们将使用Python的NumPy进行高效计算，结合Matplotlib进行高质量可视化。
前言：探索无限细节的几何宇宙
分形几何（Fractal Geometry）作为20世纪数学最激动人心的发现之一，打破了欧几里得几何的传统框架，揭示了自然界中普遍存在的「自相似性」与「无限细节」。从蜿蜒的海岸线到茂密的树冠，从显微镜下的雪花到宇宙星系的分布，分形无处不在，连接着数学的抽象之美与现实世界的复杂系统。
本文将带领读者走进分形几何的可视化世界，通过三个核心模块——复数迭代生成的Mandelbrot集、符号重写系统驱动的L分形树、生物医学中的肺部血管模拟——展示分形理论的数学本质与跨学科应用。我们将结合Python编程工具，利用NumPy的高效计算和Matplotlib的可视化能力，实现从数学公式到动态图像的完整转化，同时通过动画和案例揭示分形在科学研究中的实际价值。
一、基本概念：分形与L系统理论基础
1.1 分形几何的核心特征
**分形（Fractal）**由数学家本华·曼德博（Benoît Mandelbrot）于1975年提出，其定义可概括为：
- 自相似性（Self-Similarity）：局部与整体在形态、功能或信息上具有相似性（如蕨类叶片的分支结构）；
- 分数维度（Fractional Dimension）：几何对象的复杂度无法用整数维描述（如科赫曲线的分形维数为1.2618）；
- 迭代生成（Iterative Generation）：通过简单规则的重复应用创造复杂结构（如Mandelbrot集的复数迭代）。
典型分形示例：
- Mandelbrot集：通过复数迭代 ( z_{n+1} = z_n^2 + c ) 生成的无限复杂图形，边界处的每个点对应一个动态系统的长期行为；
- 科赫雪花（Koch Snowflake）：通过递归替换线段生成的分形曲线，具有无限周长和有限面积；
- 蕨类植物：自然界中典型的自相似分形，可用仿射变换的迭代函数系统（IFS）模拟。
1.2 L系统：符号重写的分形生成器
**L系统（L-System）**由生物学家阿里斯蒂德·林登迈耶（Aristid Lindenmayer）于1968年提出，最初用于模拟植物的生长过程。其核心要素包括：
- 字母表（Alphabet）：由符号组成的集合，分为绘图符号（如
F表示向前移动并绘制线段）和控制符号（如+、-表示旋转角度）； - 公理（Axiom）：初始字符串，代表分形的最基础形态；
- 规则（Rules）：符号替换的映射关系，通过迭代逐步构建复杂结构。
L系统工作流程：
- 从公理出发，根据规则逐代替换字符串中的每个符号；
- 通过「龟绘图」（Turtle Graphics）解释器将最终字符串转换为几何图形：
F：向前移动指定距离并绘制线段；+/-：向左/右旋转指定角度；[/]：保存/恢复当前画笔状态（用于处理分支结构）。
应用场景：
- 植物形态模拟（如蕨类、树木）；
- 晶体生长、神经网络建模；
- 生物医学中的支气管树、血管分支模拟。
二、Mandelbrot集生成器
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib.colors import ListedColormap
import matplotlib.animation as animation
from numba import jit
plt.rcParams["font.sans-serif"] = ["SimHei", "WenQuanYi Micro Hei", "Heiti TC"]
plt.rcParams["axes.unicode_minus"] = False # 正确显示负号
plt.rcParams["mathtext.fontset"] = "cm" # 设置数学字体
# 使用Numba加速计算
@jit(nopython=True, cache=False)
def mandelbrot(c, max_iter=100):
"""
计算Mandelbrot集的迭代次数
参数:
c: 复数
max_iter: 最大迭代次数
返回: 迭代次数 (如果发散) 或 max_iter (如果收敛)
"""
z = 0j
for i in range(max_iter):
z = z * z + c
if abs(z) > 2:
return i
return max_iter
@jit(nopython=True, cache=False)
def compute_mandelbrot(width, height, x_min, x_max, y_min, y_max, max_iter=100):
"""
计算整个Mandelbrot集
参数:
width, height: 图像尺寸
x_min, x_max: 实部范围
y_min, y_max: 虚部范围
max_iter: 最大迭代次数
返回: 迭代次数矩阵
"""
# 创建图像网格
x = np.linspace(x_min, x_max, width)
y = np.linspace(y_min, y_max, height)
iterations = np.zeros((height, width), dtype=np.int32)
for i in range(height):
for j in range(width):
c = x[j] + 1j * y[i]
iterations[i, j] = mandelbrot(c, max_iter)
# 返回迭代次数和统计信息（用于调试）
return iterations, iterations.min(), iterations.max(), np.mean(iterations)
def plot_mandelbrot(
width=800,
height=800,
x_min=-2.0,
x_max=1.0,
y_min=-1.5,
y_max=1.5,
max_iter=100,
cmap="viridis",
title="Mandelbrot集",
output_path="mandelbrot.png",
):
"""
绘制Mandelbrot集
参数:
width, height: 图像尺寸
x_min, x_max, y_min, y_max: 复数平面范围
max_iter: 最大迭代次数
cmap: 颜色映射
title: 图像标题
output_path: 输出文件路径
"""
# 接收完整返回值并提取迭代次数矩阵
result = compute_mandelbrot(width, height, x_min, x_max, y_min, y_max, max_iter)
iterations = result[0] # 提取迭代次数矩阵
# 创建图像
plt.figure(figsize=(10, 10), dpi=100)
# 使用对数缩放增强细节
log_iter = np.log(iterations + 1)
# 绘制Mandelbrot集
plt.imshow(log_iter, extent=(x_min, x_max, y_min, y_max), cmap=cmap, origin="lower")
# 设置标题和标签
plt.title(title, fontsize=18, pad=20)
plt.xlabel("实部", fontsize=14)
plt.ylabel("虚部", fontsize=14)
# 添加颜色条
cbar = plt.colorbar(label="迭代次数 (对数缩放)")
cbar.set_ticks([]) # 隐藏刻度，因为是对数缩放
# 添加数学公式
plt.text(
0.05,
0.95,
r"$z_{n+1} = z_n^2 + c$",
transform=plt.gca().transAxes,
fontsize=16,
bbox=dict(facecolor="white", alpha=0.8),
)
# 保存图像
plt.savefig(output_path, bbox_inches="tight", dpi=150)
plt.close()
print(f"Mandelbrot集图像已保存至: {output_path}")
def mandelbrot_zoom_animation(output_path="mandelbrot_zoom_fixed.gif"):
"""创建修复后的Mandelbrot集缩放动画"""
# 初始参数
width, height = 400, 400
max_iter = 200
# 创建鲜艳的颜色映射
colors = [
"#000000", # 黑
"#0000FF", # 蓝
"#00FFFF", # 青
"#00FF00", # 绿
"#FFFF00", # 黄
"#FF8000", # 橙
"#FF0000", # 红
"#FFFFFF", # 白
]
cmap = ListedColormap(colors)
# 设置更平滑的缩放路径
zoom_path = [
# 起始位置
{"center": (-0.5, 0.0), "scale": 1.0, "pause": 5},
# 逐步放大
{"center": (-0.5, 0.0), "scale": 0.3, "pause": 3},
{"center": (-0.5, 0.0), "scale": 0.1, "pause": 3},
# 移动到海马谷区域
{"center": (-0.745, 0.1), "scale": 0.03, "pause": 5},
{"center": (-0.745, 0.1), "scale": 0.01, "pause": 7},
{"center": (-0.745, 0.1), "scale": 0.003, "pause": 10},
]
fig = plt.figure(figsize=(8, 8), dpi=100)
ax = plt.gca()
# 初始化图像
init_data = np.zeros((height, width))
im = ax.imshow(init_data, cmap=cmap, origin="lower")
# 添加颜色条
cbar = fig.colorbar(im, ax=ax)
cbar.set_label("迭代深度")
plt.title("Mandelbrot集缩放动画", fontsize=16, pad=15)
plt.xlabel("实部", fontsize=12)
plt.ylabel("虚部", fontsize=12)
# 动画帧列表 (存储数据和范围)
frames = []
for zoom in zoom_path:
cx, cy = zoom["center"]
scale = zoom["scale"]
# 计算当前视图范围
x_min = cx - 1.5 * scale
x_max = cx + 1.5 * scale
y_min = cy - 1.5 * scale
y_max = cy + 1.5 * scale
# 计算Mandelbrot集
iterations = compute_mandelbrot(
width, height, x_min, x_max, y_min, y_max, max_iter
)
# 检查返回值类型
if not isinstance(iterations, np.ndarray):
print(f"警告: compute_mandelbrot返回了非数组类型: {type(iterations)}")
if isinstance(iterations, tuple):
print(
f"元组长度: {len(iterations)}, 元素类型: {[type(e) for e in iterations]}"
)
# 尝试从元组中提取数组 (假设第一个元素是数据)
if len(iterations) > 0 and isinstance(iterations[0], np.ndarray):
iterations = iterations[0]
else:
raise TypeError("无法从元组中提取有效的迭代数据")
log_iter = np.log(iterations + 1)
# 计算当前帧的合适颜色范围
vmin = np.min(log_iter)
vmax = np.max(log_iter)
# 避免vmin==vmax导致错误
if vmin == vmax:
vmax = vmin + 0.1
# 添加帧数据 (图像数据、范围和颜色范围)
frame_data = {
"data": log_iter,
"extent": (x_min, x_max, y_min, y_max),
"vmin": vmin,
"vmax": vmax,
}
# 添加多帧以暂停
for _ in range(zoom["pause"]):
frames.append(frame_data)
print(
f"生成帧：{len(frames)}，范围：({x_min:.5f}, {x_max:.5f}, {y_min:.5f}, {y_max:.5f})"
)
# 动画更新函数
def update(frame):
im.set_data(frame["data"])
im.set_extent(frame["extent"])
im.set_clim(frame["vmin"], frame["vmax"])
# 修复f-string括号问题
zoom_factor = 3 / (frame["extent"][1] - frame["extent"][0])
ax.set_title(
f"Mandelbrot集缩放动画 - 缩放因子: {zoom_factor:.0f}x", fontsize=14
)
return (im,)
# 创建动画
anim = animation.FuncAnimation(
fig, update, frames=frames, interval=100, blit=True, repeat_delay=1000
)
# 保存GIF
anim.save(output_path, writer="pillow", fps=10)
plt.close()
print(f"修复后的Mandelbrot集缩放动画已保存至: {output_path}")
三、L系统分形树生成工具
class LSystem:
"""L系统分形生成器"""
def __init__(self, rules, axiom, angle, distance, iterations):
"""
初始化L系统
参数:
rules: 替换规则字典
axiom: 初始字符串
angle: 旋转角度 (度)
distance: 初始步长
iterations: 迭代次数
"""
self.rules = rules
self.axiom = axiom
self.angle = np.radians(angle)
self.distance = distance
self.iterations = iterations
self.current_string = axiom
self.stack = []
self.positions = []
self.directions = []
# 生成分形字符串
self.generate_string()
def generate_string(self):
"""生成L系统字符串"""
for _ in range(self.iterations):
new_string = ""
for char in self.current_string:
new_string += self.rules.get(char, char)
self.current_string = new_string
def generate_points(self):
"""根据字符串生成点坐标"""
# 初始状态
x, y = 0, 0
direction = np.pi / 2 # 初始向上 (90度)
points = [(x, y)]
# 距离衰减因子
dist = self.distance
dist_factor = 0.6 # 每级分支长度衰减
for i, char in enumerate(self.current_string):
if char == 'F': # 向前移动
x += dist * np.cos(direction)
y += dist * np.sin(direction)
points.append((x, y))
elif char == '+': # 左转
direction += self.angle
elif char == '-': # 右转
direction -= self.angle
elif char == '[': # 保存状态
self.stack.append((x, y, direction, dist))
elif char == ']': # 恢复状态
if self.stack:
x, y, direction, dist = self.stack.pop()
points.append((x, y))
# 分支级别减少时增加距离
dist /= dist_factor
elif char == 'f': # 向前移动但不画线 (用于定位)
x += dist * np.cos(direction)
y += dist * np.sin(direction)
points.append((x, y))
elif char == '|': # 反转方向
direction += np.pi
return points
def plot_fractal_tree(rules, axiom, angle, distance, iterations,
title="L系统分形树", output_path="fractal_tree.png"):
"""
绘制L系统分形树
参数:
rules: 替换规则字典
axiom: 初始字符串
angle: 旋转角度 (度)
distance: 初始步长
iterations: 迭代次数
title: 图像标题
output_path: 输出文件路径
"""
# 创建L系统
lsys = LSystem(rules, axiom, angle, distance, iterations)
# 生成点坐标
points = lsys.generate_points()
# 创建图像
plt.figure(figsize=(10, 12), dpi=120)
# 绘制分形树
x, y = zip(*points)
plt.plot(x, y, 'k-', linewidth=0.8, alpha=0.7)
# 设置标题和标签
plt.title(title, fontsize=16, pad=15)
plt.axis('equal')
plt.axis('off')
# 添加L系统信息
info_text = (f"迭代次数: {iterations}\n"
f"旋转角度: {angle}°\n"
f"初始字符串: {axiom}\n"
f"替换规则: {str(rules)}")
plt.text(0.05, 0.05, info_text, transform=plt.gca().transAxes, fontsize=10,
bbox=dict(facecolor='white', alpha=0.8))
# 保存图像
plt.savefig(output_path, bbox_inches='tight')
plt.close()
print(f"分形树图像已保存至: {output_path}")
def plot_tree_growth(axiom, rules, angle, distance, max_iter=5,
output_path="tree_growth.gif"):
"""
创建分形树生长动画
"""
fig = plt.figure(figsize=(8, 10), dpi=100)
plt.axis('equal')
plt.axis('off')
# 动画更新函数
def update(iter_num):
plt.clf()
plt.axis('equal')
plt.axis('off')
# 创建当前迭代的L系统
lsys = LSystem(rules, axiom, angle, distance, iter_num)
points = lsys.generate_points()
# 绘制分形树
x, y = zip(*points)
plt.plot(x, y, 'k-', linewidth=0.8, alpha=0.7)
# 设置标题
plt.title(f"分形树生长: 迭代 {iter_num}", fontsize=14, pad=10)
return fig,
# 创建动画
anim = animation.FuncAnimation(fig, update, frames=range(1, max_iter+1),
interval=800, blit=True)
# 保存GIF
anim.save(output_path, writer='pillow', fps=2)
plt.close()
print(f"分形树生长动画已保存至: {output_path}")
四、跨学科案例：肺部血管分形模拟
def simulate_lung_vasculature(output_path="lung_vasculature.png"):
"""
模拟肺部血管分形结构
参数:
output_path: 输出文件路径
"""
# 创建图像
plt.figure(figsize=(12, 10), dpi=130)
# 初始主干血管
start_point = (0, 0)
end_point = (0, 10)
plt.plot([start_point[0], end_point[0]], [start_point[1], end_point[1]],
'r-', linewidth=3.0, alpha=0.8)
# 递归绘制血管分支
def draw_branches(start, direction, length, width, angle, depth, max_depth=7):
if depth > max_depth:
return
# 计算分支点
branch_point = (
start[0] + length * np.cos(direction),
start[1] + length * np.sin(direction)
)
# 绘制分支
plt.plot([start[0], branch_point[0]], [start[1], branch_point[1]],
color='darkred', linewidth=width, alpha=0.7)
# 递归绘制子分支 (左右两个方向)
new_length = length * 0.7 # 长度衰减
new_width = width * 0.6 # 宽度衰减
# 左分支
left_angle = direction + angle
draw_branches(branch_point, left_angle, new_length, new_width, angle, depth+1, max_depth)
# 右分支
right_angle = direction - angle
draw_branches(branch_point, right_angle, new_length, new_width, angle, depth+1, max_depth)
# 设置血管参数
initial_angle = np.pi / 2 # 初始向上
branch_angle = np.radians(35) # 分支角度
initial_length = 3.0 # 初始分支长度
initial_width = 1.5 # 初始分支宽度
# 从主干开始绘制
draw_branches(end_point, initial_angle, initial_length, initial_width, branch_angle, 1)
# 添加对称的另一侧
draw_branches(end_point, initial_angle + np.pi, initial_length, initial_width, branch_angle, 1)
# 设置标题和标签
plt.title("肺部血管分形结构模拟", fontsize=18, pad=20)
plt.axis('equal')
plt.axis('off')
# 添加分形信息
plt.text(0.05, 0.05,
"分形特性:\n- 自相似性\n- 递归分支\n- 尺度不变性",
transform=plt.gca().transAxes, fontsize=12,
bbox=dict(facecolor='white', alpha=0.8))
# 添加医学说明
plt.text(0.7, 0.15,
"医学意义:\n- 最大化气体交换面积\n- 优化血液流动效率\n- 最小化血液流动阻力",
transform=plt.gca().transAxes, fontsize=12,
bbox=dict(facecolor='white', alpha=0.8))
# 保存图像
plt.savefig(output_path, bbox_inches='tight')
plt.close()
print(f"肺部血管模拟图已保存至: {output_path}")
def generate_lung_bronchi(output_path="lung_bronchi.png"):
"""
使用L系统模拟肺部支气管分形结构
"""
# 定义L系统规则
rules = {
'F': 'FF', # 主干生长
'X': 'F-[[X]+X]+F[+FX]-X' # 分支规则
}
axiom = 'X'
angle = 25 # 分支角度
distance = 1.0 # 初始步长
iterations = 6 # 迭代次数
# 创建L系统
lsys = LSystem(rules, axiom, angle, distance, iterations)
# 生成点坐标
points = lsys.generate_points()
# 创建图像
plt.figure(figsize=(10, 12), dpi=130)
# 绘制支气管树
x, y = zip(*points)
plt.plot(x, y, 'b-', linewidth=1.5, alpha=0.7)
# 设置标题和标签
plt.title("肺部支气管分形结构 (L系统模拟)", fontsize=16, pad=15)
plt.axis('equal')
plt.axis('off')
# 添加分形参数
info_text = (f"迭代次数: {iterations}\n"
f"分支角度: {angle}°\n"
f"L系统规则: {str(rules)}")
plt.text(0.05, 0.05, info_text, transform=plt.gca().transAxes, fontsize=10,
bbox=dict(facecolor='white', alpha=0.8))
# 保存图像
plt.savefig(output_path, bbox_inches='tight')
plt.close()
print(f"支气管分形图已保存至: {output_path}")
五、使用示例
if __name__ == "__main__":
# 1. Mandelbrot集
# 标准Mandelbrot集
plot_mandelbrot(output_path="mandelbrot_standard.png")
# 高细节区域
plot_mandelbrot(x_min=-0.745, x_max=-0.735,
y_min=0.095, y_max=0.105,
max_iter=500, cmap='hot',
title="Mandelbrot集细节 (海马谷)",
output_path="mandelbrot_detail.png")
# Mandelbrot集缩放动画
mandelbrot_zoom_animation(output_path="mandelbrot_zoom.gif")
# 2. L系统分形树
# 二叉分形树
tree_rules1 = {
'F': 'FF',
'X': 'F-[[X]+X]+F[+FX]-X'
}
plot_fractal_tree(rules=tree_rules1, axiom='X', angle=25,
distance=1.0, iterations=6,
title="二叉分形树 (L系统)",
output_path="binary_tree.png")
# 雪花分形
snowflake_rules = {
'F': 'F+F--F+F'
}
plot_fractal_tree(rules=snowflake_rules, axiom='F--F--F', angle=60,
distance=1.0, iterations=5,
title="科赫雪花 (L系统)",
output_path="koch_snowflake.png")
# 分形树生长动画
tree_rules2 = {
'F': 'FF',
'X': 'F[+X]F[-X]+X'
}
plot_tree_growth(axiom='X', rules=tree_rules2, angle=20,
distance=1.0, max_iter=6,
output_path="tree_growth.gif")
# 3. 肺部血管分形模拟
# 递归模拟
simulate_lung_vasculature(output_path="lung_vasculature.png")
# L系统模拟
generate_lung_bronchi(output_path="lung_bronchi.png")
六、生成图像说明
-
Mandelbrot集 (
mandelbrot_standard.png,mandelbrot_detail.png):- 展示完整的Mandelbrot集和细节放大
- 使用对数缩放增强视觉细节
- 包含迭代公式的数学表示
- 缩放动画展示分形的无限复杂性
-
L系统分形树 (
binary_tree.png,koch_snowflake.png):- 二叉分形树的自然形态
- 科赫雪花的经典分形结构
- 包含L系统规则和参数说明
- 生长动画展示分形构建过程
-
肺部血管模拟 (
lung_vasculature.png,lung_bronchi.png):- 递归方法模拟的肺部血管结构
- L系统生成的支气管分形树
- 展示自相似性和分支模式
- 包含医学意义说明
七、教学要点
-
分形几何基础:
- 分形的定义与特征：自相似性、分数维度
- Mandelbrot集的数学原理：复数迭代
- 分形维度的计算方法
-
L系统理论:
- L系统的组成要素：字母表、公理、规则
- 字符串重写与递归生成
- 图形解释：龟绘图解释器
-
生物分形应用:
- 肺部血管的分形特性
- 支气管树的递归分支结构
- 分形在优化生物系统中的优势
-
计算优化技术:
- 使用Numba加速迭代计算
- 对数缩放增强视觉细节
- 递归算法的高效实现
- 动画展示分形生成过程
-
跨学科联系:
- 数学：复数动力学
- 计算机科学：递归算法
- 生物学：分形在生命系统中的普遍性
- 医学：肺部结构与功能优化
本节提供了强大的分形几何可视化工具，从数学的Mandelbrot集到生物的肺部血管结构。分形树生长动画展示了递归构建的奇妙过程，而肺部血管案例则体现了分形在生物系统中的自然存在和功能优化。
八、结语
分形几何不仅是数学的一个分支，更是连接科学、艺术与自然的通用语言。通过本文的工具和案例，读者将学会用代码捕捉分形的无限细节，同时理解其在生物学、医学等领域的深层意义——从微观的细胞结构到宏观的生态系统，分形始终是复杂系统高效运作的底层逻辑。希望这些工具能激发更多跨学科的探索，让分形的智慧在数据可视化、算法设计乃至科学发现中绽放光彩。
