遺伝的アルゴリズムにより三次元関数の最小値を探索してみた

変更履歴

• 240822 多峰性関数であるRastrigin function の探索結果を図5に追加しました

遺伝的アルゴリズムにより三次元関数の最小値を探索してみました。
使用したベンチマーク関数は粒子群最適化の場合と同じ下記の４つです。

• Sphere function
• Five well potential function
• Himmelblau's function
• Rosenbrock function

掲載した４本のプログラムは関数部分（function関数およびset_function_conditions関数）が違うだけで、他は全て同じです。
図１～図４に各関数での探索結果を示しました。

探索が進むと次第に親の座標が一点に集中してくるので、二つの親を選んでその中間ポイントを子供として選んでも、親と同じ座標になってしまいます。そこでプログラム中のこの部分では、わざとノイズ成分としてDXとDYを付加して、親と同じ座標になることを回避しています。

double DX = (X_MAX - X_MIN) * B_GA * (get_random() - 0.5);
double DY = (Y_MAX - Y_MIN) * B_GA * (get_random() - 0.5);
double next_x = get_random() * (max_cx - min_cx) + min_cx + DX;
double next_y = get_random() * (max_cy - min_cy) + min_cy + DY;

遺伝的アルゴリズムは急速に収束しているように見えますが、見えないところで遺伝子を大量に淘汰しているので、単純に探索効率の比較はできないと思います。ただ粒子群最適化や焼きなまし法に比べて、探索性能は単峰性関数においても多峰性関数においても、悪くないように感じます。

出力結果については全部貼ると膨大な量になってしまうので、Sphere function の最初の３回分だけを貼りました。

アルゴリズムを参考にさせていただいたサイトは既に無くなってしまったようで、残念ながらリンクを貼ることはできませんでした。

参考にさせていただベンチマーク関数のサイト

• 最適化アルゴリズムを評価するベンチマーク関数まとめ

• 焼きなまし法と粒子群最適化法の探索挙動を可視化してみた

参考にさせていただいたグラフ表示に関するサイト（実際に参考にしたのは「Visual C++の易しい使い方（23） ―三次元グラフの陰線処理―」の方ですが、今は公開されていないようです）

• 隠線処理を用いた三次元グラフの作成

誤りの指摘や感想等ありましたら、書き込んでいただけると幸いです。

図1　Sphere function の探索結果

図2　Five-well potential function の探索結果

図3　Himmelblau's function の探索結果

図4　Rosenbrock function の探索結果

図5 Rastrigin function の探索結果

sphere_ga.cpp

#define _USE_MATH_DEFINES
#include <math.h>
#include <iostream>
#include <fstream>
using namespace std;

const int ITERATION = 1000; // 試行回数

const int SWARM_SIZE = 20; // 群の大きさ

double X_MAX = 10; // graphのX軸の表示範囲の上限
double X_MIN = -10; // graphのX軸の表示範囲の下限

double Y_MAX = 10; // giterationaphのY軸の表示範囲の上限
double Y_MIN = -10; // graphのY軸の表示範囲の下限

double Z_MAX = 200; //  graphのZ軸の表示範囲の上限
double Z_MIN = 0; // graphのZ軸の表示範囲の下限

string str_formula = "x * x + y * y"; // グラフに表示する数式

double x[ITERATION][SWARM_SIZE];
double y[ITERATION][SWARM_SIZE];
double z[ITERATION][SWARM_SIZE];

double function(double x, double y)
{
	double z;

	z = x * x + y * y;

	return z;
}

void set_function_conditions()
{
	X_MAX = 10;
	X_MIN = -10;
	Y_MAX = 10;
	Y_MIN = -10;
	Z_MAX = 200;
	Z_MIN = 0;
	str_formula = "x * x + y * y";
}

int N_P = SWARM_SIZE * 2; // number of population

int CUT_LINE = SWARM_SIZE;
double A_GA = 0.3;
double B_GA = 0.1;
double MUTATION_ODDS = 0.1;
int change_counter = 0;

double get_random()
{
	return (double)rand() / ((double)RAND_MAX + 1);
}

double max1(double a, double b)
{
	double c = a;
	if (b > a)
		c = b;
	return c;
}

double min1(double a, double b)
{
	double c = a;
	if (b < a)
		c = b;
	return c;
}

class POPULATION
{
protected:
	double x, y, z;
	int* gene;
	long long weight_p;
	long long value_p;
public:
	POPULATION();
	~POPULATION();
	void calc_weight_and_value();
	long long get_value() { return value_p; }
	long long get_weight() { return weight_p; }
	void set_gene(int i, int g);
	int get_gene(int i) { return gene[i]; }
	void mutation();
	void initialize();
	void set_xyz(double x1, double y1, double z1) { x = x1; y = y1; z = z1; }
	double get_x() { return x; }
	double get_y() { return y; }
	double get_z() { return z; }
};

POPULATION** population = nullptr;

POPULATION::POPULATION()
{

}

POPULATION::~POPULATION()
{

}

void POPULATION::initialize()
{
	double x1 = get_random() * (X_MAX - X_MIN) + X_MIN;
	double y1 = get_random() * (Y_MAX - Y_MIN) + Y_MIN;
	double z1 = function(x1, y1);
	set_xyz(x1, y1, z1);
}

void POPULATION::mutation()
{
	initialize();
}

void create_object()
{
	population = new POPULATION*[N_P];
	for (int i = 0; i < N_P; i++)
		population[i] = new POPULATION;
}

void initialize_population()
{
	for (int i = 0; i < N_P; i++) 
		population[i]->initialize();
}

void ranking()
{
	for (int i = 0; i < N_P - 1; i++) {
		for (int j = i + 1; j < N_P; j++) {
			if (population[i]->get_z() > population[j]->get_z()) {
				POPULATION* tmp = population[i];
				population[i] = population[j];
				population[j] = tmp;
			}
		}
	}
}

