c++畫分形之Julia集與Mandelbrot集


  Julia集是一個在復平面上形成分形的點的集合,它最早由法國數學家Gaston Julia發現。

  Julia集合可以由下式進行反復迭代得到:f(z) = z+ 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); }

 

 

 


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM