September 18, 2011

繪製 Mandelbrot Set 的小程式

工作的空檔,偶然想起 fractals (碎形),為紀念去年因胰臟癌過世的 Benoît Mandelbrot 大師 (1924-2010),就嘗試撰寫繪製 Mandelbrot Set 的小程式,藉以體驗自我相似結構的美妙。這個程式採用最單純的圖像格式 [PPM],基本上把寬度與高度的資訊描述好,就是逐一填入像素資料。程式碼如下:(mandelbrot.c)
#include <stdio.h>
#include <stdlib.h>
#include <complex.h>
 
#define width_size      800
#define height_size     600
#define Maxval          255

static const float orig_x = width_size * 2/3;
static const float orig_y = height_size * 1/2;

typedef struct _pixel {
    unsigned char r;
    unsigned char g;
    unsigned char b;
} pixel;

static const pixel dim_gray = { 105, 105, 105 };

static unsigned char iteration(int x, int y)
{
    const int limit = Maxval + 1;
    int i;
    complex c = ((x - orig_x) / (width_size / 3)) +
                ((orig_y - y) / (height_size / 2)) * I;
    complex z = 0;

    for (i = 0; i < limit; i++) {
        /* basic formula */
        z = z * z + c;
        if (creal(z) > 2 || cimag(z) > 2)
            break;
    }
    return (unsigned char) (i == limit ? 0 : i);
}


int main()
{
    FILE *f = fopen("mandelbrot.ppm", "w+");

    /* PPM header */
    fprintf(f,
            "P6\n"      /* PPM magic number */
            "#Mandelbrot Set\n"
            "%d "       /* width, in ASCII decimal */
            "%d\n"      /* height, in ASCII decimal */
            "%d\n",     /* maximum color value, in ASCII decimal */
            width_size, height_size, Maxval);

    /* Write every pixel generated by Mandelbrot Set */
    for (int i = 0; i < height_size; i++) {
        for (int j = 0; j < width_size; j++) {
            unsigned char iter = iteration(j, i);
            if (iter) {
                pixel p = {
                    .r = iter,
                    .g = (float) abs(j - orig_x) / width_size * Maxval,
                    .b = (float) abs(i - orig_y) / height_size * Maxval };
                fwrite(&p, sizeof(pixel), 1, f);
            } else {
                fwrite(&dim_gray, sizeof(pixel), 1, f);
            }
        }
    }

    fclose(f);
    return 0;
}
程式利用 C99 提供的複數型態 (complex),很直覺地演繹 Mandelbrot set 的基本公式,相關的原理與繪製技巧可見 [How to Plot the Mandelbrot Set By Hand],編譯並執行本程式:
$ gcc -Wall -std=c99 mandelbrot.c -o mandelbrot && ./mandelbrot
當程式執行完畢,會輸出 "mandelbrot.ppm",可使用 GIMP 或桌面環境內建的瀏覽工具來觀賞,畫面如下:

程式以反覆代入方程式的方式,計算結果來得到上述圖形,並加以彩繪,基本上就是控制像素的 RGB 比例。關於背景知識,可參考非常詳盡的中文參考資料 [碎形 Fractal]。
由 jserv 發表於 September 18, 2011 11:53 PM
迴響
發表迴響









記住我的資訊?