void cross_over(int i)
{
	for (int j = CUT_LINE; j < N_P; j++) {
		int parent1 = rand() % CUT_LINE;
		int parent2 = rand() % CUT_LINE;
	
		double x1 = population[parent1]->get_x();
		double y1 = population[parent1]->get_y();
		double x2 = population[parent2]->get_x();
		double y2 = population[parent2]->get_y();
		double dx = abs(x1 - x2);
		double dy = abs(y1 - y2);
		double min_x = min(x1, x2);
		double min_y = min(y1, y2);
		double max_x = max(x1, x2);
		double max_y = max(y1, y2);
		double min_cx = min_x - A_GA * dx;
		double max_cx = max_x + A_GA * dx;
		double min_cy = min_y - A_GA * dy;
		double max_cy = max_y + A_GA * dy;
//		double next_x = get_random() * (max_cx - min_cx) + min_cx;
//		double next_y = get_random() * (max_cy - min_cy) + min_cy;
		double DX = (X_MAX - X_MIN) * B_GA * (get_random() - 0.5);
		double DY = (Y_MAX - Y_MIN) * B_GA * (get_random() - 0.5);
		double next_x = get_random() * (max_cx - min_cx) + min_cx + DX;
		double next_y = get_random() * (max_cy - min_cy) + min_cy + DY;
		next_x = max1(min1(next_x, X_MAX), X_MIN);
		next_y = max1(min1(next_y, Y_MAX), Y_MIN);
		double next_z = function(next_x, next_y);
		population[j]->set_xyz(next_x, next_y, next_z);

		if (get_random() < MUTATION_ODDS)
			population[j]->mutation();
	}
}

void init_genetic_algorithm()
{
	create_object();
	initialize_population();
	ranking();
	change_counter = 0;
	for (int j = 0; j < CUT_LINE; j++) {
		x[change_counter][j] = population[j]->get_x();
		y[change_counter][j] = population[j]->get_y();
		z[change_counter][j] = population[j]->get_z();
	}
	change_counter++;		
}

void delete_genetic_algorithm()
{
	for (int i = 0; i < N_P; i++)
		delete population[i];
	delete []population;
	population = nullptr;
}

void save_result(int i)
{
//	std::cout << std::endl << i << std::endl;
	bool b_changed = false;
	for (int j = 0; j < CUT_LINE; j++) {
		double xx = population[j]->get_x();
		double yy = population[j]->get_y();
		double zz = population[j]->get_z();
//		std::cout << "x:" << xx << " y:" << yy << " z:" << zz << std::endl;

		if ((float)xx != (float)x[change_counter - 1][j])
			b_changed = true;
 		if ((float)yy != (float)y[change_counter - 1][j])
			b_changed = true;
		if ((float)zz != (float)z[change_counter - 1][j])
			b_changed = true;
	}

	if (b_changed == true) {
		for (int j = 0; j < CUT_LINE; j++) {
			x[change_counter][j] = population[j]->get_x();
			y[change_counter][j] = population[j]->get_y();
			z[change_counter][j] = population[j]->get_z();
		}
		change_counter++;
	}		
}

void optimize_genetic_algorithm()
{
	for (int i = 0; i < ITERATION; i++) {
		cross_over(i);
		ranking();
		save_result(i);
	}
}	

void write_data();

int main()
{
	srand ((unsigned)time(NULL));

	set_function_conditions();

	init_genetic_algorithm();

	optimize_genetic_algorithm();

	write_data();

	delete_genetic_algorithm();

	return 0;
}

void write_data()
{
	ofstream fout;
	fout.open("graph_data.txt");

	fout << str_formula << endl;
	cout << str_formula << endl ;

//	fout << ITERATION << endl;
//	cout << ITERATION << endl;
	fout << "ITERATION " << change_counter << endl;
	cout << "ITERATION " << change_counter << endl;

	fout << "SWARM_SIZE " << SWARM_SIZE << endl;
	cout << "SWARM_SIZE " << SWARM_SIZE << endl;

	fout << "X_MAX " << X_MAX << endl;
	cout << "X_MAX " << X_MAX << endl;

	fout << "X_MIN " << X_MIN << endl;
	cout << "X_MIN " << X_MIN << endl;

	fout << "Y_MAX " << Y_MAX << endl;
	cout << "Y_MAX " << Y_MAX << endl;

	fout << "Y_MIN " << Y_MIN << endl;
	cout << "Y_MIN " << Y_MIN << endl;

	fout << "Z_MAX " << Z_MAX << endl;
	cout << "Z_MAX " << Z_MAX << endl;

	fout << "Z_MIN " << Z_MIN << endl;
	cout << "Z_MIN " << Z_MIN << endl;

	fout << "DATA" << endl;
	cout << "DATA" << endl;

//	for (int i = 0; i < ITERATION; i++) {
	for (int i = 0; i < change_counter; i++) {
		cout << i << endl;
		for (int j = 0; j < SWARM_SIZE; j++) {
			cout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
			fout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
		}
		cout << endl;
		fout << endl;
	}

	fout.close();
}

five_well_potential_ga.cpp

#define _USE_MATH_DEFINES
#include <math.h>
#include <iostream>
#include <fstream>
using namespace std;

const int ITERATION = 1000; // 試行回数

const int SWARM_SIZE = 20; // 群の大きさ

double X_MAX = 10; // graphのX軸の表示範囲の上限
double X_MIN = -10; // graphのX軸の表示範囲の下限

double Y_MAX = 10; // giterationaphのY軸の表示範囲の上限
double Y_MIN = -10; // graphのY軸の表示範囲の下限

double Z_MAX = 200; //  graphのZ軸の表示範囲の上限
double Z_MIN = 0; // graphのZ軸の表示範囲の下限

string str_formula = "x * x + y * y"; // グラフに表示する数式

double x[ITERATION][SWARM_SIZE];
double y[ITERATION][SWARM_SIZE];
double z[ITERATION][SWARM_SIZE];

double function(double x, double y)
{
	double z;

	z = (1.0 - 1.0 / (1.0 + 0.05 * (x * x + (y - 10) * (y - 10)))
		- 1.0 / (1.0 + 0.05 * ((x - 10) * (x - 10) + y * y))
		- 1.5 / (1.0 + 0.03 * ((x + 10) * (x + 10) + y * y))
		- 2.0 / (1.0 + 0.05 * ((x - 5) * (x - 5) + (y + 10) * (y + 10)))
		- 1.0 / (1.0 + 0.1 * ((x + 5) * (x + 5) + (y + 10) * (y + 10))))
		* (1.0 + 0.0001 * exp(log(x * x + y * y) * 1.2));

	return z;
}

