每日算法 — 使用java实现分形生成:曼德勃罗集复数迭代与逃逸时间分形着色

算法背景

分形(Fractal)是具有自相似性的几何结构,在任意缩放尺度下都呈现相似的复杂图案。曼德勃罗集(Mandelbrot Set)是最著名的分形之一,由数学家Benoit Mandelbrot于1980年提出。它通过简单的复数迭代公式 z_{n+1} = z_n^2 + c 就能生成无限复杂、边界精美的图形,被誉为”数学中最复杂的对象”。

核心思想

曼德勃罗集的定义基于复数平面上每个点 c 的迭代行为。对于给定的复数 c,从 z_0 = 0 开始迭代 z_{n+1} = z_n^2 + c:

  • 如果序列 |z_n| 始终保持有界(不超过2),则 c 属于曼德勃罗集,通常染为黑色
  • 如果序列逃逸到无穷大,则 c 不属于曼德勃罗集,根据逃逸速度着色

逃逸时间算法记录 |z_n| > 2 时的迭代次数,用此值生成颜色梯度。

数学原理

复数迭代公式展开为实数形式:

z = x + yi
z^2 = (x^2 - y^2) + (2xy)i
z^2 + c = (x^2 - y^2 + cx) + (2xy + cy)i

因此迭代规则为:

x_{n+1} = x_n^2 - y_n^2 + cx
y_{n+1} = 2*x_n*y_n + cy

逃逸时间着色

逃逸迭代次数 n 与最大迭代次数 maxIter 的比值可以映射为颜色。常用方法包括:

  1. 线性映射:color = n / maxIter
  2. 对数平滑:smooth = n + 1 - log(log(|z|)) / log(2),消除颜色条纹
  3. HSV空间:将迭代次数映射为色相,饱和度固定

Java完整实现

import javax.imageio.ImageIO;
import java.awt.image.BufferedImage;
import java.io.File;
import java.io.IOException;

/**
 * 曼德勃罗集分形生成器
 * 使用逃逸时间算法与平滑着色
 */
public class MandelbrotGenerator {

    // 图像参数
    private int width;
    private int height;
    private int maxIterations;

    // 视窗范围(复平面坐标)
    private double minX, maxX, minY, maxY;

    public MandelbrotGenerator(int width, int height, int maxIterations) {
        this.width = width;
        this.height = height;
        this.maxIterations = maxIterations;
        // 默认视窗:经典曼德勃罗集范围
        this.minX = -2.5;
        this.maxX = 1.0;
        this.minY = -1.25;
        this.maxY = 1.25;
    }

    /**
     * 设置视窗范围(支持缩放浏览)
     */
    public void setView(double minX, double maxX, double minY, double maxY) {
        this.minX = minX;
        this.maxX = maxX;
        this.minY = minY;
        this.maxY = maxY;
    }

    /**
     * 计算单个像素的逃逸迭代次数(带平滑处理)
     * @param cx 复数c的实部
     * @param cy 复数c的虚部
     * @return   平滑后的逃逸值
     */
    public double computeEscape(double cx, double cy) {
        double x = 0.0, y = 0.0;
        double x2 = 0.0, y2 = 0.0;

        int iteration = 0;
        // 优化:利用 x^2 和 y^2 避免重复计算
        while (x2 + y2 <= 4.0 && iteration < maxIterations) {
            y = 2.0 * x * y + cy;
            x = x2 - y2 + cx;
            x2 = x * x;
            y2 = y * y;
            iteration++;
        }

        if (iteration == maxIterations) {
            return -1; // 未逃逸,属于集合内部
        }

        // 对数平滑:消除颜色条纹
        double logZn = Math.log(x2 + y2) / 2.0;
        double nu = Math.log(logZn / Math.log(2.0)) / Math.log(2.0);
        return iteration + 1 - nu;
    }

    /**
     * HSV到RGB颜色转换
     * @param h 色相 [0, 360)
     * @param s 饱和度 [0, 1]
     * @param v 明度 [0, 1]
     * @return  RGB整数颜色
     */
    public int hsvToRgb(double h, double s, double v) {
        double c = v * s;
        double x = c * (1 - Math.abs((h / 60.0) % 2 - 1));
        double m = v - c;

        double r, g, b;
        if (h < 60)       { r = c; g = x; b = 0; }
        else if (h < 120) { r = x; g = c; b = 0; }
        else if (h < 180) { r = 0; g = c; b = x; }
        else if (h < 240) { r = 0; g = x; b = c; }
        else if (h < 300) { r = x; g = 0; b = c; }
        else              { r = c; g = 0; b = x; }

        int ri = (int) ((r + m) * 255);
        int gi = (int) ((g + m) * 255);
        int bi = (int) ((b + m) * 255);
        return (ri << 16) | (gi << 8) | bi;
    }

    /**
     * 将逃逸值映射为颜色
     * @param escape 平滑逃逸值
     * @return       RGB颜色
     */
    public int escapeToColor(double escape) {
        if (escape < 0) {
            return 0x000000; // 集合内部:黑色
        }

        // 周期性颜色映射,产生漂亮的条纹效果
        double t = escape / 50.0; // 周期因子
        double hue = (t * 360.0) % 360.0;
        double saturation = 0.8 + 0.2 * Math.sin(t * Math.PI);
        double value = 0.9;

        return hsvToRgb(hue, saturation, value);
    }

