マンデルブロ集合を描く[Processing]

　この間フラクタル図形を描いてみたけれど，一番の有名どころを描いてなかったのでやってみることにした。一応述べておくが私はフラクタルが好きというわけではないです。ただなんか無性に綺麗な図形を描きたくなる，いわば「Processingハイ」という状態であるというだけです。

1.マンデルブロ集合について

詳細はググるなり本で調べるなどしてください。ここでは図形を描くための最低限の情報について書き留めておきます。手元にある本も見ていますが，だいたいWikipediaにおんなじこと書いてある（Wikipediaさんすげー）．
　まずマンデルブロ集合というのは次の漸化式で表される数列Znが収束するcの集合です。ただしZn,cは複素数。

　このときcを

　とあらわすならば，点(a,b)について打点していけばマンデルブロ集合を描ける。
　しかし，Processingでは複素数は標準ではサポートしていない。自分でComplex型を定義するのも難儀なので複素数を実部と虚部に分けて計算する。
　Znを次のように実部と虚部に分ける。

そうするとXn，Ynは次のように表される。

　以下に示すプログラムでは，このようにして計算していく。ただしZnの大きさの二乗が閾値Lを超えたら計算をやめるようにしてる。そのへんはプログラムを参考にしてほしいです。

2.ソースコード

まぁ長々と説明してきたがつまりこういうこと。

/* Znを何回漸化式で計算するか */
int N = 255;
/* Znの最大値 */
int L = 255;
/* 拡大率 */
float SCALE = 3.0;

void setup() {
  size(512, 512);
}

void draw() {
  /*中心に原点を持ってくる*/
  translate(width/2, height/2);
  background(0);

  for (int a = -width/2; a <= width/2; a++ ) {
    for (int b = -height/2; b <= height/2; b++ ) {
      /* 複素数c=x+yiを定める */
      float x = SCALE * a / width;
      float y = SCALE * b / height;
      /* 漸化式に基づいて計算 */
      int r = calc(x, y);

      /* 収束の速度に応じて色を指定 */
      /* 下のプログラムは適当。stroke(赤成分緑成分,青成分) */
      stroke(r%256, r*4%256, r*16%256);
      /* 点の描画 */
      rect(a, b, 1, 1);
    }
  }
}

/* 漸化式にもとづいてZnが発散するnを計算 */
int calc(float x, float y) {
  float tx, ty;
  float zx = 0.0;
  float zy = 0.0;

  for (int i = 1; i <= N; i++ ) {
    tx = zx;
    ty = zy;

    zx = tx*tx-ty*ty+x;
    zy = 2*tx*ty+y;

    /* 発散した（閾値を超えた）場合はiを返す */
    if ( zx*zx + zy*zy > L )
      return i;
  }
  
  /* 発散しなかった場合は0を返す */
  return 0;
}

　割りといい感じに描けた。