void set_function_conditions()
{
	X_MAX = 20;
	X_MIN = -20;
	Y_MAX = 20;
	Y_MIN = -20;
	Z_MAX = 1.5;
	Z_MIN = -1.5;
	str_formula = "(1 - 1 / (1 + 0.05 * (x**2 + (y - 10)**2)) - 1 / (1 + 0.05 * ((x - 10)**2 + y**2)) - 1.5 / (1 + 0.03 * ((x + 10)**2 + y**2)) - 2 / (1 + 0.05 * ((x - 5)**2 + (y + 10)**2)) - 1 / (1 + 0.1 * ((x + 5)**2 + (y + 10)**2))) * (1 + 0.0001 * exp(log(x**2 + y**2) * 1.2))";
}

int N_P = SWARM_SIZE * 2; // number of population

int CUT_LINE = SWARM_SIZE;
double A_GA = 0.3;
double B_GA = 0.1;
double MUTATION_ODDS = 0.1;
int change_counter = 0;

double get_random()
{
	return (double)rand() / ((double)RAND_MAX + 1);
}

double max1(double a, double b)
{
	double c = a;
	if (b > a)
		c = b;
	return c;
}

double min1(double a, double b)
{
	double c = a;
	if (b < a)
		c = b;
	return c;
}

class POPULATION
{
protected:
	double x, y, z;
	int* gene;
	long long weight_p;
	long long value_p;
public:
	POPULATION();
	~POPULATION();
	void calc_weight_and_value();
	long long get_value() { return value_p; }
	long long get_weight() { return weight_p; }
	void set_gene(int i, int g);
	int get_gene(int i) { return gene[i]; }
	void mutation();
	void initialize();
	void set_xyz(double x1, double y1, double z1) { x = x1; y = y1; z = z1; }
	double get_x() { return x; }
	double get_y() { return y; }
	double get_z() { return z; }
};

POPULATION** population = nullptr;

POPULATION::POPULATION()
{

}

POPULATION::~POPULATION()
{

}

void POPULATION::initialize()
{
	double x1 = get_random() * (X_MAX - X_MIN) + X_MIN;
	double y1 = get_random() * (Y_MAX - Y_MIN) + Y_MIN;
	double z1 = function(x1, y1);
	set_xyz(x1, y1, z1);
}

void POPULATION::mutation()
{
	initialize();
}

void create_object()
{
	population = new POPULATION*[N_P];
	for (int i = 0; i < N_P; i++)
		population[i] = new POPULATION;
}

void initialize_population()
{
	for (int i = 0; i < N_P; i++) 
		population[i]->initialize();
}

void ranking()
{
	for (int i = 0; i < N_P - 1; i++) {
		for (int j = i + 1; j < N_P; j++) {
			if (population[i]->get_z() > population[j]->get_z()) {
				POPULATION* tmp = population[i];
				population[i] = population[j];
				population[j] = tmp;
			}
		}
	}
}

void cross_over(int i)
{
	for (int j = CUT_LINE; j < N_P; j++) {
		int parent1 = rand() % CUT_LINE;
		int parent2 = rand() % CUT_LINE;
	
		double x1 = population[parent1]->get_x();
		double y1 = population[parent1]->get_y();
		double x2 = population[parent2]->get_x();
		double y2 = population[parent2]->get_y();
		double dx = abs(x1 - x2);
		double dy = abs(y1 - y2);
		double min_x = min(x1, x2);
		double min_y = min(y1, y2);
		double max_x = max(x1, x2);
		double max_y = max(y1, y2);
		double min_cx = min_x - A_GA * dx;
		double max_cx = max_x + A_GA * dx;
		double min_cy = min_y - A_GA * dy;
		double max_cy = max_y + A_GA * dy;
//		double next_x = get_random() * (max_cx - min_cx) + min_cx;
//		double next_y = get_random() * (max_cy - min_cy) + min_cy;
		double DX = (X_MAX - X_MIN) * B_GA * (get_random() - 0.5);
		double DY = (Y_MAX - Y_MIN) * B_GA * (get_random() - 0.5);
		double next_x = get_random() * (max_cx - min_cx) + min_cx + DX;
		double next_y = get_random() * (max_cy - min_cy) + min_cy + DY;
		next_x = max1(min1(next_x, X_MAX), X_MIN);
		next_y = max1(min1(next_y, Y_MAX), Y_MIN);
		double next_z = function(next_x, next_y);
		population[j]->set_xyz(next_x, next_y, next_z);

		if (get_random() < MUTATION_ODDS)
			population[j]->mutation();
	}
}

void init_genetic_algorithm()
{
	create_object();
	initialize_population();
	ranking();
	change_counter = 0;
	for (int j = 0; j < CUT_LINE; j++) {
		x[change_counter][j] = population[j]->get_x();
		y[change_counter][j] = population[j]->get_y();
		z[change_counter][j] = population[j]->get_z();
	}
	change_counter++;		
}

void delete_genetic_algorithm()
{
	for (int i = 0; i < N_P; i++)
		delete population[i];
	delete []population;
	population = nullptr;
}

void save_result(int i)
{
//	std::cout << std::endl << i << std::endl;
	bool b_changed = false;
	for (int j = 0; j < CUT_LINE; j++) {
		double xx = population[j]->get_x();
		double yy = population[j]->get_y();
		double zz = population[j]->get_z();
//		std::cout << "x:" << xx << " y:" << yy << " z:" << zz << std::endl;

		if ((float)xx != (float)x[change_counter - 1][j])
			b_changed = true;
 		if ((float)yy != (float)y[change_counter - 1][j])
			b_changed = true;
		if ((float)zz != (float)z[change_counter - 1][j])
			b_changed = true;
	}

	if (b_changed == true) {
		for (int j = 0; j < CUT_LINE; j++) {
			x[change_counter][j] = population[j]->get_x();
			y[change_counter][j] = population[j]->get_y();
			z[change_counter][j] = population[j]->get_z();
		}
		change_counter++;
	}		
}

