这是本人之前用C语言代码写的Mandelbrot Set的复平面分形的可视化实现。尽管这次的作业并不难，但是我觉得仍然值得记录。
叠甲：本人不才，代码水平有限，只是单纯的记录。

项目介绍

项目目的

用C语言实现mandelbrot集，将其通过静态图片、动态视频等进行可视化展现，使用多线程并行计算加速。

实现功能

  1. mandelbrot集的总体静态图片展现
  2. mandelbrot集的局部动态视频展现
  3. 多线程并行计算加速

使用方法

  1. make video或者make后输入./mandelbrot来编译运行程序，可生成mandelbrot集bmp静态图片和mp4动态视频
  2. make doc生成项目报告pdf文档
  3. make clean-all清除除了动态视频和说明文档外所有生成的文件
  4. make clean-rest清除动态视频和说明文档

依赖关系

  1. ffmpeg用于生成动态视频
  2. math库
  3. OpenMP库

代码实现

.h文件

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
# include <stdio.h>
# include <stdlib.h>
# include <math.h>
# include <complex.h>
# include <ctype.h>
# include <unistd.h>
/**
 * @brief 结构体mandelbrot
 * @param x_min 定义区域的x轴最小值，复平面左边界
 * @param x_max 定义区域的x轴最大值，复平面右边界
 * @param y_min 定义区域的y轴最小值，复平面下边界
 * @param y_max 定义区域的y轴最大值，复平面上边界
 * @param max_iter 最大迭代次数
 * @param width 输出图像的宽度或列数
 * @param height 输出图像的高度或行数
 */
typedef struct
{
    double x_min;
    double x_max;
    double y_min;
    double y_max;
    int max_iter;
    int width;
    int height;
    int **iter_count;
}mandelbrot;
/**
 * @brief 初始化mandelbrot结构体的函数
 */
void mandelbrot_init(mandelbrot *m, double x_min, double x_max, double y_min, double y_max, int max_iter, int width, int height);
/**
 * @brief 生成并保存mandelbrot图像为bmp文件的函数
 * @param m mandelbrot结构体指针
 * @param filename 输出图像的文件名
 */
void mandelbrot_save_bmp(const mandelbrot *m,const char *filename);

/**
 * @brief 释放mandelbrot结构体的内存的函数
 */
void mandelbrot_free(mandelbrot *m);

.c文件

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
# include <stdio.h>
# include <stdlib.h>
# include <math.h>
# include <complex.h>
# include "mandelbrot.h"
# define M_PI 3.14159265358979323846
/**
 * @brief 文件头结构体
 * @param bftype 位图的类型
 * @param bfsize 位图的大小
 * @param bfreserved1 保留位1，必须为0
 * @param bfreserved2 保留位2，必须为0
 * @param bfoffbits 位图的数据偏移量
 */
#pragma pack(push, 1)
struct BMPfileheader
{
    unsigned short bftype;
    unsigned int bfsize;
    unsigned short bfreserved1;
    unsigned short bfreserved2;
    unsigned int bfoffbits;
};
/**
 * @brief 信息头结构体
 * @param bisize 信息头的大小
 * @param biwidth 位图的宽度
 * @param biheight 位图的高度
 * @param bipplanes 颜色平面数
 * @param bitcount 位深度，每个像素占用的位数 
 * @param bicompression 压缩类型
 * @param bisizeimage 位图数据大小
 * @param bixpelspermeter 水平分辨率
 * @param biypelspermeter 垂直分辨率
 * @param bicolor_used 实际使用的颜色数
 * @param bicolor_important 重要颜色数
 */
struct BMPinfoheader
{
    unsigned int bisize;
    int biwidth;
    int biheight;
    unsigned short bipplanes;
    unsigned short bitcount;
    unsigned int bicompression;
    unsigned int bisizeimage;
    int bixpelspermeter;
    int biypelspermeter;
    unsigned int bicolor_used;
    unsigned int bicolor_important;
};
#pragma pack(pop)

void mandelbrot_init(mandelbrot *m, double x_min, double x_max, double y_min, double y_max, int max_iter, int width, int height)
{
    m->x_min = x_min;
    m->x_max = x_max;
    m->y_min = y_min;
    m->y_max = y_max;
    m->max_iter = max_iter;
    m->width = width;
    m->height = height;
    m->iter_count = (int **)malloc(sizeof(int*) * m->width); //分配行内存
    for(int i = 0; i < m->width; i++)
    {
        m->iter_count[i] = (int*)malloc(sizeof(int) * m->height); //分配列内存
    }
    double dx = (m->x_max - m->x_min) / m->width; //迭代的水平变化量，即列宽度
    double dy = (m->y_max - m->y_min) / m->height; //迭代的垂直变化量，即行高度
    #pragma omp parallel for schedule(dynamic) //OpenMP多线程并行计算
    for(int i = 0; i < m->width; i++)
    {
        for(int j = 0; j < m->height; j++)
        {
            complex z = 0 + 0 * I; //初始化Z为0
            complex c = m->x_min + i * dx + (m->y_min + j * dy) * I; //计算当前点对应的复数值
            m->iter_count[i][j] = m->max_iter; //初始化迭代次数
            for(int iter_num = 0; iter_num < m->max_iter; iter_num++)
            {
                z = z * z + c; //迭代过程
                if(cabs(z) > 2.0) 
                {
                m->iter_count[i][j] = iter_num; //如果发散，更新迭代次数
                break; //退出迭代
                }
            }
        }
    }
}

