#include"evops.h"
namespace EVO {
	void shuffle(int* a, int n)
	{
		for (int i = n; i > 0; i--)
		{
			int k = rand() % i;
			int temp = a[i - 1];
			a[i - 1] = a[k];
			a[k] = temp;
		}
	}

	bool  comp(const Account induval1, const Account induval2)
	{
		return induval1.energy < induval2.energy;
	}

	bool convergent(ConInfo* ifo1, ConInfo* ifo2)
	{ // return fabs( (ifo1->min_energy) - (ifo2->min_energy))< precision ;
		return 0;
	}

	Population::Population()
	{
		this->Generation = this->POPULATION_NUM = this->maxiterator = 0;
		this->init_volume = 0; this->max_formula_atom = 12;
		this->maxcaltime = 2;
		FILE* getpop = fopen("evops.in", "r");
		if (getpop == NULL)
			perror("evops.in is not exist!");
		char popnum[80];
		int flag1 = 0, flag2 = 0;
		while (!feof(getpop))
		{
			char temp[80];
			fgets(temp, 80, getpop);
			string value = GetValueAfterEqual(temp);
			if (strstr(temp, "POPULATION") && !flag1 && !value.empty())
			{
				this->POPULATION_NUM = atoi(value.c_str());
				flag1 = 1;
			}
			if (strstr(temp, "MAX_ITERATOR") && !flag2 && !value.empty())
			{
				this->maxiterator = atoi(value.c_str());
				flag2 = 1;
			}
			if (strstr(temp, "VOLUME") && !value.empty())
			{
				this->init_volume = atof(value.c_str());
			}
			if (strstr(temp, "FormulaAtomNum") && !value.empty())
			{
				this->max_formula_atom = atoi(value.c_str());
			}
			if (strstr(temp, "MAXCALTIME") && !value.empty())
			{
				this->maxcaltime = atof(value.c_str());
			}
			if (strstr(temp, "SingleEnergy") && !value.empty())
			{
				int iat = 0, cnt = 0, flag = 0;
				vector<string> string_single_energy(1);
				for (int i = 0; i < strlen(temp); i++)
				{
					if (temp[i] == '=')
						iat = 1;
					if (iat)
					{
						if (temp[i] != '=' && temp[i] != ' ' && temp[i] != '\t' && temp[i] != '\n' && temp[i] != '\r')
						{
							string_single_energy[cnt] += temp[i];
							flag = 1;
						}
						else if (temp[i] == ' ' && flag)
						{
							cnt++;
							string_single_energy.resize(string_single_energy.size() + 1);
							flag = 0;
						}
					}
				}
				for (int i = 0; i < string_single_energy.size(); i++)
				{
					double tmp = atof(string_single_energy[i].c_str());
					this->single_energy.push_back(tmp);
				}
			}
			if (strstr(temp, "ATOM_NUM") && !value.empty())
			{
				vector<string> Atom(1);
				int iat = 0, cnt = 0, flag = 0;
				for (int i = 0; i < strlen(temp); i++)
				{
					if (temp[i] == '=')
						iat = 1;
					if (iat)
					{
						if (temp[i] != '=' && temp[i] != ' ' && temp[i] != '\t' && temp[i] != '\n' && temp[i] != '\r')
						{
							Atom[cnt] += temp[i];
							flag = 1;
						}
						else if (temp[i] == ' ' && flag)
						{
							cnt++;
							Atom.resize(Atom.size() + 1);
							flag = 0;
						}
					}
				}
				for (int i = 0; i < Atom.size(); i++)
				{
					int tmp = atoi(Atom[i].c_str());
					if (tmp != 0)
						this->init_atom.push_back(tmp);
				}
			}
			if (strstr(temp, "ELEMENT_TYPE"))
			{
				int iat = 0, cnt = 0, flag = 0;
				this->init_element.resize(1);
				for (int i = 0; i < strlen(temp); i++)
				{
					if (temp[i] == '=')
						iat = 1;
					if (iat)
					{
						if (temp[i] != '=' && temp[i] != ' ' && temp[i] != '\t' && temp[i] != '\n' && temp[i] != '\r')
						{
							this->init_element[cnt] += temp[i];
							flag = 1;
						}
						else if (temp[i] == ' ' && flag)
						{
							cnt++;
							this->init_element.resize(this->init_element.size() + 1);
							flag = 0;
						}
					}
				}
				if (*(this->init_element.end() - 1) == "")
					this->init_element.pop_back();
			}
		}
		fclose(getpop);
		population.resize(POPULATION_NUM);
		population_child.resize(POPULATION_NUM);
		vector<Account>::iterator ite = population.begin();
		vector<Account>::iterator ite_end = population.end();
		for (; ite != ite_end; ++ite)
		{
			strcpy(ite->name, "none"); ite->energy = 50; ite->volume = 50; ite->method = 0; ite->formation_energy = 0;
		}
	}

