マンデルブロ集合の描画結果をもっと早く綺麗にみるために

前のエントリで、マンデルブロ集合の計算量の多さを隠蔽するために画像を縦に分割して描画して早く計算結果を見られるようにした。このやりかたはわかりやすいものの、もっといい方法が無いか考えた。3x3のマスを考える。左上から右下に描画していく初めのやり方はプログラムしやすい。計算結果をみるのに時間がかかる。縦に分割すると、ある行はまったく描画されない。

縦と横両方共に分けて塗ってみたらどうか。つまり、描画を2回に分けて、初めは×字のマスを塗って、次に残った菱形のマスを描画する。大まかな輪郭がこれまでより綺麗に表示され、細部は後から描画される。よさそうだ

今回は点を打つRGBの比を変えた。

実際よさげ。
描画中/描画後

#include<GL/glut.h>
#include<stdio.h>
#include<math.h> // fabs, fmod

#define INF (1e8)

//double scale=5.;
const double scale=4.0;
const int iter = 60; // この値まで繰り替えしてINFより小さいならマンデルブロ集合とみなす
const double diff = 0.00390625;// 2^-8
int w, h;

#define EPS (1e-8)
#define FLOAT_EQ(x,y) (fabs(x - y) < EPS)

void render(double additional_diff)
{	
	int i;
	double a, b;
	double ca, cb;
	double za, zb;
	double tmp;

	// 行の描画
	for(b = -scale; b<scale;  b+=diff){
		double s = -scale + (FLOAT_EQ(fmod(b, diff*2), 0) ? additional_diff: diff - additional_diff);

		for(a = s; a<scale; a+=diff*2.0){
			// c = a + biに関してテスト
			ca = a;
			cb = b;
			za = zb = 0.0;
			for(i=0; i<iter; i++){
				tmp = za*za - zb*zb + ca;
				zb = 2*za*zb + cb;
				za = tmp;
				if(fabs(za) > INF || fabs(zb) > INF){ // 発散したら
					double x = (double)i /iter;
					glColor3f(x, x/2, x*2.5);
					glVertex2d(a, b);
					break;
				}
			}
			if(i==iter){
				glColor3f(0, 0, 0); // マンデルブロ集合である
				glVertex2d(a, b);
			}
		}
	}
}

void mandelbrot()
{
	render(0);
	render(diff);
}

void display()
{
	glClear(GL_COLOR_BUFFER_BIT);

	glMatrixMode(GL_MODELVIEW);
	glLoadIdentity();
	gluPerspective(20.0, (double)w/h, 1, 100);
	glTranslated(0.5, 0, -scale*2);

	glBegin(GL_POINT);
	mandelbrot();
	glEnd();
}

void reshape(int ww, int hh)
{
	w = ww;
	h = hh;
}

int main(int argc, char **argv)
{
	glutInit(&argc, argv);
	w = h = 400;
	glutInitWindowSize(w, h);
	glutCreateWindow("mandelbrot");

	glClearColor(0,0,0,0);

	glutReshapeFunc(reshape);
	glutDisplayFunc(display);

	glutMainLoop();
	return 0;
}