PR-20201212

LinearAlgebra/GA.cpp

@@ -0,0 +1,111 @@
+#include"linearvector.h"
+#include"matrix.h"
+#include"QuetionA.h"
+#include"GA.h"
+
+/*
+source: The Matrix object you want to solve;
+size: Number of individuals in the species (200 in the article);
+max_generation: Maximum iteration times (1000 in the article);
+cross_pos: Possibility of crossing (0.6 in the article);
+cross_rate: Degree of crossing (0.66 in the article);
+mutate_pos: Possibility of mutation (0.05 in the article);
+mutate_rate: Dependence of mutation degree on the generation times (0.8 in the article);
+*/
+Eigenvalue::Eigenvalue(const Matrix& source, const int& size, const int& max_generation, const double& cross_pos, const double& cross_rate, const double& mutate_pos, const double& mutate_rate) {
+	this->source = source;
+	this->size = size;
+	this->max_generation = max_generation;
+	this->cross_pos = cross_pos;
+	this->cross_rate = cross_rate;
+	this->mutate_pos = mutate_pos;
+	this->mutate_rate = mutate_rate;
+}
+
+void Eigenvalue::initiallize(std::mt19937 gen, std::uniform_real_distribution<double> u) {
+	generation = 0;
+	parents.clear(); children.clear();
+	eigenpolyval.clear(); sum_cost.clear();
+	upper_bound = 0.0;
+	for (int i = 0; i < source.get_row(); i++) {
+		for (int j = 0; j < source.get_column(); j++) {
+			upper_bound += std::fabs(source(i,j));
+		}
+	}
+	lower_bound = -upper_bound;
+	answers.clear();
+	for (int i = 0; i < size; i++) {
+		parents.push_back(u(gen)*(upper_bound-lower_bound)+lower_bound);
+	}
+}
+double Eigenvalue::solve_eigenpolyval(double lam) {
+	return std::fabs((source-lam*identity_matrix(source.get_column())).determinant());
+}
+void Eigenvalue::solve_total() {
+	double max_eigenpolyval = 0.0;
+	eigenpolyval.clear(); sum_cost.clear();
+	for (int i = 0; i < size; i++) {
+		eigenpolyval.push_back(solve_eigenpolyval(parents[i]));
+		if (eigenpolyval[i] > max_eigenpolyval) max_eigenpolyval = eigenpolyval[i];
+	}
+	sum_cost.push_back(max_eigenpolyval - eigenpolyval[0]);
+	for (int i = 1; i < size; i++) {
+		sum_cost.push_back(sum_cost[i-1] + max_eigenpolyval - eigenpolyval[i]);
+	}
+}
+int Eigenvalue::choose(double proportion) {
+	for (int i = 0; i < size; i++) {
+		if (sum_cost[i] > proportion * sum_cost[size - 1]) return i;
+	}
+}
+void Eigenvalue::mutate(int idx, bool opt, double ran) {
+	double tmp = children[idx];
+	if (opt) {
+		children[idx] += (upper_bound - tmp) * (1.0 - std::pow(ran, mutate_rate * (double)generation / max_generation));
+	}
+	else {
+		children[idx] -= (tmp - lower_bound) * (1.0 - std::pow(ran, mutate_rate * (double)generation / max_generation));
+	}
+}
+void Eigenvalue::crossover(int ida, int idb, bool opt) {
+	if (opt) {
+		children.push_back(cross_rate * parents[idb] + (1.0 - cross_rate) * parents[ida]);
+		children.push_back(cross_rate * parents[ida] + (1.0 - cross_rate) * parents[idb]);
+	}
+	else {
+		children.push_back(parents[ida]);
+		children.push_back(parents[idb]);
+	}
+}
+void Eigenvalue::produce_next(std::mt19937 gen, std::uniform_real_distribution<double> u) {
+	children.clear();
+	generation++;
+	solve_total();
+	int ida, idb;
+	for (int i = 0; i < (size >> 1); i++) {
+		ida = choose(u(gen));
+		idb = choose(u(gen));
+		crossover(ida, idb, u(gen)<cross_pos);
+	}
+	for (int i = 0; i < size; i++) if (u(gen) < mutate_pos) {
+		mutate(i,u(gen)<0.5,u(gen));
+	}
+	std::swap(parents, children);
+}
+void Eigenvalue::solve_eigenvalue() {
+	std::mt19937 gen(time(NULL));
+	std::uniform_real_distribution<double> u(0.0, 1.0);
+	if (source.get_column() != source.get_row()) {
+		std::cerr << "[ERROR]Non-square matrices have no eigenvalue!!" << std::endl;
+		return;
+	}
+	initiallize(gen,u);
+	for (int i = 0; i < max_generation; i++) {
+		produce_next(gen, u);
+	}
+	std::sort(parents.begin(), parents.end());
+
+	for (int i = 0; i < size; i++) {
+		if (solve_eigenpolyval(parents[i]) < eps2&&(answers.empty()||parents[i]-*(answers.end()-1)>eps))	answers.push_back(parents[i]);
+	}
+}

