Julia集是一個在復平面上形成分形的點的集合,它最早由法國數學家Gaston Julia發現。
Julia集合可以由下式進行反復迭代得到:f(z) = z2 + c, 其中z是復平面某一點,c是一個復常數。把這個公式反復迭代,最終會得到一個復數C,然后根據C的模的大小,把這個點映射成不同的顏色,漂亮的Julia集分形就出來了。可以參閱M67的這篇文章,其中有詳細的介紹。
下面的幾幅圖是我用c++和opencv畫出來的,一張500*500的圖片,迭代15次,在我的i5電腦上跑約不到10秒,速度稍慢。圖中不同的C值對應着不同的Julia集,配色比較爛,大家湊合着看。
c = 0.73i

c = -0.3128 + 0.756i

c = 0.45

c = 0.3

當把復數C替換成該點的坐標時,著名的Mandelbrot集就出現了,圖中每一個像素點都是一個Julia集。可以把Mandlbrot集是理解為是Julia集的一個縮影,圖中不同的顏色表示該點所在的Julia集的發散速度。這張我圖畫的很不滿意,因為迭代次數太少,精度不夠,而且顏色也非常難看。在我的x86系統上迭代15次就已經有一些點會overflow了,嘆氣!期待高人指點

這里盜用M67的一段話:“生成 Mandelbrot 集的算法和生成 Julia 集的算法完全一樣,只是這一次我們固定的是初始值,而把 c 當作了變量。Mandelbrot 集內的每一個點就對應了一個連通的 Julia 集,Mandelbrot 集合外的點則對應了不連通的 Julia 集,並且很容易想到,越靠近 Mandelbrot 集的邊界,對應的 Julia 集形狀就越詭異。因此, Mandelbrot 集還有另外一種解讀方法:它就是 Julia 集的縮略圖!完全沒有比喻的意思,它真的就是 Julia 集的縮略圖”。此處的初始值就是每個像素點所在坐標。

代碼如下,c和d對應着上文中復數C的實部和虛部,改變它們就可以生成不同的Julia集
其中有兩行注釋掉的語句,把它們和上兩行替換掉,程序就可以畫出Mandelbrot集了。顏色的配置我寫的很麻煩,因為自己手工設置的有一種親切感。
#include "stdafx.h" #include <cv.h> #include <highgui.h> #include <cxcore.h> const int icount = 15; //迭代次數 const float c = -0.3128; //實部 const float d = 0.756; //虛部
double m_real, m_image; //Mandelbro集
CvScalar colortab[21]; class Complex { public: double real; double image; Complex(double r=0, double i=0){real = r, image = i;} }; Complex operator+(const Complex& a, const Complex &b) { Complex c; c.real = a.real + b.real; c.image = a.image + b.image; return c; } Complex operator*(const Complex& a, const Complex &b) { Complex c; c.real = a.real * b.real - a.image * b.image; c.image = a.image * b.real + a.real * b.image; return c; } double Model(Complex a) { return sqrtf(a.real * a.real + a.image * a.image); } double Iteration(Complex a, int n) { if(n==0) return Model(a); else { Complex temp = a*a ; temp.real += c; temp.image += d;
// temp.real += m_real; 把這兩句代替前面的兩句就是mandelbrot集了
// temp.image += m_image;
return Iteration(temp, n-1);
}
} CvScalar dye(double dist) { if(dist < 1.0/4096) return colortab[0]; else if(dist < 1.0/1024) return colortab[1]; else if(dist < 1.0/256) return colortab[2]; else if(dist < 1.0/64) return colortab[3]; else if(dist < 1.0/16) return colortab[4]; else if(dist < 1.0/4) return colortab[5]; else if(dist < 1) return colortab[6]; else if(dist < 4) return colortab[7]; else if(dist < 16) return colortab[8]; else if(dist < 64) return colortab[9]; else if(dist < 256) return colortab[10]; else if(dist < 1024) return colortab[11]; else if(dist < 4096) return colortab[12]; else if(dist < 16384) return colortab[13]; else if(dist < 65536) return colortab[14]; else if(dist < 262144) return colortab[15]; else if(dist < 1048576) return colortab[16]; else if(dist < 4194304) return colortab[17]; else if(dist < 16777216) return colortab[18]; else if(dist < 67108864) return colortab[19]; else return colortab[20]; } int _tmain(int argc, _TCHAR* argv[]) { colortab[0] = CV_RGB(28, 28, 28); colortab[1] = CV_RGB(130, 130, 130); colortab[2] = CV_RGB(85, 26, 139); colortab[3] = CV_RGB(224, 102, 255); colortab[4] = CV_RGB(255, 187, 255); colortab[5] = CV_RGB(0,0,205); colortab[6] = CV_RGB(72, 118, 255); colortab[7] = CV_RGB(0, 191, 255); colortab[8] = CV_RGB(0, 255, 255); colortab[9] = CV_RGB(0, 255, 127); colortab[10] = CV_RGB(0, 255, 0); colortab[11] = CV_RGB(50, 205, 50); colortab[12] = CV_RGB(173, 255, 47); colortab[13] = CV_RGB(255, 185, 15); colortab[14] = CV_RGB(255, 215, 0); colortab[15] = CV_RGB(255, 255, 0); colortab[16] = CV_RGB(255, 69, 0); colortab[17] = CV_RGB(255, 140, 0); colortab[18] = CV_RGB(255, 211, 155); colortab[19] = CV_RGB(255, 231, 186); colortab[20] = CV_RGB(255, 239, 213); IplImage* img = cvCreateImage(cvSize(500,500), 8, 3); for (int Y=0; Y<img->height; Y++) { for (int X=0; X<img->width; X++) { float x = (X - img->width/2) / 200.0; float y = (Y - img->height/2) / 200.0; m_real = x;
m_image = y;
Complex a(x,y); float dist = Iteration(a, icount); CvScalar color = dye(dist); cvSet2D(img, Y, X, color); } } cvNamedWindow("Julia"); cvShowImage("Julia", img); cvWaitKey(0); //保存圖片 char* path = "C:\\Users\\Administrator\\Desktop\\Julia.bmp"; cvSaveImage(path, img); cvReleaseImage(&img); }