void mandelbrot_save_bmp(const mandelbrot *m,const char *filename)
{
    FILE *fp = fopen(filename,"wb");
    if(fp == NULL) return;  
    int width_size = m->width * 3; //每个像素点需要3个字节rgb
    int filler = (4 - (width_size % 4)) % 4; //填充字节数，需要是4的倍数
    width_size += filler; 
    #pragma pack(push, 1)
    struct BMPfileheader fileheader;
    struct BMPinfoheader infoheader;
    fileheader.bftype = 0x4d42; //BMP文件的文件标识
    fileheader.bfsize = 54 + width_size * m->height; //14 + 40 +  像素点大小和
    fileheader.bfreserved1 = 0; //必须设为0
    fileheader.bfreserved2 = 0; //必须设为0
    fileheader.bfoffbits = 54; //偏移量，14 + 40
    infoheader.bisize = 40; //信息头的大小
    infoheader.biwidth = m->width; 
    infoheader.biheight = m->height; 
    infoheader.bipplanes = 1; //固定为1
    infoheader.bitcount = 24; //常见值24
    infoheader.bicompression = 0; //无压缩
    infoheader.bisizeimage = width_size * m->height; //像素点大小和
    infoheader.bixpelspermeter = 0;
    infoheader.biypelspermeter = 0;
    infoheader.bicolor_used = 0; //默认
    infoheader.bicolor_important = 0; //颜色同等重要
    #pragma pack(pop)
    //写入文件头和信息头
    fwrite(&fileheader, sizeof(fileheader), 1, fp);
    fwrite(&infoheader, sizeof(infoheader), 1, fp);
    for(int j = m->height - 1; j >= 0; j--) {
        for(int i = 0; i < m->width; i++) {
            int iter = m->iter_count[i][j];
            //通过迭代次数计算颜色分配
            double ratio = (iter == m->max_iter) ? 0.0 : (double)iter / m->max_iter;
            ratio = pow(ratio, 0.3); 
            unsigned char red = 255 * 0.5 * (sin(2.0 * M_PI * ratio)) + 128;
            unsigned char green = 255 * 0.5 * (sin(2.0 * M_PI * ratio + 2.0 * M_PI / 3)) + 128;
            unsigned char blue = 255 * 0.5 * (sin(2.0 * M_PI * ratio + 4.0 * M_PI / 3)) + 128;
            //写入像素点的颜色数据，BRG格式
            fwrite(&blue, 1, 1, fp);
            fwrite(&green, 1, 1, fp);
            fwrite(&red, 1, 1, fp);
        } 
        for(int k = 0; k < filler; k++) fputc(0, fp);     
    }
    fclose(fp);
}

void mandelbrot_free(mandelbrot *m)
{
    for(int i = 0; i <m->width; i++) {free(m->iter_count[i]);}
    free(m->iter_count);
}

测试文件

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
#include <stdio.h>
#include "mandelbrot.h"

int main() {
    mandelbrot m;
    int max_iter = 1000;     // 最大迭代次数
    int width = 800;      // 图像宽度  
    int height = 600;       // 图像高度
    int frames = 90;      // 帧数
    double center_x = -0.743643887;    //缩放中心点横坐标
    double center_y = 0.131825911;     //缩放中心点纵坐标
    double zoom = 0.95;    //缩放倍数
    double x0_min = -1.5;
    double x0_max = 1.5;
    double y0_min = -1.0;
    double y0_max = 1.0;
        // 初始化mandelcrot集
        mandelbrot_init(&m, -1.5, 1.5, -1.0, 1.0, max_iter, width, height);
        // 保存mandelbrot图像为BMP文件
        const char *filename = "mandelbrot.bmp";
        mandelbrot_save_bmp(&m, filename);
        printf("已将mandelbrot集图像保存为 %s\n", filename);
        // 释放先前分配的行列内存
        mandelbrot_free(&m);
        // 生成mandelbrot集的动画
        for(int i = 0; i < frames; i++)
        {
            int current_max_iter = max_iter + i * 100; 
            double decay = pow(zoom, i);  // 缩放
            double x_min = center_x - (center_x - x0_min) * decay;
            double x_max = center_x + (x0_max - center_x) * decay;
            double y_min = center_y - (center_y - y0_min) * decay;
            double y_max = center_y + (y0_max - center_y) * decay;    
            mandelbrot_init(&m, x_min, x_max, y_min, y_max, current_max_iter, width, height);
            char filename[100];
            sprintf(filename, "mandelbrot%03d.bmp", i);
            mandelbrot_save_bmp(&m, filename);
            mandelbrot_free(&m);
        }
    return 0;
}

运行结果

Makefile和report.tex就不放了

静态图片

mandelbrot.bmp

动态视频