LinearAlgebra/GA.h

@@ -0,0 +1,45 @@
+#pragma once
+#ifndef EVOLUTION_H
+#define EVOLUTION_H
+#include"linearvector.h"
+#include"matrix.h"
+#include"QuetionA.h"
+#include<iostream>
+#include<cmath>
+#include<algorithm>
+#include<vector>
+#include<random>
+#include<ctime>
+
+class Eigenvalue {
+private:
+	Matrix source;
+	Matrix tmp;
+	int size,generation,max_generation;
+	std::vector<double> parents;
+	std::vector<double> children;
+	std::vector<double> eigenpolyval;
+	std::vector<double> sum_cost;
+	double cross_pos, cross_rate, mutate_pos, mutate_rate;
+	double eps = 1e-3, eps2=1e-2;
+public:
+
+	double upper_bound, lower_bound;
+	std::vector<double> answers;
+
+	Eigenvalue(const Matrix&, const int&, const int&, const double&, const double&, const double&, const double&);
+	~Eigenvalue()=default;
+
+	void initiallize(std::mt19937, std::uniform_real_distribution<double>);
+	double solve_eigenpolyval(double x);
+	void solve_total();
+	int choose(double);
+	void mutate(int, bool, double);
+	void crossover(int, int, bool);
+	void produce_next(std::mt19937, std::uniform_real_distribution<double>);
+	void solve_eigenvalue();
+};
+
+
+#endif
+

LinearAlgebra/LinearAlgebra.vcxproj

@@ -140,13 +140,15 @@
     </Link>
   </ItemDefinitionGroup>
   <ItemGroup>
+    <ClCompile Include="GA.cpp" />
     <ClCompile Include="linearvector.cpp" />
     <ClCompile Include="matrix.cpp" />
     <ClCompile Include="main.cpp" />
     <ClCompile Include="QuetionA.cpp" />
     <ClCompile Include="QuetionB.cpp" />
   </ItemGroup>
   <ItemGroup>
+    <ClInclude Include="GA.h" />
     <ClInclude Include="linearvector.h" />
     <ClInclude Include="matrix.h" />
     <ClInclude Include="QuetionA.h" />

LinearAlgebra/LinearAlgebra.vcxproj.filters

@@ -26,6 +26,9 @@
     <ClCompile Include="QuetionB.cpp">
       <Filter>源文件</Filter>
     </ClCompile>
+    <ClCompile Include="GA.cpp">
+      <Filter>源文件</Filter>
+    </ClCompile>
   </ItemGroup>
   <ItemGroup>
     <ClInclude Include="matrix.h">
@@ -40,5 +43,8 @@
     <ClInclude Include="QuetionB.h">
       <Filter>头文件</Filter>
     </ClInclude>
+    <ClInclude Include="GA.h">
+      <Filter>头文件</Filter>
+    </ClInclude>
   </ItemGroup>
 </Project>

LinearAlgebra/main.cpp

@@ -1,12 +1,21 @@
 #include<iostream>
 #include"QuetionA.h"
 #include"QuetionB.h"
+#include"GA.h"
 int main()
 {
-	std::cout<<"////////////////////////////////////////QUETION1////////////////////////////////////////"<<std::endl;
-	QuestionA::main();
-	std::cout<<std::endl<<std::endl;
-	std::cout<<"////////////////////////////////////////QUETION2////////////////////////////////////////"<<std::endl;
-	QuestionB::main();
+//	std::cout<<"////////////////////////////////////////QUETION1////////////////////////////////////////"<<std::endl;
+//	QuestionA::main();
+//	std::cout<<std::endl<<std::endl;
+//	std::cout<<"////////////////////////////////////////QUETION2////////////////////////////////////////"<<std::endl;
+//	QuestionB::main();
+	double a[] = { 1,0,0,0,2,0,0,0,3 };
+	Matrix A = Matrix((unsigned int)3, (unsigned int)3, a);
+	Eigenvalue eigenv = Eigenvalue(A, 500, 5000, 0.6, 0.66, 0.05, 0.8);
+	eigenv.solve_eigenvalue();
+	printf("%lld\n", eigenv.answers.size());
+	for (int i = 0; i < eigenv.answers.size(); i++) {
+		printf("%.4lf\n",eigenv.answers[i]);
+	}
 	return 0;
 }