	void Population::first_generation(bool R_flag)
	{
		FILE* fp_olog = fopen("olog", "at+");
		if (Generation != 0)
		{
			fprintf(fp_olog, " error:this is not first generation\n");
			fclose(fp_olog);
		}
		else
		{
			setvbuf(stdout, NULL, _IOLBF, 0);
			system("rm -rf ./poscar");
			mkdir("poscar", S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
			chdir("./poscar");
			int number = 0;
			if (!R_flag)
				gene_first_population();
			else
				gene_first_population_Random(this->init_atom, this->init_element, this->POPULATION_NUM, this->max_formula_atom, this->init_volume);
			DIR* dir_ptr;
			struct dirent* direntp;
			if ((dir_ptr = opendir(".")) == NULL)
			{
				printf("cannot open ./poscar \n");
				return;
			}
			else
			{
				while ((direntp = readdir(dir_ptr)) != NULL)
				{
					string str1(direntp->d_name);
					if (str1 == "." || str1 == "..") continue;
					number++;
				}
				closedir(dir_ptr);
			}
			chdir("..");
			Generation = 1;
			FILE* g = fopen("generation", "w");
			fprintf(g, "%d", Generation);
			fclose(g);
			if (number != POPULATION_NUM)
			{
				printf("The number of files in lib is not equal to POPULATION_NUM\n");
				exit(1);
			}
			FILE* N = fopen("pos_name", "w");
			for (int i = 0; i < POPULATION_NUM; i++)
			{
				sprintf(population[i].name, "1st_%d", i);
				fprintf(N, "%s\n", population[i].name);
			}
			fclose(N);
		}
		printf("Generate the first generation structures!\n");
	}

	inline void Population::Sort()
	{
		fprintf(stdout, "SORTING\n");
		sort(population.begin(), population.end(), comp);
		fprintf(stdout, "SORT DONE\n");
	}

	std::vector<unsigned int> readtype(const char* file, double vec[3][3], double xyz[MAX_NATOM][3])
	{
		FILE* fp = fopen(file, "r");
		char   title[MAX_FNAME];                   // title of the system
		double latt;                     // Scaling Lattice parameter 
		int    ifix;                     // Select 0 or normal 1
		int    iflg;                     // Direct 0 or Cartes 1
		int    nant[2];                      // nant[0];nant[1] number of atoms and types
		int    typenum[MAX_NELEM];                 // number of atoms for each ion type
		int    elemnum[MAX_NELEM];                 // atomic number of each element
		char   elemsym[MAX_NELEM][3];              // atomic symbol of each element
		char(*fix)[3];
		fix = new char[MAX_NATOM][3];                  // Selective fix on each ato
		readposcar(fp, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
		if (fp == NULL)
		{
			FILE* fp_compare = fopen("compare_olog", "at+");
			fprintf(fp_compare, "the %s is open wrong!\n", file);
			fclose(fp_compare);
		}
		std::vector<unsigned int>type;
		type.reserve(nant[0]);
		for (int i = 0; i < typenum[0]; i++)
			type.push_back(elemnum[0]);
		for (int j = typenum[0]; j < nant[0]; j++)
			type.push_back(elemnum[1]);
		fclose(fp);
		delete[] fix;
		return type;
	}

	int Population::COMPARE(char const* file1, char const* file2, double tolerance1, double tolerance2)
		//tolerance1 : Tolerance for position / length comparisons.Same units as cell matrix. (default = 0.05)
		//tolerance2 : Angular tolerance for lattice comparisons(default = 0.25 degrees) :
	{
		// file is empty
		FILE* fp1 = fopen(file1, "r");
		FILE* fp2 = fopen(file2, "r");
		char buf[1024];
		if (fgets(buf, 1024, fp1) == NULL || fgets(buf, 1024, fp2) == NULL)
			return 0;
		fclose(fp1);
		fclose(fp2);
		double(*xyz1)[3];
		xyz1 = new double[MAX_NATOM][3];
		double(*xyz2)[3];
		xyz2 = new double[MAX_NATOM][3];
		double cell1[3][3];
		double cell2[3][3];
		std::vector<unsigned int>type1 = readtype(file1, cell1, xyz1);
		std::vector<unsigned int>type2 = readtype(file2, cell2, xyz2);
		std::vector<XcVector> pos1;
		pos1.reserve(type1.size());
		for (int i = 0; i < type1.size(); i++)
		{
			pos1.push_back(XcVector(xyz1[i][0], xyz1[i][1], xyz1[i][2]));
		}
		std::vector<XcVector> pos2;
		pos2.reserve(type2.size());
		for (int i = 0; i < type2.size(); i++)
		{
			pos2.push_back(XcVector(xyz2[i][0], xyz2[i][1], xyz2[i][2]));
		}
		//bool match = XtalComp::compare(cell1, type1, pos1, cell2, type2, pos2, NULL, 0.05, 0.25);
		bool match = XtalComp::compare(cell1, type1, pos1, cell2, type2, pos2, NULL, tolerance1, tolerance2);
		delete[] xyz1;
		delete[] xyz2;
		if (!match)
			return 0;//the structures are different
		else
			return 1;//the structures are same
	}

	void Population::Pso(char current[], char last[], char pbest[], char gbest[], double w, int Generation, int num)
	{
		double(*xyz_pb)[3];
		char(*fix_pb)[3];
		char   title_pb[MAX_NCHAR];        // Title of system
		double latt_pb;                    // Scaling factor
		int    ifix_pb;                    // Select 0 or normal 1
		int    iflg_pb;                    // Direct 0 or Cartes 1
		int    nant_pb[2];                 // number of atoms and types
		int    typenum_pb[MAX_NELEM];      // number of atoms of each type
		int    elemnum_pb[MAX_NELEM];      // atomic number
		char   elemsym_pb[MAX_NELEM][3];   // atomic symbol
		double vec_pb[3][3];               // lattice vector
		xyz_pb = new double[MAX_NATOM][3];	// atomic coordinate
		fix_pb = new char[MAX_NATOM][3];       // Selective fix on each atom
		double(*xyz_gb)[3];
		char(*fix_gb)[3];
		char   title_gb[MAX_NCHAR];        // Title of system
		double latt_gb;                    // Scaling factor
		int    ifix_gb;                    // Select 0 or normal 1
		int    iflg_gb;                    // Direct 0 or Cartes 1
		int    nant_gb[2];                 // number of atoms and types
		int    typenum_gb[MAX_NELEM];      // number of atoms of each type
		int    elemnum_gb[MAX_NELEM];      // atomic number
		char   elemsym_gb[MAX_NELEM][3];   // atomic symbol
		double vec_gb[3][3];               // lattice vector
		xyz_gb = new double[MAX_NATOM][3];	// atomic coordinate
		fix_gb = new char[MAX_NATOM][3];       // Selective fix on each atom
		double(*xyz_last)[3];
		char(*fix_last)[3];
		char   title_last[MAX_NCHAR];        // Title of system
		double latt_last;                    // Scaling factor
		int    ifix_last;                    // Select 0 or normal 1
		int    iflg_last;                    // Direct 0 or Cartes 1
		int    nant_last[2];                 // number of atoms and types
		int    typenum_last[MAX_NELEM];      // number of atoms of each type
		int    elemnum_last[MAX_NELEM];      // atomic number
		char   elemsym_last[MAX_NELEM][3];   // atomic symbol
		double vec_last[3][3];               // lattice vector
		xyz_last = new double[MAX_NATOM][3];	// atomic coordinate
		fix_last = new char[MAX_NATOM][3];       // Selective fix on each atom
		double(*xyz_cur)[3];
		xyz_cur = new double[MAX_NATOM][3];	// atomic coordinate
		FILE* fp_pb = fopen(pbest, "r");
		readposcar(fp_pb, title_pb, latt_pb, ifix_pb, iflg_pb, nant_pb, typenum_pb, elemnum_pb, elemsym_pb, vec_pb, xyz_pb, fix_pb);
		fclose(fp_pb);
		FILE* fp_gb = fopen(gbest, "r");
		readposcar(fp_gb, title_gb, latt_gb, ifix_gb, iflg_gb, nant_gb, typenum_gb, elemnum_gb, elemsym_gb, vec_gb, xyz_gb, fix_gb);
		fclose(fp_gb);
		FILE* fp_last = fopen(last, "r");
		readposcar(fp_last, title_last, latt_last, ifix_last, iflg_last, nant_last, typenum_last, elemnum_last, elemsym_last, vec_last, xyz_last, fix_last);
		fclose(fp_last);
		const double c1 = 2;
		const double c2 = 2;
		double r1 = rand() % 10000 / 10000.1;
		double r2 = rand() % 10000 / 10000.1;
		double(*last_v)[3];
		last_v = new double[nant_last[0]][3];
		double(*current_v)[3];
		current_v = new double[nant_last[0]][3];
		char last_vf[100];
		char current_vf[100];
		sprintf(last_vf, "../velocity_%d", Generation - 1);
		sprintf(current_vf, "../velocity_buf");
		if (Generation == 1)
		{
			for (int i = 0; i < nant_last[0]; i++)
				for (int j = 0; j < 3; j++)
					last_v[i][j] = rand() % 10000 / 10001.0;
		}
		else
		{
			FILE* fp = fopen(last_vf, "r");
			char line[1024];
			int start_line = num * nant_last[0];
			int cnt = 0;
			int iat = 0;
			while (!feof(fp))
			{
				fgets(line, 1024, fp);
				if (cnt == start_line)
				{
					sscanf(line, "%lf %lf %lf", &last_v[iat][0], &last_v[iat][1], &last_v[iat][2]);
					if (iat < nant_last[0] - 1)
					{
						start_line++;
						iat++;
					}
				}
				cnt++;
			}
			fclose(fp);
		}
		for (int i = 0; i < nant_last[0]; i++)
		{
			for (int j = 0; j < 3; j++)
			{
				current_v[i][j] = w * last_v[i][j] + c1 * r1 * (xyz_pb[i][j] - xyz_last[i][j]) + c2 * r2 * (xyz_gb[i][j] - xyz_last[i][j]);
				xyz_cur[i][j] = xyz_last[i][j] + current_v[i][j];
			}
		}
		FILE* fp1 = fopen(current_vf, "w");
		for (int i = 0; i < nant_last[0]; i++)
			fprintf(fp1, "%lf %lf %lf\n", current_v[i][0], current_v[i][1], current_v[i][2]);
		fclose(fp1);
		strcpy(title_last, current);
		FILE* fp_current = fopen(current, "w");
		savposcar(fp_current, title_last, latt_last, ifix_last, iflg_last, nant_last, typenum_last, elemnum_last, elemsym_last, vec_last, xyz_cur, fix_last);
		fclose(fp_current);
		delete[] xyz_cur;
		delete[] xyz_gb;
		delete[] xyz_last;
		delete[] xyz_pb;
		delete[] fix_gb;
		delete[] fix_last;
		delete[] fix_pb;
		delete[] last_v;
		delete[] current_v;
	}

	void box_fix(const char file[], double layer)
	{
		FILE* fp = fopen(file, "r");
		char   title[MAX_NCHAR];        // Title of system
		double latt;                    // Scaling factor
		int    ifix;                    // Select 0 or normal 1
		int    iflg;                    // Direct 0 or Cartes 1
		int    nant[2];                 // number of atoms and types
		int    typenum[MAX_NELEM];      // number of atoms of each type
		int    elemnum[MAX_NELEM];      // atomic number
		char   elemsym[MAX_NELEM][3];   // atomic symbol
		double vec[3][3];               // lattice vector
		double(*xyz)[3] = new double[MAX_NATOM][3];  // atomic coordinate
		char(*fix)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom
		readposcar(fp, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
		fclose(fp);
		direct_to_carts(iflg, vec, xyz, nant);
		double x_max = INT_MIN, x_min = INT_MAX, y_max = INT_MIN, y_min = INT_MAX, z_max = INT_MIN, z_min = INT_MAX;
		double n_vec[3][3] = { 0 };
		for (int i = 0; i < nant[0]; i++)
		{
			x_max = max(x_max, xyz[i][0]);
			y_max = max(y_max, xyz[i][1]);
			z_max = max(z_max, xyz[i][2]);
			x_min = min(x_min, xyz[i][0]);
			y_min = min(y_min, xyz[i][1]);
			z_min = min(z_min, xyz[i][2]);
		}
		for (int i = 0; i < nant[0]; i++)
		{
			xyz[i][0] -= x_min - layer;
			xyz[i][1] -= y_min - layer;
			xyz[i][2] -= z_min - layer;
		}
		n_vec[0][0] = x_max - x_min + 2 * layer;
		n_vec[1][1] = y_max - y_min + 2 * layer;
		n_vec[2][2] = z_max - z_min + 2 * layer;
		carts_to_direct(iflg, n_vec, xyz, nant);
		FILE* fp1 = fopen(file, "w");
		savposcar(fp1, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, n_vec, xyz, fix);
		fclose(fp1);
		delete[] xyz;
		delete[] fix;
	}

	void Population::Renew(char label[])
	{
		// get current generation
		FILE* fp2 = fopen("generation", "r");
		fscanf(fp2, "%d", &Generation);
		fclose(fp2);
		//sort induvals from the  most stable (energy mininum ) to the most unstable
		//creat parent pool we chose a member probability of (Emax-E)/(Emax-Emin) 
		//read energy and volume -> population[i].energy population[i].volume
		FILE* fp = fopen("V_E", "r");
		char buf[1024];
		for (int j = 0; j < POPULATION_NUM; j++)
		{
			fgets(buf, 1024, fp);
			sscanf(buf, "%lf	%lf	%lf", &population[j].energy, &population[j].formation_energy, &population[j].volume);
			//save stable structure to pool
			if (population[j].formation_energy < 0)
			{
				char target_file[1024], source_file[1024];
				sprintf(source_file, "./generation_%d/%d/CONTCAR", Generation, j);
				sprintf(target_file, "./structure_pool/%d_%d", Generation, j);
				copy_file(target_file, source_file);
			}
		}
		fclose(fp);
		FILE* fp1 = fopen("V_E", "w");
		fclose(fp1);
		// read current posname -> population[i].name
		FILE* fp4 = fopen("pos_name", "r");
		char name2[POPULATION_NUM][20];
		for (int i = 0; i < POPULATION_NUM; i++)
		{
			fscanf(fp4, "%s", name2[i]);
			strcpy(population[i].name, name2[i]);
		}
		fclose(fp4);
		// PSO algorithm
		vector<Account> pbest(POPULATION_NUM);
		double gbest_energy;
		char gbest_name[20];
		double pbest_energy;
		char pbest_name[20];
		if (!strcmp(label, "PSO"))
		{
			if (Generation == 1)
			{
				FILE* fp_pbe = fopen("pbest", "w");
				for (int i = 0; i < POPULATION_NUM; i++)
				{
					fprintf(fp_pbe, "%s %lf\n", population[i].name, population[i].energy);
					pbest[i].energy = population[i].energy;
					strcpy(pbest[i].name, population[i].name);
				}
				fclose(fp_pbe);
			}
			else
			{
				FILE* fp_pbe = fopen("pbest", "r");
				for (int i = 0; i < POPULATION_NUM; i++)
				{
					fscanf(fp_pbe, "%s %lf", pbest_name, &pbest_energy);
					if (population[i].energy < pbest_energy)
					{
						pbest[i].energy = population[i].energy;
						strcpy(pbest[i].name, population[i].name);
					}
					else
					{
						pbest[i].energy = pbest_energy;
						strcpy(pbest[i].name, pbest_name);
					}
				}
				fclose(fp_pbe);
				FILE* fp_pbe2 = fopen("pbest", "w");
				for (int i = 0; i < POPULATION_NUM; i++)
					fprintf(fp_pbe2, "%s %lf\n", pbest[i].name, pbest[i].energy);
				fclose(fp_pbe2);
			}
			Sort();
			if (Generation == 1)
			{
				FILE* fp_gbe = fopen("gbest", "w");
				fprintf(fp_gbe, "%s %lf\n", population[0].name, population[0].energy);
				fclose(fp_gbe);
				gbest_energy = population[0].energy;
				strcpy(gbest_name, population[0].name);
			}
			else
			{
				FILE* fp_gbe = fopen("gbest", "r");
				fscanf(fp_gbe, "%s%lf", gbest_name, &gbest_energy);
				fclose(fp_gbe);
				FILE* fp_gbe2 = fopen("gbest", "w");
				if (population[0].energy < gbest_energy)
				{
					gbest_energy = population[0].energy;
					strcpy(gbest_name, population[0].name);
				}
				fprintf(fp_gbe2, "%s %lf\n", gbest_name, gbest_energy);
				fclose(fp_gbe2);
			}
		}
		else
			Sort();
		//record all structures to exclude similarity strucutres
		vector<string> compare_name;
		compare_name.reserve((Generation + 1) * POPULATION_NUM);
		char name[20];
		for (int j = 1; j < Generation + 1; j++)
		{
			for (int i = 0; i < POPULATION_NUM; i++)
			{
				if (j == 1)
					sprintf(name, "%dst_%d", j, i);
				else
					sprintf(name, "%dth_%d", j, i);
				compare_name.push_back(name);
			}
		}
		//update Generation
		Generation = Generation + 1;
		FILE* fp3 = fopen("generation", "w");
		fprintf(fp3, "%d", Generation);
		fclose(fp3);
		chdir("./poscar");
		char POSCAR_Child[20];
		int i = 0;
		int cnt = 5;//the num of generate similarity structure times
		if (!strcmp(label, "EA"))
		{
			//printf min energy and volume
			double Emin = population[0].energy;
			double Emax = population[POPULATION_NUM - 1].energy;
			FILE* fp_min = fopen("../min_energy_structure", "at+");
			fprintf(fp_min, "generation: %d\n  name: %s  energy:  %lf   volume:  %lf\n", Generation - 1, population[0].name, population[0].energy, population[0].volume);
			fclose(fp_min);
			while (i < POPULATION_NUM)
			{
			again:
				int elite = max(POPULATION_NUM / 10, 1);
				int R = 0;
				char olog[1024] = { NULL };
				sprintf(olog, "../olog_%d", Generation);
				if (i == 0)
				{
					FILE* fp_title = fopen(olog, "at+");
					fprintf(fp_title, "**************************Generation = %d**************************\n", Generation);
					fclose(fp_title);
				}
				sprintf(POSCAR_Child, "%dth_%d", Generation, i);
				// We choose 1/10 of the lowest energy individuals as elite individuals to be retained to the next generation
				if (i < elite)
				{
					copy_file(POSCAR_Child, population[i].name);
					FILE* _fp1 = fopen(population[i].name, "r");
					POSCAR pos;
					readposcar(_fp1, pos);
					fclose(_fp1);
					strcpy(pos.title, POSCAR_Child);
					FILE* _fp2 = fopen(population[i].name, "w");
					savposcar(_fp2, pos);
					fclose(_fp2);
					FILE* fp_elite = fopen(olog, "at+");
					fprintf(fp_elite, "**********************************\n");
					fprintf(fp_elite, "%s is the elite inherited from %s\n", POSCAR_Child, population[i].name);
					strcpy(population_child[i].name, POSCAR_Child);
					fclose(fp_elite);
					i++;
					continue;
				}
				//We will eliminate the structure that highest(1/10)
				else if (i >= POPULATION_NUM - elite)
				{
					gene_Random_structrue(init_atom, init_element, max_formula_atom, POSCAR_Child, init_volume);
					FILE* fp_elim = fopen(olog, "at+");
					fprintf(fp_elim, "**********************************\n");
					fprintf(fp_elim, "%s is eliminated as unqualified, %s will be generated randomly!\n", population[i].name, POSCAR_Child);
					fclose(fp_elim);
				}
				else
				{
					// we chose method randomly to change and get new structues we can chose method another way 
					int method = rand() % 10;//???? the ratio
					int k1 = rand() % (POPULATION_NUM - elite);
					int k2 = rand() % (POPULATION_NUM - elite);
					double r = rand() % 10000 / 10000.1;
					while (r > (Emax - population[k1].energy) / (Emax - Emin) && population[k1].formation_energy < 0)
					{
						r = rand() % 10000 / 10000.1;
						k1 = rand() % POPULATION_NUM;
					}
					while (r > (Emax - population[k2].energy) / (Emax - Emin) && population[k1].formation_energy < 0)
					{
						r = rand() % 10000 / 10000.1;
						k2 = rand() % POPULATION_NUM;
					}
					if (method < 6)
					{
						FILE* variant = fopen(olog, "at+");
						fprintf(variant, "**********************************\n");
						fprintf(variant, ">>>parent of %s is %s\n", POSCAR_Child, population[k1].name);
						fclose(variant);
						if (!Variant_singer_parent(population[k1].name, POSCAR_Child, method, Generation)) // parent file is empty
						{
							FILE* fp_gen = fopen(olog, "at+");
							fprintf(fp_gen, "%s may be an empty file!, %s will be generated randomly!\n", population[k1].name, POSCAR_Child);
							fclose(fp_gen);
							gene_Random_structrue(init_atom, init_element, max_formula_atom, POSCAR_Child, init_volume);
						}
					}
					else
					{
						FILE* variant = fopen(olog, "at+");
						fprintf(variant, "**********************************\n");
						fprintf(variant, ">>>parent of %s is %s and %s\n", POSCAR_Child, population[k1].name, population[k2].name);
						fclose(variant);
						if (Variant_two_parent(population[k1].name, population[k2].name, POSCAR_Child, Generation))// parent file is empty
						{
							FILE* fp_gen = fopen(olog, "at+");
							fprintf(fp_gen, "%s or %s may be an empty file!, %s will be generated randomly!\n", population[k1].name, population[k2].name, POSCAR_Child);
							fclose(fp_gen);
							gene_Random_structrue(init_atom, init_element, max_formula_atom, POSCAR_Child, init_volume);
						}
					}
				}
				strcpy(population_child[i].name, POSCAR_Child);
				for (int k = 0; k < compare_name.size(); k++)
				{
					R += COMPARE(compare_name[k].c_str(), population_child[i].name, 0.05, 0.25);
					if (R != 0 && cnt != 0)
					{
						FILE* variant = fopen(olog, "at+");
						fprintf(variant, ">>>the %s is repeated with %s\n", population_child[i].name, compare_name[k].c_str());
						fclose(variant);
						cnt--;
						goto again;
					}
				}
				compare_name.push_back(population_child[i].name);
				i++;
			}
		}
		else if (!strcmp(label, "PSO"))
		{
			char current_vf[50];
			while (i < POPULATION_NUM)
			{
			again2:
				int R = 0;
				sprintf(POSCAR_Child, "%dth_%d", Generation, i);
				double w = 0.9 - (0.9 - 0.4) / maxiterator * (Generation - 1);
				Pso(POSCAR_Child, name2[i], pbest[i].name, gbest_name, w, Generation - 1, i);
				strcpy(population_child[i].name, POSCAR_Child);
				for (int k = 0; k < compare_name.size(); k++)
				{
					R += COMPARE(compare_name[k].c_str(), population_child[i].name, 0.05, 0.25);
					if (R != 0 && cnt != 0)
					{
						cnt--;
						goto again2;
					}
				}
				sprintf(current_vf, "../velocity_%d", Generation - 1);
				copy_file(current_vf, "../velocity_buf");
				compare_name.push_back(population_child[i].name);
				i++;
			}
		}
		for (int j = 0; j < POPULATION_NUM; j++)
		{
			strcpy(population[j].name, population_child[j].name);
		}
		chdir("..");
		FILE* fp5 = fopen("pos_name", "w");
		for (int j = 0; j < POPULATION_NUM; j++)
			fprintf(fp5, "%s\n", population[j].name);
		fclose(fp5);
	}

	void Population::JudgeVaspComplete()
	{
		timeval start_time, end_time;
		int i = 0;
		while (i < POPULATION_NUM)
		{
			int reflag = 0, end_flag = 0;
			char dir[10];
			sprintf(dir, "./%d", i);
			chdir(dir);
			FILE* fp = fopen("OUTCAR", "r");
			if (fp == NULL)
			{
				if (i != 0)
					reflag = 1;
				chdir("../");
				continue;
			}
			else if (reflag == 1 || i == 0)
				gettimeofday(&start_time, 0);
			char buf[1024];
			while (fgets(buf, 1024, fp) != NULL)
				if (strstr(buf, "Voluntary context switches"))
					end_flag = 1;
			gettimeofday(&end_time, 0);
			double timeuse = 1000000 * (end_time.tv_sec - start_time.tv_sec) + end_time.tv_usec - start_time.tv_usec;
			if (end_flag == 0 && timeuse / 1000000 / 60 / 60 <= this->maxcaltime)
			{
				chdir("../");
				fclose(fp);
				continue;
			}
			chdir("../");
			fclose(fp);
			i++;
		}
	}

	void Population::derive()
	{
		FILE* fp = fopen("../../V_E", "at+");
		double energy = get_energy();
		double volume = get_volume();
		//calculate formation energy 
		double formation = energy;
		int total_num = 0;
		for (int i = 0; i < this->init_atom.size(); i++)
		{
			formation -= single_energy[i] * init_atom[i];
			total_num += init_atom[i];
		}
		formation /= total_num;
		fprintf(fp, "%lf  %lf	%lf\n", energy, formation, volume);
		fclose(fp);
		remove("CHG");
		remove("DOSCAR");
		remove("EIGENVAL");
		remove("PCDAT");
		remove("vasprun.xml");
		remove("WAVECAR");
		remove("XDATCAR");
		remove("IBZKPT");
		remove("PROCAR");
	}

	void Population::Record()
	{
		FILE* fp_olog = fopen("olog", "at+");
		FILE* G = fopen("generation", "r");
		fscanf(G, "%d", &Generation);
		FILE* fp = fopen("Record", "at+");
		if (fp == NULL)
		{
			fprintf(fp_olog, "Record cannot open\n");
			return;
		}
		fprintf(fp, "generation:	%d \n", Generation);
		//		population.name
		FILE* fp1 = fopen("pos_name", "r");
		char name2[POPULATION_NUM];
		for (int i = 0; i < POPULATION_NUM; i++)
		{
			fscanf(fp1, "%s", &name2[i]);
			strcpy(population[i].name, &name2[i]);
		}
		fclose(fp1);
		// 		population.energy and population.volume
		FILE* fp2 = fopen("V_E", "r");
		char buf[1024];
		for (int j = 0; j < POPULATION_NUM; j++)
		{
			fgets(buf, 1024, fp2);
			sscanf(buf, "%lf	%lf	%lf", &population[j].energy, &population[j].formation_energy, &population[j].volume);
		}
		fclose(fp2);
		for (int k = 0; k < POPULATION_NUM; k++)
			fprintf(fp, "name: %s number:%d  energy: %lfev formation_energy: %lfev	volume: %lf\n", population[k].name, k, population[k].energy, population[k].formation_energy, population[k].volume);
		fclose(fp);
		fclose(fp_olog);
	}

	void Population::control()
	{
		FILE* fp_con = fopen("control.txt", "r");
		if (fp_con == NULL)
		{
			FILE* fp_con = fopen("command_olog", "at+");
			fprintf(fp_con, "the control file is wrong!\n");
			fclose(fp_con);
		}
		//count lines
		char buf[MAX_NLINE][MAX_NTEXT];
		strcpy(buf[0], "EVOPS");
		int lines = 0;
		while (fgets(buf[lines + 1], MAX_NTEXT, fp_con) != NULL)
		{
			int j = 0;
			while (1)
			{
				if (buf[lines + 1][j] == ' ')
				{
					buf[lines + 1][j] = '\0';
					break;
				}
				else
					j++;
			}
			lines++;
		}
		char* TEXT[MAX_NLINE];
		for (int i = 0; i < MAX_NLINE; i++)
		{
			TEXT[i] = buf[i];
		}
		fclose(fp_con);
		if (!strcmp("mf", buf[1]))
			MultiFile_Variant(lines + 1, TEXT);
		else
			Command(lines + 1, TEXT);
	}

	void Population::Command(int argc, char* argv[])
	{
		system("rm ../buf");
		if (!strcmp("Exchange", argv[1]) || !strcmp("Strain", argv[1]) || !strcmp("Ripple", argv[1]) ||
			!strcmp("Randmove", argv[1]) || !strcmp("Slip", argv[1]) || !strcmp("Twist", argv[1]) ||
			!strcmp("Rotation", argv[1]) || !strcmp("Stripple", argv[1]) ||
			!strcmp("Permustrain", argv[1]) || !strcmp("Rotstrain", argv[1]))
		{
			char   title[MAX_NCHAR];        // Title of system
			double latt;                    // Scaling factor
			int    ifix;                    // Select 0 or normal 1
			int    iflg;                    // Direct 0 or Cartes 1
			int    nant[2];                 // number of atoms and types
			int    typenum[MAX_NELEM];      // number of atoms of each type
			int    elemnum[MAX_NELEM];      // atomic number
			char   elemsym[MAX_NELEM][3];   // atomic symbol
			double vec[3][3];               // lattice vector
			double(*xyz)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*fix)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom
			char   title_e[MAX_NCHAR];        // Title of system
			double latt_e;                    // Scaling factor
			int    ifix_e;                    // Select 0 or normal 1
			int    iflg_e;                    // Direct 0 or Cartes 1
			int    nant_e[2];                 // number of atoms and types
			int    typenum_e[MAX_NELEM];      // number of atoms of each type
			int    elemnum_e[MAX_NELEM];      // atomic number
			char   elemsym_e[MAX_NELEM][3];   // atomic symbol
			double vec_e[3][3];               // lattice vector
			double(*xyz_e)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*fix_e)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom
			FILE* fp_olog = fopen("command_olog", "at+");
			fprintf(fp_olog, "*************************************\n");
			fclose(fp_olog);
			FILE* fp_parent;
			fp_parent = fopen(argv[2], "r");
			if (fp_parent == NULL)
			{
				FILE* fp_olog = fopen("command_olog", "at+");
				fprintf(fp_olog, "the %s is open wrong\n", argv[2]);
				fclose(fp_olog);
			}
			readposcar(fp_parent, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
			fclose(fp_parent);
			FILE* fp_olog_ = fopen("command_olog", "at+");
			fprintf(fp_olog_, ">>>the parent file is %s\n", argv[2]);
			fclose(fp_olog_);
			if (!strcmp("Exchange", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Exchange(nant, typenum, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < (atoi(argv[4]) * 2 + 5))
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					int* Num1 = (int*)malloc(sizeof(int) * atoi(argv[4]));
					int* Num2 = (int*)malloc(sizeof(int) * atoi(argv[4]));
					for (int i = 0; i < atoi(argv[4]); i++)
					{
						Num1[i] = atoi(argv[5 + 2 * i]);
						Num2[i] = atoi(argv[6 + 2 * i]);
					}
					if (atoi(argv[4]) > nant[0] / 2)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, ">>>warning! Exchange atoms too many times!");
						fclose(fp_olog);
					}
					Exchange_command(nant, typenum, xyz, xyz_e, atoi(argv[4]), Num1, Num2);
					free(Num1);
					free(Num2);
					Num1 = NULL;
					Num2 = NULL;
				}
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
			}
			if (!strcmp("Strain", argv[1]))
			{
				char mode[5];
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					int RandE_Strain_Mode = rand() % 7;
					switch (RandE_Strain_Mode)
					{
					case 0:strcpy(mode, "xx"); break;
					case 1:strcpy(mode, "yy"); break;
					case 2:strcpy(mode, "zz"); break;
					case 3:strcpy(mode, "xy"); break;
					case 4:strcpy(mode, "xz"); break;
					case 5:strcpy(mode, "yz"); break;
					case 6:strcpy(mode, "xyz"); break;
					default:break;
					}
					Strain(mode, vec, vec_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					strcpy(mode, argv[4]);
					if ((!strcmp(mode, "xx") && argc < 7) || (!strcmp(mode, "yy") && argc < 7) ||
						(!strcmp(mode, "zz") && argc < 7) || (!strcmp(mode, "xy") && argc < 7) ||
						(!strcmp(mode, "xz") && argc < 7) || (!strcmp(mode, "yz") && argc < 7) ||
						(!strcmp(mode, "xyz") && argc < 12))
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					if (!strcmp(mode, "xx"))
						Strain_command(mode, vec, vec_e, atof(argv[5]), 0, 0, 0, 0, 0);
					if (!strcmp(mode, "yy"))
						Strain_command(mode, vec, vec_e, 0, atof(argv[5]), 0, 0, 0, 0);
					if (!strcmp(mode, "zz"))
						Strain_command(mode, vec, vec_e, 0, 0, atof(argv[5]), 0, 0, 0);
					if (!strcmp(mode, "xy"))
						Strain_command(mode, vec, vec_e, 0, 0, 0, atof(argv[5]), 0, 0);
					if (!strcmp(mode, "xz"))
						Strain_command(mode, vec, vec_e, 0, 0, 0, 0, atof(argv[5]), 0);
					if (!strcmp(mode, "yz"))
						Strain_command(mode, vec, vec_e, 0, 0, 0, 0, 0, atof(argv[5]));
					if (!strcmp(mode, "xyz"))
						Strain_command(mode, vec, vec_e, atof(argv[5]), atof(argv[6]),
							atof(argv[7]), atof(argv[8]), atof(argv[9]), atof(argv[10]));
				}
				double v = volume(vec);
				double v_e = volume(vec_e);
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec_e[i][j] *= pow(v / v_e, 1 / 3);
					}
				}
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec[i][j] = vec_e[i][j];
					}
				}
			}
			if (!strcmp("Ripple", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Ripple(nant, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < 8)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Ripple_command(nant, xyz, xyz_e, atoi(argv[4]), atof(argv[5]), atof(argv[6]));
					//The reference value of temp is as follows
					/*case 0: x = 0, y = 1, z = 2; axis = 'z'; break;
					case 1: x = 0, y = 2, z = 1; axis = 'y'; break;
					case 2: x = 1, y = 0, z = 2; axis = 'z'; break;
					case 3: x = 1, y = 2, z = 0; axis = 'x'; break;
					case 4: x = 2, y = 0, z = 1; axis = 'y'; break;
					case 5: x = 2, y = 1, z = 0; axis = 'x'; break;*/
				}
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
			}
			if (!strcmp("Randmove", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Randmove(nant, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				if (argc < 7)
				{
					FILE* fp_olog = fopen("command_olog", "at+");
					fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
					fclose(fp_olog);
				}
				Randmove_command(nant, xyz, xyz_e, atoi(argv[4]), atof(argv[5]));
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
			}
			if (!strcmp("Slip", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Slip(nant, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < 11)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Slip_command(nant, xyz, xyz_e, atof(argv[4]), atof(argv[5]),
						atof(argv[6]), atof(argv[7]), atof(argv[8]), atof(argv[9]));
				}
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
			}
			if (!strcmp("Twist", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Twist(nant, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < 7)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Twist_command(nant, xyz, xyz_e, argv[4][0], atof(argv[5]));
				}
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
			}
			if (!strcmp("Rotation", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Rotation(vec, vec_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < 8)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Rotation_command(vec, vec_e, atof(argv[4]), atof(argv[5]), atof(argv[6]));
				}
				double v = volume(vec);
				double v_e = volume(vec_e);
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec_e[i][j] *= pow(v / v_e, 1 / 3);
					}
				}
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec[i][j] = vec_e[i][j];
					}
				}
			}
			if (!strcmp("Stripple", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Stripple(nant, vec, vec_e, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < 14)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Stripple_command(nant, vec, vec_e, xyz, xyz_e, atof(argv[4]), atof(argv[5]),
						atof(argv[6]), atof(argv[7]), atof(argv[8]), atof(argv[9]),
						atoi(argv[10]), atof(argv[11]), atof(argv[12]));
				}
				double v = volume(vec);
				double v_e = volume(vec_e);
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec_e[i][j] *= pow(v / v_e, 1 / 3);
					}
				}
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec[i][j] = vec_e[i][j];
					}
				}
			}
			if (!strcmp("Permustrain", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Permustrain(nant, typenum, vec, vec_e, xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < (atoi(argv[10]) * 2 + 11))
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					int* Num1 = (int*)malloc(sizeof(int) * atoi(argv[10]));
					int* Num2 = (int*)malloc(sizeof(int) * atoi(argv[10]));
					for (int i = 0; i < atoi(argv[10]); i++)
					{
						Num1[i] = atoi(argv[11 + 2 * i]);
						Num2[i] = atoi(argv[12 + 2 * i]);
					}
					if (atoi(argv[4]) > nant[0] / 2)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, ">>>warning! Exchange atoms too many times!");
						fclose(fp_olog);
					}
					Permustrain_command(nant, typenum, vec, vec_e, xyz, xyz_e, atof(argv[4]), atof(argv[5]),
						atof(argv[6]), atof(argv[7]), atof(argv[8]),
						atof(argv[9]), atoi(argv[10]), Num1, Num2);
					free(Num1);
					free(Num2);
					Num1 = NULL;
					Num2 = NULL;
				}
				double v = volume(vec);
				double v_e = volume(vec_e);
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec_e[i][j] *= pow(v / v_e, 1 / 3);
					}
				}
				for (int i = 0; i < nant[0]; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						xyz[i][j] = xyz_e[i][j];
					}
				}
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec[i][j] = vec_e[i][j];
					}
				}
			}
			if (!strcmp("Rotstrain", argv[1]))
			{
				if (!strcmp("Random", argv[4]) || 'R' == argv[4][0])
				{
					Rotstrain(xyz, xyz_e);
					system("cat ../buf >> command_olog");
					system("rm ../buf");
				}
				else
				{
					if (argc < 10)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Rotstrain_command(vec, vec_e, atof(argv[4]), atof(argv[5]),
						atof(argv[6]), argv[7], atof(argv[8]));
				}
				double v = volume(vec);
				double v_e = volume(vec_e);
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec_e[i][j] *= pow(v / v_e, 1 / 3);
					}
				}
				for (int i = 0; i < 3; i++)
				{
					for (int j = 0; j < 3; j++)
					{
						vec[i][j] = vec_e[i][j];
					}
				}
			}
			FILE* fp_child;
			fp_child = fopen(argv[3], "w");
			if (fp_child == NULL)
			{
				FILE* fp_olog = fopen("command_olog", "at+");
				fprintf(fp_olog, "the %s is open wrong\n", argv[3]);
				fclose(fp_olog);
			}
			FILE* fp_olog2 = fopen("command_olog", "at+");
			fprintf(fp_olog2, ">>>the child file is %s\n", argv[3]);
			fprintf(fp_olog2, "*************************************\n");
			fclose(fp_olog2);
			strcpy(title, argv[3]);
			savposcar(fp_child, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
			fclose(fp_child);
			delete[] xyz;
			delete[] fix;
		}
		else if (!strcmp("Crossover", argv[1]))
		{
			char   title[MAX_NCHAR];        // Title of system
			double latt;                    // Scaling factor
			int    ifix;                    // Select 0 or normal 1
			int    iflg;                    // Direct 0 or Cartes 1
			int    nant[2];                 // number of atoms and types
			int    typenum[MAX_NELEM];      // number of atoms of each type
			int    elemnum[MAX_NELEM];      // atomic number
			char   elemsym[MAX_NELEM][3];   // atomic symbol
			double vec[3][3];               // lattice vector
			double(*xyz)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*fix)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom

			char   _title[MAX_NCHAR];        // Title of system
			double _latt;                    // Scaling factor
			int    _ifix;                    // Select 0 or normal 1
			int    _iflg;                    // Direct 0 or Cartes 1
			int    _nant[2];                 // number of atoms and types
			int    _typenum[MAX_NELEM];      // number of atoms of each type
			int    _elemnum[MAX_NELEM];      // atomic number
			char   _elemsym[MAX_NELEM][3];   // atomic symbol
			double _vec[3][3];               // lattice vector
			double(*_xyz)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*_fix)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom

			char   title_e[MAX_NCHAR];        // Title of system
			double latt_e;                    // Scaling factor
			int    ifix_e;                    // Select 0 or normal 1
			int    iflg_e;                    // Direct 0 or Cartes 1
			int    nant_e[2];                 // number of atoms and types
			int    typenum_e[MAX_NELEM];      // number of atoms of each type
			int    elemnum_e[MAX_NELEM];      // atomic number
			char   elemsym_e[MAX_NELEM][3];   // atomic symbol
			double vec_e[3][3];               // lattice vector
			double(*xyz_e)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*fix_e)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom

			FILE* fp_olog = fopen("command_olog", "at+");
			fprintf(fp_olog, "*************************************\n");
			fclose(fp_olog);
			FILE* fp_parent;
			fp_parent = fopen(argv[2], "r");
			if (fp_parent == NULL)
			{
				FILE* fp_olog = fopen("command_olog", "at+");
				fprintf(fp_olog, "the %s is open wrong\n", argv[2]);
				fclose(fp_olog);
			}
			readposcar(fp_parent, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
			fclose(fp_parent);
			FILE* _fp_parent;
			_fp_parent = fopen(argv[3], "r");
			if (_fp_parent == NULL)
			{
				FILE* fp_olog = fopen("command_olog", "at+");
				fprintf(fp_olog, "the %s is open wrong\n", argv[3]);
				fclose(fp_olog);
			}
			readposcar(_fp_parent, _title, _latt, _ifix, _iflg, _nant, _typenum, _elemnum, _elemsym, _vec, _xyz, _fix);
			fclose(_fp_parent);
			if (argc < 7)
			{
				FILE* fp_olog = fopen("command_olog", "at+");
				fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
				fclose(fp_olog);
			}
			FILE* fp_olog_ = fopen("command_olog", "at+");
			fprintf(fp_olog_, ">>>the parent file is %s and %s\n", argv[2], argv[3]);
			fclose(fp_olog_);
			if (nant[0] != _nant[0])return;
			if (!strcmp("Random", argv[5]) || 'R' == argv[5][0])
			{
				Crossover(latt, _latt, latt_e, nant, typenum, vec, _vec, vec_e, xyz, _xyz, xyz_e);
				system("cat ../buf >> command_olog");
				system("rm ../buf");
			}
			else
				Crossover_command(latt, _latt, latt_e, nant, typenum, vec, _vec, vec_e, xyz, _xyz, xyz_e, argv[5][0], atof(argv[6]));
			strcpy(title_e, argv[4]);
			FILE* fp_child;
			fp_child = fopen(argv[4], "w");
			savposcar(fp_child, title_e, latt_e, ifix, iflg, nant, typenum, elemnum, elemsym, vec_e, xyz_e, fix);
			fclose(fp_child);
			FILE* fp_olog2 = fopen("command_olog", "at+");
			fprintf(fp_olog2, ">>>the child file is %s\n", argv[4]);
			fprintf(fp_olog2, "*************************************\n");
			fclose(fp_olog2);
			delete[] xyz;
			delete[] fix;
			delete[] _xyz;
			delete[] _fix;
			delete[] xyz_e;
			delete[] fix_e;
		}
	}

	void Population::MultiFile_Variant(int argc, char* argv[])
	{
		system("rm ../buf");
		char(*OpenFileName)[MAX_FNAME] = new char[MAX_NATOM][MAX_FNAME];
		char(*SaveFileName)[MAX_FNAME] = new char[MAX_NATOM][MAX_FNAME];
		char(*buf)[MAX_FNAME] = new char[MAX_NATOM][MAX_FNAME];
		for (int i = 0; i < atoi(argv[3]); i++)
		{
			sprintf(buf[i], "%s%d", argv[2], i + 1);
		}
	again:
		for (int i = 0; i < atoi(argv[3]); i++)
		{
			sprintf(OpenFileName[i], "%s%d", argv[2], i + 1);
			FILE* FP = fopen("command_olog", "at+");
			fprintf(FP, "*************************************\n");
			fclose(FP);
			FILE* fp_parent = fopen(OpenFileName[i], "r");
			if (fp_parent == NULL)
			{
				FILE* FP_olog = fopen("command_olog", "at+");
				fprintf(FP_olog, "the %s is open wrong\n", OpenFileName[i]);
				fclose(FP_olog);
			}
			FILE* fp_olog_ = fopen("command_olog", "at+");
			fprintf(fp_olog_, ">>>the parent file is %s\n", OpenFileName[i]);
			fclose(fp_olog_);
			char   title[MAX_NCHAR];        // Title of system
			double latt;                    // Scaling factor
			int    ifix;                    // Select 0 or normal 1
			int    iflg;                    // Direct 0 or Cartes 1
			int    nant[2];                 // number of atoms and types
			int    typenum[MAX_NELEM];      // number of atoms of each type
			int    elemnum[MAX_NELEM];      // atomic number
			char   elemsym[MAX_NELEM][3];   // atomic symbol
			double vec[3][3];               // lattice vector
			double(*xyz)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*fix)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom
			char   title_e[MAX_NCHAR];        // Title of system
			double latt_e;                    // Scaling factor
			int    ifix_e;                    // Select 0 or normal 1
			int    iflg_e;                    // Direct 0 or Cartes 1
			int    nant_e[2];                 // number of atoms and types
			int    typenum_e[MAX_NELEM];      // number of atoms of each type
			int    elemnum_e[MAX_NELEM];      // atomic number
			char   elemsym_e[MAX_NELEM][3];   // atomic symbol
			double vec_e[3][3];               // lattice vector
			double(*xyz_e)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*fix_e)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom

			char   _title[MAX_NCHAR];        // Title of system
			double _latt;                    // Scaling factor
			int    _ifix;                    // Select 0 or normal 1
			int    _iflg;                    // Direct 0 or Cartes 1
			int    _nant[2];                 // number of atoms and types
			int    _typenum[MAX_NELEM];      // number of atoms of each type
			int    _elemnum[MAX_NELEM];      // atomic number
			char   _elemsym[MAX_NELEM][3];   // atomic symbol
			double _vec[3][3];               // lattice vector
			double(*_xyz)[3] = new double[MAX_NATOM][3];  // atomic coordinate
			char(*_fix)[3] = new char[MAX_NATOM][3];  // Selective fix on each atom

			readposcar(fp_parent, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
			fclose(fp_parent);
			int MAX_LENGTH = 5;//EVOPS MF INPUT N OUTPUT + VARIANT + PARAMENT + ...
			for (int k = 5; k < argc; k++)
			{
				if (!strcmp("Exchange", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Exchange(nant, typenum, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += atoi(argv[k + 1]) * 2 + 2;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[i]);
							fclose(fp_olog);
						}
						int* Num1 = (int*)malloc(sizeof(int) * atoi(argv[k + 1]));
						int* Num2 = (int*)malloc(sizeof(int) * atoi(argv[k + 1]));
						for (int l = 0; l < atoi(argv[k + 1]); l++)
						{
							Num1[l] = atoi(argv[k + 2 + 2 * l]);
							Num2[l] = atoi(argv[k + 3 + 2 * l]);
						}
						if (atoi(argv[k + 1]) > nant[0] / 2)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, ">>>warning! Exchange atoms too many times!");
							fclose(fp_olog);
						}
						Exchange_command(nant, typenum, xyz, xyz_e, atoi(argv[k + 1]), Num1, Num2);
						free(Num1);
						free(Num2);
						Num1 = NULL;
						Num2 = NULL;
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
				}
				if (!strcmp("Strain", argv[k]))
				{
					char mode[5];
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						int RandE_Strain_Mode = rand() % 7;
						switch (RandE_Strain_Mode)
						{
						case 0:strcpy(mode, "xx"); break;
						case 1:strcpy(mode, "yy"); break;
						case 2:strcpy(mode, "zz"); break;
						case 3:strcpy(mode, "xy"); break;
						case 4:strcpy(mode, "xz"); break;
						case 5:strcpy(mode, "yz"); break;
						case 6:strcpy(mode, "xyz"); break;
						default:break;
						}
						Strain(mode, vec, vec_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						strcpy(mode, argv[k + 1]);
						if ((!strcmp(mode, "xx")) || (!strcmp(mode, "yy")) || (!strcmp(mode, "zz")) ||
							(!strcmp(mode, "xy")) || (!strcmp(mode, "xz")) || (!strcmp(mode, "yz")))
							MAX_LENGTH += 3;
						else if ((!strcmp(mode, "xyz")))
							MAX_LENGTH += 8;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						if (!strcmp(mode, "xx"))
							Strain_command(mode, vec, vec_e, atof(argv[k + 2]), 0, 0, 0, 0, 0);
						if (!strcmp(mode, "yy"))
							Strain_command(mode, vec, vec_e, 0, atof(argv[k + 2]), 0, 0, 0, 0);
						if (!strcmp(mode, "zz"))
							Strain_command(mode, vec, vec_e, 0, 0, atof(argv[k + 2]), 0, 0, 0);
						if (!strcmp(mode, "xy"))
							Strain_command(mode, vec, vec_e, 0, 0, 0, atof(argv[k + 2]), 0, 0);
						if (!strcmp(mode, "xz"))
							Strain_command(mode, vec, vec_e, 0, 0, 0, 0, atof(argv[k + 2]), 0);
						if (!strcmp(mode, "yz"))
							Strain_command(mode, vec, vec_e, 0, 0, 0, 0, 0, atof(argv[k + 2]));
						if (!strcmp(mode, "xyz"))
							Strain_command(mode, vec, vec_e, atof(argv[k + 2]), atof(argv[k + 3]),
								atof(argv[k + 4]), atof(argv[k + 5]), atof(argv[k + 6]), atof(argv[k + 7]));
					}
					double v = volume(vec);
					double v_e = volume(vec_e);
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec_e[i][j] *= pow(v / v_e, 1 / 3);
						}
					}
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec[i][j] = vec_e[i][j];
						}
					}
				}
				if (!strcmp("Ripple", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Ripple(nant, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += 4;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						Ripple_command(nant, xyz, xyz_e, atoi(argv[k + 1]), atof(argv[k + 2]), atof(argv[k + 3]));
						//The reference value of temp is as follows
						/*case 0: x = 0, y = 1, z = 2; axis = 'z'; break;
						case 1: x = 0, y = 2, z = 1; axis = 'y'; break;
						case 2: x = 1, y = 0, z = 2; axis = 'z'; break;
						case 3: x = 1, y = 2, z = 0; axis = 'x'; break;
						case 4: x = 2, y = 0, z = 1; axis = 'y'; break;
						case 5: x = 2, y = 1, z = 0; axis = 'x'; break;*/
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
				}
				if (!strcmp("Randmove", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Randmove(nant, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					MAX_LENGTH += 3;
					if (argc < MAX_LENGTH)
					{
						FILE* fp_olog = fopen("command_olog", "at+");
						fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
						fclose(fp_olog);
					}
					Randmove_command(nant, xyz, xyz_e, atoi(argv[k + 1]), atof(argv[k + 2]));
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
				}
				if (!strcmp("Slip", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Slip(nant, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += 7;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						Slip_command(nant, xyz, xyz_e, atof(argv[k + 1]), atof(argv[k + 2]),
							atof(argv[k + 3]), atof(argv[k + 4]), atof(argv[k + 5]), atof(argv[k + 6]));
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
				}
				if (!strcmp("Twist", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Twist(nant, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += 3;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						Twist_command(nant, xyz, xyz_e, argv[k + 1][0], atof(argv[k + 2]));
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
				}
				if (!strcmp("Rotation", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Rotation(vec, vec_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += 4;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						Rotation_command(vec, vec_e, atof(argv[k + 1]), atof(argv[k + 2]), atof(argv[k + 3]));
					}
					double v = volume(vec);
					double v_e = volume(vec_e);
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec_e[i][j] *= pow(v / v_e, 1 / 3);
						}
					}
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec[i][j] = vec_e[i][j];
						}
					}
				}
				if (!strcmp("Stripple", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Stripple(nant, vec, vec_e, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += 10;
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						Stripple_command(nant, vec, vec_e, xyz, xyz_e, atof(argv[k + 1]), atof(argv[k + 2]),
							atof(argv[k + 3]), atof(argv[k + 4]), atof(argv[k + 5]), atof(argv[k + 6]),
							atoi(argv[k + 7]), atof(argv[k + 8]), atof(argv[k + 9]));
					}
					double v = volume(vec);
					double v_e = volume(vec_e);
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec_e[i][j] *= pow(v / v_e, 1 / 3);
						}
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec[i][j] = vec_e[i][j];
						}
					}
				}
				if (!strcmp("Permustrain", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Permustrain(nant, typenum, vec, vec_e, xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						MAX_LENGTH += 8 + 2 * atoi(argv[k + 7]);
						if (argc < MAX_LENGTH)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						int* Num1 = (int*)malloc(sizeof(int) * atoi(argv[k + 7]) / 2);
						int* Num2 = (int*)malloc(sizeof(int) * atoi(argv[k + 7]) / 2);
						for (int l = 0; l < atoi(argv[k + 7]); l++)
						{
							Num1[l] = atoi(argv[k + 8 + 2 * l]);
							Num2[l] = atoi(argv[k + 9 + 2 * l]);
						}
						if (atoi(argv[k + 7]) > nant[0] / 2)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, ">>>warning! Exchange atoms too many times!");
							fclose(fp_olog);
						}
						Permustrain_command(nant, typenum, vec, vec_e, xyz, xyz_e, atof(argv[k + 1]), atof(argv[k + 2]),
							atof(argv[k + 3]), atof(argv[k + 4]), atof(argv[k + 5]),
							atof(argv[k + 6]), atoi(argv[k + 7]), Num1, Num2);
						free(Num1);
						free(Num2);
						Num1 = NULL;
						Num2 = NULL;
					}
					double v = volume(vec);
					double v_e = volume(vec_e);
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec_e[i][j] *= pow(v / v_e, 1 / 3);
						}
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec[i][j] = vec_e[i][j];
						}
					}
				}
				if (!strcmp("Rotstrain", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						Rotstrain(xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						if (argc < 6)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s operator is missing parameter!\n", argv[1]);
							fclose(fp_olog);
						}
						Rotstrain_command(vec, vec_e, atof(argv[k + 1]), atof(argv[k + 2]),
							atof(argv[k + 3]), argv[k + 4], atof(argv[k + 5]));
					}
					double v = volume(vec);
					double v_e = volume(vec_e);
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec_e[i][j] *= pow(v / v_e, 1 / 3);
						}
					}
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec[i][j] = vec_e[i][j];
						}
					}
				}
				if (!strcmp("Crossover", argv[k]))
				{
					if (!strcmp("Random", argv[k + 1]) || 'R' == argv[k + 1][0])
					{
						int RandFileNum = rand() % (atoi(argv[3])) + 1;
						char RandFile[MAX_FNAME];
						sprintf(RandFile, "%s%d", argv[2], RandFileNum);
						FILE* FP_RandFile = fopen(RandFile, "r");
						if (FP_RandFile == NULL)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s is open wrong\n", RandFile);
							fclose(fp_olog);
						}
						readposcar(FP_RandFile, _title, _latt, _ifix, _iflg, _nant, _typenum, _elemnum, _elemsym, _vec, _xyz, _fix);
						fclose(FP_RandFile);
						FILE* fp_olog_ = fopen("command_olog", "at+");
						fprintf(fp_olog_, ">>>the other random parent file is %s\n", RandFile);
						fclose(fp_olog_);
						Crossover(latt, _latt, latt_e, nant, typenum, vec, _vec, vec_e, xyz, _xyz, xyz_e);
						system("cat ../buf >> command_olog");
						system("rm ../buf");
					}
					else
					{
						FILE* FP_Parent_ = fopen(argv[k + 1], "r");
						if (FP_Parent_ == NULL)
						{
							FILE* fp_olog = fopen("command_olog", "at+");
							fprintf(fp_olog, "the %s is open wrong\n", argv[k + 1]);
							fclose(fp_olog);
						}
						readposcar(FP_Parent_, _title, _latt, _ifix, _iflg, _nant, _typenum, _elemnum, _elemsym, _vec, _xyz, _fix);
						fclose(FP_Parent_);
						FILE* fp_olog_ = fopen("command_olog", "at+");
						fprintf(fp_olog_, ">>>the other parent file is %s\n", argv[k + 1]);
						fclose(fp_olog_);
						Crossover_command(latt, _latt, latt_e, nant, typenum, vec, _vec, vec_e, xyz, _xyz, xyz_e, argv[k + 2][0], atof(argv[k + 3]));
					}
					double v = volume(vec);
					double v_e = volume(vec_e);
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec_e[i][j] *= pow(v / v_e, 1 / 3);
						}
					}
					for (int i = 0; i < nant[0]; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							xyz[i][j] = xyz_e[i][j];
						}
					}
					for (int i = 0; i < 3; i++)
					{
						for (int j = 0; j < 3; j++)
						{
							vec[i][j] = vec_e[i][j];
						}
					}
				}
			}
			sprintf(SaveFileName[i], "%s%d", argv[4], i + 1);
			FILE* fp_child = fopen(SaveFileName[i], "w");
			if (fp_child == NULL)
			{
				FILE* fp_olog = fopen("command_olog", "at+");
				fprintf(fp_olog, "the %s is open wrong\n", SaveFileName[i]);
				fclose(fp_olog);
			}
			FILE* fp_olog2 = fopen("command_olog", "at+");
			fprintf(fp_olog2, ">>>the child file is %s\n", SaveFileName[i]);
			fprintf(fp_olog2, "*************************************\n");
			fclose(fp_olog2);
			strcpy(title, argv[3]);
			savposcar(fp_child, title, latt, ifix, iflg, nant, typenum, elemnum, elemsym, vec, xyz, fix);
			fclose(fp_child);
			for (int k = 0; k < atoi(argv[3]); k++)//all input files
			{
				if (COMPARE(SaveFileName[i], buf[k], 0.05, 0.25))
				{
					FILE* fp = fopen("command_olog", "r");
					fprintf(fp, "structure is repeated!\n");
					fclose(fp);
					goto again;
				}
			}
			for (int l = 0; l < i; l++)//current output files
			{
				if (COMPARE(SaveFileName[i], SaveFileName[l], 0.05, 0.25))
				{
					FILE* fp = fopen("command_olog", "r");
					fprintf(fp, "structure is repeated!\n");
					fclose(fp);
					goto again;
				}
			}
			delete[] xyz;
			delete[] fix;
			delete[] _xyz;
			delete[] _fix;
			delete[] xyz_e;
			delete[] fix_e;
		}
		delete[] OpenFileName;
		delete[] SaveFileName;
		delete[] buf;
	}


	//EVOPS --mcomp file1 10 savefile2 p1 p2
	int Population::MultiFile_COMPARE(int argc, char* argv[])
	{
		FILE* fp_compare = fopen("compare_olog", "at+");
		if (argc != 7)
			fprintf(fp_compare, "the parament is wrong!");
		char(*File1)[MAX_FNAME] = new char[MAX_NATOM][MAX_FNAME];
		char(*SaveFile)[MAX_FNAME] = new char[MAX_NATOM][MAX_FNAME];
		for (int i = 0; i < atoi(argv[3]); i++)
		{
			int count = 0;
			sprintf(File1[i], "%s%d", argv[2], i + 1);//1-n
			sprintf(SaveFile[i], "%s%d", argv[4], i + 1);
			for (int j = 0; j < i; j++)
			{
				if (COMPARE(File1[i], File1[j], atof(argv[5]), atof(argv[6])))
				{
					fprintf(fp_compare, "the %s and %s are the same!\n", File1[i], File1[j]);
					count += 1;
					break;
				}
			}
			if (!count)
				copy_file(SaveFile[i], File1[i]);
		}
		fclose(fp_compare);
		delete[] File1;
		delete[] SaveFile;
	}
}