void optimize_genetic_algorithm()
{
	for (int i = 0; i < ITERATION; i++) {
		cross_over(i);
		ranking();
		save_result(i);
	}
}	

void write_data();

int main()
{
	srand ((unsigned)time(NULL));

	set_function_conditions();

	init_genetic_algorithm();

	optimize_genetic_algorithm();

	write_data();

	delete_genetic_algorithm();

	return 0;
}

void write_data()
{
	ofstream fout;
	fout.open("graph_data.txt");

	fout << str_formula << endl;
	cout << str_formula << endl ;

//	fout << ITERATION << endl;
//	cout << ITERATION << endl;
	fout << "ITERATION " << change_counter << endl;
	cout << "ITERATION " << change_counter << endl;

	fout << "SWARM_SIZE " << SWARM_SIZE << endl;
	cout << "SWARM_SIZE " << SWARM_SIZE << endl;

	fout << "X_MAX " << X_MAX << endl;
	cout << "X_MAX " << X_MAX << endl;

	fout << "X_MIN " << X_MIN << endl;
	cout << "X_MIN " << X_MIN << endl;

	fout << "Y_MAX " << Y_MAX << endl;
	cout << "Y_MAX " << Y_MAX << endl;

	fout << "Y_MIN " << Y_MIN << endl;
	cout << "Y_MIN " << Y_MIN << endl;

	fout << "Z_MAX " << Z_MAX << endl;
	cout << "Z_MAX " << Z_MAX << endl;

	fout << "Z_MIN " << Z_MIN << endl;
	cout << "Z_MIN " << Z_MIN << endl;

	fout << "DATA" << endl;
	cout << "DATA" << endl;

//	for (int i = 0; i < ITERATION; i++) {
	for (int i = 0; i < change_counter; i++) {
		cout << i << endl;
		for (int j = 0; j < SWARM_SIZE; j++) {
			cout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
			fout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
		}
		cout << endl;
		fout << endl;
	}

	fout.close();
}

himmelblau_ga.cpp

#define _USE_MATH_DEFINES
#include <math.h>
#include <iostream>
#include <fstream>
using namespace std;

const int ITERATION = 1000; // 試行回数

const int SWARM_SIZE = 20; // 群の大きさ

double X_MAX = 10; // graphのX軸の表示範囲の上限
double X_MIN = -10; // graphのX軸の表示範囲の下限

double Y_MAX = 10; // giterationaphのY軸の表示範囲の上限
double Y_MIN = -10; // graphのY軸の表示範囲の下限

double Z_MAX = 200; //  graphのZ軸の表示範囲の上限
double Z_MIN = 0; // graphのZ軸の表示範囲の下限

string str_formula = "x * x + y * y"; // グラフに表示する数式

double x[ITERATION][SWARM_SIZE];
double y[ITERATION][SWARM_SIZE];
double z[ITERATION][SWARM_SIZE];

double function(double x, double y)
{
	double z;

	z = pow((x * x + y - 11), 2) + pow((x + y * y - 7), 2);

	return z;
}

void set_function_conditions()
{
	X_MAX = 5;
	X_MIN = -5;
	Y_MAX = 5;
	Y_MIN = -5;
	Z_MAX = 500;
	Z_MIN = 0;
	str_formula = "(x**2 + y - 11)**2 + (x + y**2 - 7)**2";
}

int N_P = SWARM_SIZE * 2; // number of population

int CUT_LINE = SWARM_SIZE;
double A_GA = 0.3;
double B_GA = 0.1;
double MUTATION_ODDS = 0.1;
int change_counter = 0;

double get_random()
{
	return (double)rand() / ((double)RAND_MAX + 1);
}

double max1(double a, double b)
{
	double c = a;
	if (b > a)
		c = b;
	return c;
}

double min1(double a, double b)
{
	double c = a;
	if (b < a)
		c = b;
	return c;
}

class POPULATION
{
protected:
	double x, y, z;
	int* gene;
	long long weight_p;
	long long value_p;
public:
	POPULATION();
	~POPULATION();
	void calc_weight_and_value();
	long long get_value() { return value_p; }
	long long get_weight() { return weight_p; }
	void set_gene(int i, int g);
	int get_gene(int i) { return gene[i]; }
	void mutation();
	void initialize();
	void set_xyz(double x1, double y1, double z1) { x = x1; y = y1; z = z1; }
	double get_x() { return x; }
	double get_y() { return y; }
	double get_z() { return z; }
};

POPULATION** population = nullptr;

POPULATION::POPULATION()
{

}

POPULATION::~POPULATION()
{

}

void POPULATION::initialize()
{
	double x1 = get_random() * (X_MAX - X_MIN) + X_MIN;
	double y1 = get_random() * (Y_MAX - Y_MIN) + Y_MIN;
	double z1 = function(x1, y1);
	set_xyz(x1, y1, z1);
}

void POPULATION::mutation()
{
	initialize();
}

void create_object()
{
	population = new POPULATION*[N_P];
	for (int i = 0; i < N_P; i++)
		population[i] = new POPULATION;
}

void initialize_population()
{
	for (int i = 0; i < N_P; i++) 
		population[i]->initialize();
}

void ranking()
{
	for (int i = 0; i < N_P - 1; i++) {
		for (int j = i + 1; j < N_P; j++) {
			if (population[i]->get_z() > population[j]->get_z()) {
				POPULATION* tmp = population[i];
				population[i] = population[j];
				population[j] = tmp;
			}
		}
	}
}