    /**
     * 生成曼德勃罗集图像
     * @return BufferedImage对象
     */
    public BufferedImage generate() {
        BufferedImage image = new BufferedImage(width, height, BufferedImage.TYPE_INT_RGB);

        double dx = (maxX - minX) / width;
        double dy = (maxY - minY) / height;

        // 逐像素计算
        for (int py = 0; py < height; py++) {
            double cy = maxY - py * dy; // 注意:图像y轴向下
            for (int px = 0; px < width; px++) {
                double cx = minX + px * dx;
                double escape = computeEscape(cx, cy);
                int color = escapeToColor(escape);
                image.setRGB(px, py, color);
            }
        }

        return image;
    }

    /**
     * 多线程并行生成(利用现代CPU多核性能)
     */
    public BufferedImage generateParallel() {
        BufferedImage image = new BufferedImage(width, height, BufferedImage.TYPE_INT_RGB);
        double dx = (maxX - minX) / width;
        double dy = (maxY - minY) / height;

        int numThreads = Runtime.getRuntime().availableProcessors();
        Thread[] threads = new Thread[numThreads];

        for (int t = 0; t < numThreads; t++) {
            final int threadId = t;
            threads[t] = new Thread(() -> {
                for (int py = threadId; py < height; py += numThreads) {
                    double cy = maxY - py * dy;
                    for (int px = 0; px < width; px++) {
                        double cx = minX + px * dx;
                        double escape = computeEscape(cx, cy);
                        image.setRGB(px, py, escapeToColor(escape));
                    }
                }
            });
            threads[t].start();
        }

        for (Thread thread : threads) {
            try { thread.join(); } catch (InterruptedException e) { e.printStackTrace(); }
        }

        return image;
    }

    /**
     * 放大到指定坐标,探索分形细节
     */
    public void zoomTo(double centerX, double centerY, double zoomFactor) {
        double rangeX = (maxX - minX) / zoomFactor;
        double rangeY = (maxY - minY) / zoomFactor;
        setView(centerX - rangeX / 2, centerX + rangeX / 2,
                centerY - rangeY / 2, centerY + rangeY / 2);
    }

    public static void main(String[] args) throws IOException {
        // 参数设置
        int width = 800;
        int height = 800;
        int maxIter = 500;

        MandelbrotGenerator generator = new MandelbrotGenerator(width, height, maxIter);

        System.out.println("=== 生成经典曼德勃罗集 ===");
        long start = System.currentTimeMillis();
        BufferedImage image = generator.generateParallel();
        long end = System.currentTimeMillis();
        System.out.println("生成耗时: " + (end - start) + "ms");

        // 保存图像
        ImageIO.write(image, "png", new File("mandelbrot.png"));
        System.out.println("图像已保存至 mandelbrot.png");

        // 放大到经典"海马谷"区域
        System.out.println("\n=== 放大到海马谷区域 ===");
        generator.zoomTo(-0.745, 0.105, 100.0);
        generator.maxIterations = 2000; // 高放大需要更多迭代
        start = System.currentTimeMillis();
        BufferedImage zoomed = generator.generateParallel();
        end = System.currentTimeMillis();
        System.out.println("放大生成耗时: " + (end - start) + "ms");
        ImageIO.write(zoomed, "png", new File("mandelbrot_zoom.png"));
        System.out.println("放大图像已保存至 mandelbrot_zoom.png");

        // 性能统计
        System.out.println("\n=== 性能统计 ===");
        System.out.println("总像素: " + (width * height));
        System.out.println("最大迭代: " + maxIter);
        System.out.println("平均每像素耗时: " +
            String.format("%.4f", (end - start) * 1000000.0 / (width * height)) + " ns");
    }
}

复杂度分析

指标 复杂度 说明
时间复杂度 O(W * H * M) W=宽度, H=高度, M=最大迭代次数
空间复杂度 O(W * H) 存储输出图像
迭代优化 O(M) 每次 利用x²和y²缓存避免重复计算

优化技巧

  1. 周期检测:若 z_n = z_m 则进入循环,可提前终止
  2. 边界检查:利用主心形和周期2 bulb 的解析边界快速判断内部点
  3. SIMD向量化:利用AVX/SSE指令同时计算多个像素
  4. GPU加速:OpenCL/CUDA实现每秒数十亿次迭代

其他经典分形

Julia集

Julia集与曼德勃罗集密切相关,固定复数 c,变化初始 z_0 为像素坐标:

public double computeJulia(double zx, double zy, double cx, double cy) {
    double x = zx, y = zy;
    for (int i = 0; i < maxIterations; i++) {
        double xNew = x * x - y * y + cx;
        y = 2 * x * y + cy;
        x = xNew;
        if (x * x + y * y > 4.0) return i;
    }
    return -1;
}

牛顿分形

基于牛顿迭代法求根,用不同根的区域着色:

z_{n+1} = z_n - f(z_n) / f'(z_n)

对于 f(z) = z^3 - 1,可得到精美的三色分形图案。

扩展应用

  1. 地形生成:Perlin噪声结合分形维度生成自然地形
  2. 纹理合成:分形布朗运动生成木纹、云层纹理
  3. 网络拓扑:分形网络模型研究互联网拓扑结构
  4. 金融市场:分形分析股票价格波动特征

总结

曼德勃罗集展示了简单规则如何产生无限复杂的美。逃逸时间算法是该分形的核心渲染技术,对数平滑着色能显著提升视觉效果。通过并行计算,现代CPU可在秒级生成高清分形图像。理解分形生成原理,是探索混沌理论与复杂系统的绝佳入口。