メインコンテンツへスキップ
見出し画像

マンデルブロ集合×反復関数系

    どうも、108Hassiumです。

    先日(2ヶ月くらい前)、こんなフラクタル図形を発見しました。

    画像

    形はマンデルブロ集合に似ていますが、細部は全然違った見た目になっています。

    この記事ではこのフラクタル図形の描画方法や関連するほかの図形を紹介します。

    反復関数系

    フラクタル図形の有名な描画方法として、反復関数系(IFS:Itrated Function System)というものがあります。

    まず、簡単な例を1つ紹介します。

    最初に、以下のように関数を3つ用意します。

    • f1(z)=z2f_1(z)=\frac{z}{2}

    • f2(z)=z+12f_2(z)=\frac{z+1}{2}

    • f3(z)=z+i2f_3(z)=\frac{z+i}{2}

    次に、初期値として適当な複素数を1つ選びます。

    今回は0を初期値とします。

    そして、以下のような手順で複素平面上に点をプロットしていきます。

    • 初期値をプロットする。

    • ある点α\alphaをプロットしたら、f1(α)f_1(\alpha)、f2(α)f_2(\alpha)、f3(α)f_3(\alpha)もプロットする

    今回の場合だと、
    00、
    f1(0)f_1(0)、f2(0)f_2(0)、f3(0)f_3(0)、
    f1(f1(0))f_1(f_1(0))、f2(f1(0))f_2(f_1(0))、f3(f1(0))f_3(f_1(0))、f1(f2(0))f_1(f_2(0))、f2(f2(0))f_2(f_2(0))、f3(f2(0))f_3(f_2(0))、f1(f3(0))f_1(f_3(0))、f2(f3(0))f_2(f_3(0))、f3(f3(0))f_3(f_3(0))・・・
    というように点が増えていきます。

    実際にコンピュータを使って描画するときは、このプロセスを適当なところで打ち切って有限時間の計算で終わらせます。

    Processingというプログラミング言語を使うと、以下のように実装できます。

    void setup(){
      size(2000,2000);
      background(0);
      stroke(255);
      f(0,0,15);
    }
    
    void f(float x,float y,int n){
      point(x*2000,y*2000);
      if(0<n){
        f(x/2,y/2,n-1);
        f(x/2+0.5,y/2,n-1);
        f(x/2,y/2+0.5,n-1);
      }
    }

    このコードを実行できるソフトウェアは、以下のページからダウンロードできます。

    実行するとこのような画像が表示されます。

    画像

    ※コードの意味や画像の保存方法は解説しないので知りたい人はググってください。

    こんな感じで、反復関数系によるフラクタル図形は関数の計算を繰り返すことで生成されます。

    先程のサンプルでは1次関数を使いましたが、他にもいろいろな関数を使ってフラクタル図形を描画することができます。

    画像
    ☝z^2-0.2+0.2i、z^2-0.2-0.2i
    画像
    ☝z^2-0.1+0.25i、z^2-0.1-0.25i
    画像
    ☝z^2+0.3+0.4i、z^2+0.3-0.4i
    画像
    ☝z^2-0.5、z^2+0.2+0.2i


    一方、マンデルブロ集合の描画にはescape time algorithmという方法が使われています。

    ざっくり説明するとこんな感じです。

    1. 複素平面上から1点ccを選ぶ。

    2. z0=0z_0=0、zn+1=zn2+cz_{n+1}=z_n^2+cという数列を、znz_nの絶対値が閾値(2以上の定数)を超えるまで計算する。閾値を超えない場合、適当な大きさのnnで計算を止める。

    3. 発散するまでの計算回数に応じてccを彩色する。

    void setup(){
      size(2000,2000);
      background(0);
      noStroke();
      double x,y,px,py,cx,cy;
      boolean o;
      for(int a=0;a<2000;a++){
        for(int b=0;b<2000;b++){
          cx=(double)a/500.0-2.0;
          cy=(double)b/500.0-2.0;
          x=0;
          y=0;
          o=true;
          for(int c=0;c<500&&o;c++){
            px=x;
            py=y;
            x=px*px-py*py+cx;
            y=2.0*px*py+cy;
            if(x*x+y*y>4){
              o=false;
              fill(cr(c*7),cr(c*8),cr(c*9));
              rect(a,b,1,1);
            }
          }
        }
      }
    }
    
    int cr(int c){
      return (c%256)*(255-(c%256))/65;
    }
    画像

    三重for文の2層目を以下のように書き換えることで、ジュリア集合も描画できます。

          ...
        for(int b=0;b<2000;b++){
          cx=0.3;
          cy=0.5;
          x=(double)a/500.0-2.0;
          y=(double)b/500.0-2.0;
          o=true;
          for(int c=0;c<500&&o;c++){
            ...

    マンデルブロ集合×IFS

    さて冒頭の画像ですが、大体以下のような計算をしています。

    1. 複素平面上から1点ccを選ぶ。

    2. f1(z)=z2+cf_1(z)=z^2+c、f2(z)=2zf_2(z)=2z、z0=0z_0=0という組み合わせでIFSの計算をして、発散した回数(=絶対値が閾値を超えた回数)を求める。

    3. 発散回数に応じてccを彩色する。

    Processingのコードは以下のようになります。

    void setup(){
      size(2000,2000);
      background(0);
      noStroke();
      float n=0;
      for(int a=0;a<2000;a++){
        for(int b=0;b<2000;b++){
          n=f(0,0,(float)a/500.0-2.0,(float)b/500.0-2.0,25);
          fill(cr(sqrt(n)*1.0),cr(sqrt(n)*2.0),cr(sqrt(n)*3.0));
          rect(a,b,1,1);
        }
      }
    }
    
    int f(float x,float y,float cx,float cy,int n){
      if(n==0){
        return 0;
      }else if(x*x+y*y>4){
        return 1;
      }else{
        return f(x*x-y*y+cx,2*x*y+cy,cx,cy,n-1)+f(2*x,2*y,cx,cy,n-1);
      }
    }
    
    float cr(float c){
      return (c%256)*(255-(c%256))/65;
    }

    ちなみに、実は冒頭の画像に似たものを既にTwitter上に投稿しているのですが、コードは若干異なりfの定義が以下のようになっています。

    int f(float x,float y,float cx,float cy,int n){
      if(n==0){
        if(x*x+y*y>100){
          return 0;
        }else{
          return 1;
        }
      }else{
        if(x*x+y*y>100){
          return 0;
        }else{
          return f(x*x-y*y+cx,2*x*y+cy,cx,cy,n-1)+f(2*x,2*y,cx,cy,n-1);
        }
      }
    }

    先程のものとは異なり、発散した回数ではなく収束した回数を数えています。

    生成される画像は、模様が若干変わる(主に原点付近)だけで細部の構造は大体同じです。

    画像

    ギャラリー

    以下、関連画像集です。

    画像

    IFS版マンデルブロ集合の、-0.16+1.03i付近の拡大図です。

    普通のマンデルブロ集合だと飛び地のある場所ですが、この画像でもそれっぽい形が見えます。

    画像
    ☝z^2+0.25,2z

    先程のコードの、

    n=f(0,0,(float)a/500.0-2.0,(float)b/500.0-2.0,25);

    この部分を、

    n=f((float)a/500.0-2.0,(float)b/500.0-2.0,0.25,0,25);

    このように書き換えることでISF版ジュリア集合を描画できます。

    画像
    ☝z^2+0.3+0.5i,2z
    画像
    ☝z^2-0.1+0.9i,2z
    画像
    ☝z^2-0.75,2z
    画像
    ☝z^2-1,2z
    画像
    ☝z^2-1.75,2z

    IFS版ジュリア集合は普通のジュリア集合と似た形になるようですが、枝分かれ等の細かい特徴は反映されないこともあるようです。

    画像
    ☝z^3+c,2z
    画像
    ☝z^3+0.4,2z
    画像
    ☝z^3+i,2z

    z3+cz^3+cも、普通のマンデルブロ集合・ジュリア集合と同じように描画できるようです。

    なお、どんな関数でもescape time algorithmと同様に描画できるわけではなく、例えばcz2−1+1\frac{c}{z^2-1}+1とかは2z2zとの組み合わせだと上手くいかないようです。

    画像
    ☝z^2+c,4z
    画像
    ☝z^2+c,4z(拡大)
    画像
    ☝z^3+c,3z
    画像
    ☝z^3+c,3z(拡大)

    1次関数の係数を大きくすると、何故か内側と外側の2層に分離するようです。

    画像
    ☝z^2+0.3+0.5i,4z
    画像
    ☝z^2-0.1+0.7i,4z

    z2+c,2zz^2+c,2zの場合と違い、z2+c,4zz^2+c,4zの組み合わせだとジュリア集合は通常の形に近くなるようです。

    画像
    ☝con(z)^2+c,2z

    トリコーンのIFS版です。

    元のトリコーンでは何もなかった中央部分の周りにもギザギザした模様ができています。

    画像
    ☝con(z)^2-1.0+0.1i,2z
    画像
    ☝(Abs(x)+iAbs(y))^2+c,2z

    バーニングシップ・フラクタルのIFSバージョンです。

    画像
    ☝z^2+c,z^2

    「マンデルブロ集合とIFSを組み合わせる」という事を思いついてから、最初に試したパターンがこれです。

    z2+c,2zz^2+c,2zとと比べるとあまりマンデルブロ集合っぽくないです。

    画像
    z^2+1,z^2
    画像
    ☝z^2-2,z^2
     
     
    数学関係の記事を書きます。毎週日曜更新予定。

    あなたへのおすすめ