void cross_over(int i)
{
	for (int j = CUT_LINE; j < N_P; j++) {
		int parent1 = rand() % CUT_LINE;
		int parent2 = rand() % CUT_LINE;
	
		double x1 = population[parent1]->get_x();
		double y1 = population[parent1]->get_y();
		double x2 = population[parent2]->get_x();
		double y2 = population[parent2]->get_y();
		double dx = abs(x1 - x2);
		double dy = abs(y1 - y2);
		double min_x = min(x1, x2);
		double min_y = min(y1, y2);
		double max_x = max(x1, x2);
		double max_y = max(y1, y2);
		double min_cx = min_x - A_GA * dx;
		double max_cx = max_x + A_GA * dx;
		double min_cy = min_y - A_GA * dy;
		double max_cy = max_y + A_GA * dy;
//		double next_x = get_random() * (max_cx - min_cx) + min_cx;
//		double next_y = get_random() * (max_cy - min_cy) + min_cy;
		double DX = (X_MAX - X_MIN) * B_GA * (get_random() - 0.5);
		double DY = (Y_MAX - Y_MIN) * B_GA * (get_random() - 0.5);
		double next_x = get_random() * (max_cx - min_cx) + min_cx + DX;
		double next_y = get_random() * (max_cy - min_cy) + min_cy + DY;
		next_x = max1(min1(next_x, X_MAX), X_MIN);
		next_y = max1(min1(next_y, Y_MAX), Y_MIN);
		double next_z = function(next_x, next_y);
		population[j]->set_xyz(next_x, next_y, next_z);

		if (get_random() < MUTATION_ODDS)
			population[j]->mutation();
	}
}

void init_genetic_algorithm()
{
	create_object();
	initialize_population();
	ranking();
	change_counter = 0;
	for (int j = 0; j < CUT_LINE; j++) {
		x[change_counter][j] = population[j]->get_x();
		y[change_counter][j] = population[j]->get_y();
		z[change_counter][j] = population[j]->get_z();
	}
	change_counter++;		
}

void delete_genetic_algorithm()
{
	for (int i = 0; i < N_P; i++)
		delete population[i];
	delete []population;
	population = nullptr;
}

void save_result(int i)
{
//	std::cout << std::endl << i << std::endl;
	bool b_changed = false;
	for (int j = 0; j < CUT_LINE; j++) {
		double xx = population[j]->get_x();
		double yy = population[j]->get_y();
		double zz = population[j]->get_z();
//		std::cout << "x:" << xx << " y:" << yy << " z:" << zz << std::endl;

		if ((float)xx != (float)x[change_counter - 1][j])
			b_changed = true;
 		if ((float)yy != (float)y[change_counter - 1][j])
			b_changed = true;
		if ((float)zz != (float)z[change_counter - 1][j])
			b_changed = true;
	}

	if (b_changed == true) {
		for (int j = 0; j < CUT_LINE; j++) {
			x[change_counter][j] = population[j]->get_x();
			y[change_counter][j] = population[j]->get_y();
			z[change_counter][j] = population[j]->get_z();
		}
		change_counter++;
	}		
}

void optimize_genetic_algorithm()
{
	for (int i = 0; i < ITERATION; i++) {
		cross_over(i);
		ranking();
		save_result(i);
	}
}	

void write_data();

int main()
{
	srand ((unsigned)time(NULL));

	set_function_conditions();

	init_genetic_algorithm();

	optimize_genetic_algorithm();

	write_data();

	delete_genetic_algorithm();

	return 0;
}

void write_data()
{
	ofstream fout;
	fout.open("graph_data.txt");

	fout << str_formula << endl;
	cout << str_formula << endl ;

//	fout << ITERATION << endl;
//	cout << ITERATION << endl;
	fout << "ITERATION " << change_counter << endl;
	cout << "ITERATION " << change_counter << endl;

	fout << "SWARM_SIZE " << SWARM_SIZE << endl;
	cout << "SWARM_SIZE " << SWARM_SIZE << endl;

	fout << "X_MAX " << X_MAX << endl;
	cout << "X_MAX " << X_MAX << endl;

	fout << "X_MIN " << X_MIN << endl;
	cout << "X_MIN " << X_MIN << endl;

	fout << "Y_MAX " << Y_MAX << endl;
	cout << "Y_MAX " << Y_MAX << endl;

	fout << "Y_MIN " << Y_MIN << endl;
	cout << "Y_MIN " << Y_MIN << endl;

	fout << "Z_MAX " << Z_MAX << endl;
	cout << "Z_MAX " << Z_MAX << endl;

	fout << "Z_MIN " << Z_MIN << endl;
	cout << "Z_MIN " << Z_MIN << endl;

	fout << "DATA" << endl;
	cout << "DATA" << endl;

//	for (int i = 0; i < ITERATION; i++) {
	for (int i = 0; i < change_counter; i++) {
		cout << i << endl;
		for (int j = 0; j < SWARM_SIZE; j++) {
			cout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
			fout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
		}
		cout << endl;
		fout << endl;
	}

	fout.close();
}

rosenbrock_ga.cpp

#define _USE_MATH_DEFINES
#include <math.h>
#include <iostream>
#include <fstream>
using namespace std;

const int ITERATION = 1000; // 試行回数

const int SWARM_SIZE = 20; // 群の大きさ

double X_MAX = 10; // graphのX軸の表示範囲の上限
double X_MIN = -10; // graphのX軸の表示範囲の下限

double Y_MAX = 10; // giterationaphのY軸の表示範囲の上限
double Y_MIN = -10; // graphのY軸の表示範囲の下限

double Z_MAX = 200; //  graphのZ軸の表示範囲の上限
double Z_MIN = 0; // graphのZ軸の表示範囲の下限

string str_formula = "x * x + y * y"; // グラフに表示する数式

double x[ITERATION][SWARM_SIZE];
double y[ITERATION][SWARM_SIZE];
double z[ITERATION][SWARM_SIZE];

double function(double x, double y)
{
	double z;

	z = 100 * (y - x * x) * (y  - x * x) + (1 - x) * (1 - x);

	return z;
}

void set_function_conditions()
{
	X_MAX = 2.048;
	X_MIN = -2.048;
	Y_MAX = 2.048;
	Y_MIN = -2.048;
	Z_MAX = 4000;
	Z_MIN = 0;
	str_formula = "100 * (y - x * x)**2 + (1 - x)**2";
}

int N_P = SWARM_SIZE * 2; // number of population

int CUT_LINE = SWARM_SIZE;
double A_GA = 0.3;
double B_GA = 0.1;
double MUTATION_ODDS = 0.1;
int change_counter = 0;

double get_random()
{
	return (double)rand() / ((double)RAND_MAX + 1);
}

double max1(double a, double b)
{
	double c = a;
	if (b > a)
		c = b;
	return c;
}

double min1(double a, double b)
{
	double c = a;
	if (b < a)
		c = b;
	return c;
}

class POPULATION
{
protected:
	double x, y, z;
	int* gene;
	long long weight_p;
	long long value_p;
public:
	POPULATION();
	~POPULATION();
	void calc_weight_and_value();
	long long get_value() { return value_p; }
	long long get_weight() { return weight_p; }
	void set_gene(int i, int g);
	int get_gene(int i) { return gene[i]; }
	void mutation();
	void initialize();
	void set_xyz(double x1, double y1, double z1) { x = x1; y = y1; z = z1; }
	double get_x() { return x; }
	double get_y() { return y; }
	double get_z() { return z; }
};

POPULATION** population = nullptr;

POPULATION::POPULATION()
{

}

POPULATION::~POPULATION()
{

}

void POPULATION::initialize()
{
	double x1 = get_random() * (X_MAX - X_MIN) + X_MIN;
	double y1 = get_random() * (Y_MAX - Y_MIN) + Y_MIN;
	double z1 = function(x1, y1);
	set_xyz(x1, y1, z1);
}

void POPULATION::mutation()
{
	initialize();
}

void create_object()
{
	population = new POPULATION*[N_P];
	for (int i = 0; i < N_P; i++)
		population[i] = new POPULATION;
}

void initialize_population()
{
	for (int i = 0; i < N_P; i++) 
		population[i]->initialize();
}

void ranking()
{
	for (int i = 0; i < N_P - 1; i++) {
		for (int j = i + 1; j < N_P; j++) {
			if (population[i]->get_z() > population[j]->get_z()) {
				POPULATION* tmp = population[i];
				population[i] = population[j];
				population[j] = tmp;
			}
		}
	}
}

void cross_over(int i)
{
	for (int j = CUT_LINE; j < N_P; j++) {
		int parent1 = rand() % CUT_LINE;
		int parent2 = rand() % CUT_LINE;
	
		double x1 = population[parent1]->get_x();
		double y1 = population[parent1]->get_y();
		double x2 = population[parent2]->get_x();
		double y2 = population[parent2]->get_y();
		double dx = abs(x1 - x2);
		double dy = abs(y1 - y2);
		double min_x = min(x1, x2);
		double min_y = min(y1, y2);
		double max_x = max(x1, x2);
		double max_y = max(y1, y2);
		double min_cx = min_x - A_GA * dx;
		double max_cx = max_x + A_GA * dx;
		double min_cy = min_y - A_GA * dy;
		double max_cy = max_y + A_GA * dy;
//		double next_x = get_random() * (max_cx - min_cx) + min_cx;
//		double next_y = get_random() * (max_cy - min_cy) + min_cy;
		double DX = (X_MAX - X_MIN) * B_GA * (get_random() - 0.5);
		double DY = (Y_MAX - Y_MIN) * B_GA * (get_random() - 0.5);
		double next_x = get_random() * (max_cx - min_cx) + min_cx + DX;
		double next_y = get_random() * (max_cy - min_cy) + min_cy + DY;
		next_x = max1(min1(next_x, X_MAX), X_MIN);
		next_y = max1(min1(next_y, Y_MAX), Y_MIN);
		double next_z = function(next_x, next_y);
		population[j]->set_xyz(next_x, next_y, next_z);

		if (get_random() < MUTATION_ODDS)
			population[j]->mutation();
	}
}

void init_genetic_algorithm()
{
	create_object();
	initialize_population();
	ranking();
	change_counter = 0;
	for (int j = 0; j < CUT_LINE; j++) {
		x[change_counter][j] = population[j]->get_x();
		y[change_counter][j] = population[j]->get_y();
		z[change_counter][j] = population[j]->get_z();
	}
	change_counter++;		
}

void delete_genetic_algorithm()
{
	for (int i = 0; i < N_P; i++)
		delete population[i];
	delete []population;
	population = nullptr;
}

void save_result(int i)
{
//	std::cout << std::endl << i << std::endl;
	bool b_changed = false;
	for (int j = 0; j < CUT_LINE; j++) {
		double xx = population[j]->get_x();
		double yy = population[j]->get_y();
		double zz = population[j]->get_z();
//		std::cout << "x:" << xx << " y:" << yy << " z:" << zz << std::endl;

		if ((float)xx != (float)x[change_counter - 1][j])
			b_changed = true;
 		if ((float)yy != (float)y[change_counter - 1][j])
			b_changed = true;
		if ((float)zz != (float)z[change_counter - 1][j])
			b_changed = true;
	}

	if (b_changed == true) {
		for (int j = 0; j < CUT_LINE; j++) {
			x[change_counter][j] = population[j]->get_x();
			y[change_counter][j] = population[j]->get_y();
			z[change_counter][j] = population[j]->get_z();
		}
		change_counter++;
	}		
}

void optimize_genetic_algorithm()
{
	for (int i = 0; i < ITERATION; i++) {
		cross_over(i);
		ranking();
		save_result(i);
	}
}	

void write_data();

int main()
{
	srand ((unsigned)time(NULL));

	set_function_conditions();

	init_genetic_algorithm();

	optimize_genetic_algorithm();

	write_data();

	delete_genetic_algorithm();

	return 0;
}

void write_data()
{
	ofstream fout;
	fout.open("graph_data.txt");

	fout << str_formula << endl;
	cout << str_formula << endl ;

//	fout << ITERATION << endl;
//	cout << ITERATION << endl;
	fout << "ITERATION " << change_counter << endl;
	cout << "ITERATION " << change_counter << endl;

	fout << "SWARM_SIZE " << SWARM_SIZE << endl;
	cout << "SWARM_SIZE " << SWARM_SIZE << endl;

	fout << "X_MAX " << X_MAX << endl;
	cout << "X_MAX " << X_MAX << endl;

	fout << "X_MIN " << X_MIN << endl;
	cout << "X_MIN " << X_MIN << endl;

	fout << "Y_MAX " << Y_MAX << endl;
	cout << "Y_MAX " << Y_MAX << endl;

	fout << "Y_MIN " << Y_MIN << endl;
	cout << "Y_MIN " << Y_MIN << endl;

	fout << "Z_MAX " << Z_MAX << endl;
	cout << "Z_MAX " << Z_MAX << endl;

	fout << "Z_MIN " << Z_MIN << endl;
	cout << "Z_MIN " << Z_MIN << endl;

	fout << "DATA" << endl;
	cout << "DATA" << endl;

//	for (int i = 0; i < ITERATION; i++) {
	for (int i = 0; i < change_counter; i++) {
		cout << i << endl;
		for (int j = 0; j < SWARM_SIZE; j++) {
			cout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
			fout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
		}
		cout << endl;
		fout << endl;
	}

	fout.close();
}

rastrigin_ga.cpp

#define _USE_MATH_DEFINES
#include <math.h>
#include <iostream>
#include <fstream>
using namespace std;

const int ITERATION = 1000; // 試行回数

const int SWARM_SIZE = 20; // 群の大きさ

double X_MAX = 10; // graphのX軸の表示範囲の上限
double X_MIN = -10; // graphのX軸の表示範囲の下限

double Y_MAX = 10; // giterationaphのY軸の表示範囲の上限
double Y_MIN = -10; // graphのY軸の表示範囲の下限

double Z_MAX = 200; //  graphのZ軸の表示範囲の上限
double Z_MIN = 0; // graphのZ軸の表示範囲の下限

string str_formula = "x * x + y * y"; // グラフに表示する数式

double x[ITERATION][SWARM_SIZE];
double y[ITERATION][SWARM_SIZE];
double z[ITERATION][SWARM_SIZE];

double function(double x, double y)
{
	double z;

	z = 20.0 + x * x - 10.0 * cos(2.0 * M_PI * x) + y * y - 10.0 * cos(2.0 * M_PI * y);

	return z;
}

void set_function_conditions()
{
	X_MAX = 10;
	X_MIN = -10;
	Y_MAX = 10;
	Y_MIN = -10;
	Z_MAX = 200;
	Z_MIN = 0;
	str_formula = "20.0 + x * x - 10.0 * cos(2.0 * pi * x) + y * y - 10.0 * cos(2.0 * pi * y)";
}

int N_P = SWARM_SIZE * 2; // number of population

int CUT_LINE = SWARM_SIZE;
double A_GA = 0.3;
double B_GA = 0.1;
double MUTATION_ODDS = 0.1;
int change_counter = 0;

double get_random()
{
	return (double)rand() / ((double)RAND_MAX + 1);
}

double max1(double a, double b)
{
	double c = a;
	if (b > a)
		c = b;
	return c;
}

double min1(double a, double b)
{
	double c = a;
	if (b < a)
		c = b;
	return c;
}

class POPULATION
{
protected:
	double x, y, z;
	int* gene;
	long long weight_p;
	long long value_p;
public:
	POPULATION();
	~POPULATION();
	void calc_weight_and_value();
	long long get_value() { return value_p; }
	long long get_weight() { return weight_p; }
	void set_gene(int i, int g);
	int get_gene(int i) { return gene[i]; }
	void mutation();
	void initialize();
	void set_xyz(double x1, double y1, double z1) { x = x1; y = y1; z = z1; }
	double get_x() { return x; }
	double get_y() { return y; }
	double get_z() { return z; }
};

POPULATION** population = nullptr;

POPULATION::POPULATION()
{

}

POPULATION::~POPULATION()
{

}

void POPULATION::initialize()
{
	double x1 = get_random() * (X_MAX - X_MIN) + X_MIN;
	double y1 = get_random() * (Y_MAX - Y_MIN) + Y_MIN;
	double z1 = function(x1, y1);
	set_xyz(x1, y1, z1);
}

void POPULATION::mutation()
{
	initialize();
}

void create_object()
{
	population = new POPULATION*[N_P];
	for (int i = 0; i < N_P; i++)
		population[i] = new POPULATION;
}

void initialize_population()
{
	for (int i = 0; i < N_P; i++) 
		population[i]->initialize();
}

void ranking()
{
	for (int i = 0; i < N_P - 1; i++) {
		for (int j = i + 1; j < N_P; j++) {
			if (population[i]->get_z() > population[j]->get_z()) {
				POPULATION* tmp = population[i];
				population[i] = population[j];
				population[j] = tmp;
			}
		}
	}
}

void cross_over(int i)
{
	for (int j = CUT_LINE; j < N_P; j++) {
		int parent1 = rand() % CUT_LINE;
		int parent2 = rand() % CUT_LINE;
	
		double x1 = population[parent1]->get_x();
		double y1 = population[parent1]->get_y();
		double x2 = population[parent2]->get_x();
		double y2 = population[parent2]->get_y();
		double dx = abs(x1 - x2);
		double dy = abs(y1 - y2);
		double min_x = min(x1, x2);
		double min_y = min(y1, y2);
		double max_x = max(x1, x2);
		double max_y = max(y1, y2);
		double min_cx = min_x - A_GA * dx;
		double max_cx = max_x + A_GA * dx;
		double min_cy = min_y - A_GA * dy;
		double max_cy = max_y + A_GA * dy;
//		double next_x = get_random() * (max_cx - min_cx) + min_cx;
//		double next_y = get_random() * (max_cy - min_cy) + min_cy;
		double DX = (X_MAX - X_MIN) * B_GA * (get_random() - 0.5);
		double DY = (Y_MAX - Y_MIN) * B_GA * (get_random() - 0.5);
		double next_x = get_random() * (max_cx - min_cx) + min_cx + DX;
		double next_y = get_random() * (max_cy - min_cy) + min_cy + DY;
		next_x = max1(min1(next_x, X_MAX), X_MIN);
		next_y = max1(min1(next_y, Y_MAX), Y_MIN);
		double next_z = function(next_x, next_y);
		population[j]->set_xyz(next_x, next_y, next_z);

		if (get_random() < MUTATION_ODDS)
			population[j]->mutation();
	}
}

void init_genetic_algorithm()
{
	create_object();
	initialize_population();
	ranking();
	change_counter = 0;
	for (int j = 0; j < CUT_LINE; j++) {
		x[change_counter][j] = population[j]->get_x();
		y[change_counter][j] = population[j]->get_y();
		z[change_counter][j] = population[j]->get_z();
	}
	change_counter++;		
}

void delete_genetic_algorithm()
{
	for (int i = 0; i < N_P; i++)
		delete population[i];
	delete []population;
	population = nullptr;
}

void save_result(int i)
{
//	std::cout << std::endl << i << std::endl;
	bool b_changed = false;
	for (int j = 0; j < CUT_LINE; j++) {
		double xx = population[j]->get_x();
		double yy = population[j]->get_y();
		double zz = population[j]->get_z();
//		std::cout << "x:" << xx << " y:" << yy << " z:" << zz << std::endl;

		if ((float)xx != (float)x[change_counter - 1][j])
			b_changed = true;
 		if ((float)yy != (float)y[change_counter - 1][j])
			b_changed = true;
		if ((float)zz != (float)z[change_counter - 1][j])
			b_changed = true;
	}

	if (b_changed == true) {
		for (int j = 0; j < CUT_LINE; j++) {
			x[change_counter][j] = population[j]->get_x();
			y[change_counter][j] = population[j]->get_y();
			z[change_counter][j] = population[j]->get_z();
		}
		change_counter++;
	}		
}

void optimize_genetic_algorithm()
{
	for (int i = 0; i < ITERATION; i++) {
		cross_over(i);
		ranking();
		save_result(i);
	}
}	

void write_data();

int main()
{
	srand ((unsigned)time(NULL));

	set_function_conditions();

	init_genetic_algorithm();

	optimize_genetic_algorithm();

	write_data();

	delete_genetic_algorithm();

	return 0;
}

void write_data()
{
	ofstream fout;
	fout.open("graph_data.txt");

	fout << str_formula << endl;
	cout << str_formula << endl ;

//	fout << ITERATION << endl;
//	cout << ITERATION << endl;
	fout << "ITERATION " << change_counter << endl;
	cout << "ITERATION " << change_counter << endl;

	fout << "SWARM_SIZE " << SWARM_SIZE << endl;
	cout << "SWARM_SIZE " << SWARM_SIZE << endl;

	fout << "X_MAX " << X_MAX << endl;
	cout << "X_MAX " << X_MAX << endl;

	fout << "X_MIN " << X_MIN << endl;
	cout << "X_MIN " << X_MIN << endl;

	fout << "Y_MAX " << Y_MAX << endl;
	cout << "Y_MAX " << Y_MAX << endl;

	fout << "Y_MIN " << Y_MIN << endl;
	cout << "Y_MIN " << Y_MIN << endl;

	fout << "Z_MAX " << Z_MAX << endl;
	cout << "Z_MAX " << Z_MAX << endl;

	fout << "Z_MIN " << Z_MIN << endl;
	cout << "Z_MIN " << Z_MIN << endl;

	fout << "DATA" << endl;
	cout << "DATA" << endl;

//	for (int i = 0; i < ITERATION; i++) {
	for (int i = 0; i < change_counter; i++) {
		cout << i << endl;
		for (int j = 0; j < SWARM_SIZE; j++) {
			cout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
			fout << x[i][j] << " " << y[i][j] << " " << z[i][j] << endl;
		}
		cout << endl;
		fout << endl;
	}

	fout.close();
}

出力結果 (sphere_ga)

x * x + y * y
ITERATION 78
SWARM_SIZE 20
X_MAX 10
X_MIN -10
Y_MAX 10
Y_MIN -10
Z_MAX 200
Z_MIN 0
DATA
-0.509644 0.926514 1.11816
1.26709 1.24573 3.15735
0.991211 2.72217 8.3927
2.15637 1.9574 8.48135
-3.52966 -1.24573 14.0104
0.19104 4.00696 16.0922
2.38953 -3.35327 16.9543
-0.177002 -4.17664 17.4756
0.278931 -4.79248 23.0457
2.2876 4.68506 27.1829
5.08972 2.05383 30.1235
-3.30872 4.6814 32.8631
1.79565 5.70496 35.7709
-5.24597 3.90137 42.7409
4.1095 5.29297 44.9035
-6.9812 2.47498 54.8627
7.35046 -2.62085 60.8982
8.0249 -0.209351 64.4429
8.15735 1.01868 67.58
-7.04407 4.61243 70.8934

-0.509644 0.926514 1.11816
1.54323 -0.386026 2.53057
1.26709 1.24573 3.15735
2.04824 -0.79537 4.8279
0.991211 2.72217 8.3927
2.15637 1.9574 8.48135
1.02831 2.93215 9.65489
1.38677 3.23235 12.3712
0.734293 3.51365 12.8849
-3.52966 -1.24573 14.0104
0.19104 4.00696 16.0922
-0.194648 4.05091 16.4478
2.38953 -3.35327 16.9543
-0.177002 -4.17664 17.4756
-0.956997 -4.08895 17.6354
2.06133 -3.66054 17.6487
2.45162 3.76245 20.1665
0.278931 -4.79248 23.0457
4.86389 0.236206 23.7132
2.2876 4.68506 27.1829

-0.509644 0.926514 1.11816
0.810067 -1.06241 1.78493
1.22274 0.762198 2.07604
0.73347 -1.24746 2.09414
1.54323 -0.386026 2.53057
1.19641 1.05328 2.54081
1.62358 0.518242 2.90459
1.26709 1.24573 3.15735
-0.915695 1.71772 3.78907
1.54106 1.35039 4.19841
-1.15405 1.82138 4.64924
2.04824 -0.79537 4.8279
2.05438 0.946472 5.1163
1.52393 -1.95427 6.14151
-0.793891 2.7244 8.0526
0.991211 2.72217 8.3927
2.15637 1.9574 8.48135
1.02831 2.93215 9.65489
1.38677 3.23235 12.3712
0.734293 3.51365 12.8849