//   *****************************************************************************************************************************************************
//   *																																					 *
//   *   Vehicle Routing Programming																													 *
//   *   Compiled by Microsoft Visual Studio 2008																										 *
//   *   Written by Yang Zuo in 2009																													 *
//   *   Copyright by Uncertainty Theory Laboratory																										 *
//   *																																					 *
//   *****************************************************************************************************************************************************

#include<iostream>
#include<vector>
#include<iomanip>
#include<cmath>
#include <fstream>
#include"UTLab.h"

#define pi 3.1415
#define popsize 3 //the number of Chromosome
#define pc 0.4
#define pm 0.3
#define N 1000
#define node  7   //the number of Nodess
#define vehicle 3 //the number of Vehicle
using namespace std;//using the standard namespace std

//   *****************************************************************************************************************************************************
//   *																																					 *
//   *		The introduction of	variables.																												 *	
//   *																																					 *
//   *****************************************************************************************************************************************************
double cr0;                    //the uncertainty measure  alpha.
int x[popsize][node+1];        //the decision variable X
int y[popsize][vehicle+1];     //the decision variable Y
double t[popsize][vehicle+1];  //the decision variable T

int x_tmp[popsize][node+1];         //temp storaging space for the decision variable  X
int y_tmp[popsize][vehicle+1];      //temp storaging space for the decision variable  Y
double t_tmp[popsize][vehicle+1];   //temp storaging space for the decision variable  T

vector<double> dis(99,0);
vector<  vector<double>  >  nod ( node+1, dis);
vector<  vector<  vector<double>   >   >   Ttime( popsize, nod);  
//  Ttime[i][j] is a 99-dimension vector storaging  the uncertainty distribution of the time when the vehicle  arrives at the  node X[i][j]


int d[node+1][node+1]=     
{ 0, 1, 2, 3, 4, 5, 6, 7, 
  1, 0, 1, 2, 3, 4, 5, 6, 
  2, 1, 0, 1, 2, 3, 4, 5, 
  3, 2, 1, 0, 1, 2, 3, 4, 
  4, 3, 2, 1, 0, 1, 2, 3, 
  5, 4, 3, 2, 1, 0, 1, 2, 
  6, 5, 4, 3, 2, 1, 0, 1, 
  7, 6, 5, 4, 3, 2, 1, 0, 
};        
//The distances between nodes.

int tw[node+1][2]=
{0,    0,
 7,    9,
 7,    9,
 14,   17,
 14,   17,
 14,   17,
 19,   21,
 19,   21,
};
//The time windows of nodes.  

int tt[node+1][node+1][2]=   
{
	0,0,    2,1,    4,1,   6,1,   8,1,   10,1,  12,1,  14,1,    
    2,1,  	0,0,    2,1,   4,1,   6,1,   8,1,   10,1,  12,1,   
	4,1,    2,1,  	0,0,   2,1,   4,1,   6,1,   8,1,   10,1,    
	6,1,   	4,1,    2,1,   0,0,   2,1,   4,1,   6,1,   8,1,   
	8,1,   	6,1,   	4,1,   2,1,   0,0,   2,1,   4,1,   6,1,  
    10,1,   8,1,   	6,1,   4,1,   2,1,   0,0,   2,1,   4,1, 
	12,1,   10,1,   8,1,   6,1,   4,1,   2,1,  	0,0,   2,1, 
	14,1,   12,1,   10,1,  8,1,   6,1,   4,1,   2,1,   0,0, 
};   
// The travelling times between Node i and Node j  is an normal uncertian variable N(2|i-j|,1).



//   *****************************************************************************************************************************************************
//   *																																					 *
//   *		The introduction of	functions																												 *	
//   *																																					 *
//   *****************************************************************************************************************************************************

               
void norm_distr(vector<double> &distr, double a,  double b); //Generating the uncertainty distribution of the normal uncertain variable N(a,b) vai 99 method
void TtimeCreate(int m); //Generating the uncertainty distribution of the times of Chromosome m
double cr(int m);       //Computing the uncertaity measure alpha for Chromosome m.
  
//The following four functions are some operations of vector.
void vec_plus1( vector<double> &X, double Y);             //Plusing  Vector X and Number Y.
void vec_plus2( vector<double> &X, vector<double> &Y);    //Plusing Vector X and Vevtor Y.
void vec_max(vector<double> &X, double Y  );              //Taking large between Vector X and Number Y.

void GAfordis(double cr0);                // Genetic algorithms
void mutation(int a, double cr0);		  //Mutation operation
void crossover(int a, int b, double cr0);  //Crossover Operation
bool feasible(int m,double cr0);          //Testing the feasibility of Chromosome m.

double cr(int m);         //Computing the smallest alpha for the m-th chromosome.
void bmap(int a, int b);  // temp a--> original b duplicate chromosome.
void map(int a, int b);   // original a-->temp b duplicate chromomsome.
int distance(int m);      //Computing the total travel distance of the m-th chromosome 
double min(double a, double b);
double max(double a, double b);

int rdint(int a,int b); //Generating a random number between a and b
void init(double cr0);  //Iniliating the decison varialbes X, Y, T.
void chrom(int m);      //Outputing the chromosome m.
void output();          //Outputing the result.
void rresult();         //Outputing the resulut into result.txt


//   *****************************************************************************************************************************************************
//   *																																					 *
//   *		The function  main.         																												 *	
//   *																																					 *
//   *****************************************************************************************************************************************************


void main()
{  
	
	
    cout<<"Please input the uncertain measure alpha :   ";
	cin>>cr0;
	cout<<"Press Enter to start. "<<endl;
	getchar();
	  
	init(cr0);    	
	GAfordis(cr0);


    ofstream result ("result.txt");
	rresult();
    output();

	
	getchar();


}


//   *****************************************************************************************************************************************************
//   *																																					 *
//   *		The definitions of  functions.         																										 *	
//   *																																					 *
//   *****************************************************************************************************************************************************

void rresult(){
	int i, k;
	double temp;
	ofstream result ("result.txt");
	if(result.is_open()){

		result<<"The Result is :"<<endl;
		result<<"distance="<<distance(0)<<endl;

		result<<"x=";
	for(i=1;  i<=node; i++){
			result<<x[0][i]<<"  ";
		}
		result<<endl;
		result<<"y=";
		for(i=1; i<=vehicle; i++){
			result<< y[0][i]<<"  ";
		}
		result<<endl;

		TtimeCreate(0);
		for(i=1; i<=vehicle; i++){
			
			temp=tw[ x[0][y[0][i-1]+1] ][1] -  (  Ttime[0][y[0][i-1]+1][int(cr0/0.01 -1)]  -  t[0][i] );
			t[0][i]=temp;
		}
		result<<"t=";
		for(i=1; i<=vehicle; i++){			
			result<<int(t[0][i])<<":"<< int ((t[0][i]-int(t[0][i]))*60);
			if(t[0][i]-int(t[0][i]))
				result<<"  ";
			else
				result<<"0"<<"  ";
		}
		result<<endl;


		for( i=1; i<=vehicle; i++){
			if(y[0][i]!=y[0][i-1]){
				
				result<< "0  -  ";
				for( k=y[0][i-1]+1; k<=y[0][i]; k++){
					result<<x[0][k]<<"  -  ";
				}
				result<<"0"<<endl;
			}
		}
	
	}


	result.close();
}



void output(){
	int i, k;
	double temp;
		cout<<"The Result is :"<<endl;

		cout<<"distance="<<distance(0)<<endl;

		cout<<"x=";
		for(i=1;  i<=node; i++){
			cout<<x[0][i]<<"  ";
		}
		cout<<endl;
		cout<<"y=";
		for(i=1; i<=vehicle; i++){
			cout<< y[0][i]<<"  ";
		}
		cout<<endl;

		TtimeCreate(0);
		for(i=1; i<=vehicle; i++){
			
			temp=tw[ x[0][y[0][i-1]+1] ][1] -  (  Ttime[0][y[0][i-1]+1][int(cr0/0.01 -1)]  -  t[0][i] );
			t[0][i]=temp;
		}
		cout<<"t=";
		for(i=1; i<=vehicle; i++){	
		
			cout<<int(t[0][i])<<":"<<int(  (t[0][i]-int(t[0][i]))*60   );
			if(t[0][i]-int(t[0][i])>0)
				cout<<"  ";
			else
				cout<<"0"<<"  ";
		}
		cout<<endl;


		for( i=1; i<=vehicle; i++){
			if(y[0][i]!=y[0][i-1]){
				
				cout<< "0  -  ";
				for( k=y[0][i-1]+1; k<=y[0][i]; k++){
					cout<<x[0][k]<<"  -  ";
				}
				cout<<"0"<<endl;
			}
		}
}

void vec_plus1( vector<double> &X, double Y){
	vector<double>::iterator   iter=X.begin();   
	
	for(  iter=X.begin(); iter!=X.end(); iter++)
		*iter=*iter+Y;	
}

void vec_plus2( vector<double> &X, vector<double> &Y ){
    
	if(  (X.size())==(Y.size()) ){
		for(int i=0; i<X.size(); i++)
			X[i]=X[i]+Y[i];
	}
	else
		cout<<"Can't plus the two vector."<<endl;
}


void vec_max(	vector<double> &X, double Y  ){
		for(int i=0; i < X.size(); i++)
			X[i]=  (X[i]> Y)? X[i] : Y;
}
		


int rdint(int a,int b)
{
	double x;
	int temp;
	x=myu(a,b+0.99);
	temp=int(floor(x));
	return(temp);
}


void  norm_distr(vector<double> &distr, double a,  double b){  

	distr.clear();
	for(int i=0;i<99;i++){
		distr.push_back( a  -log(100/double(i+1) - 1)*b/pi ) ;
		
	}

}

void chrom(int a)
{
	int distance(int a);
	double cr(int a);
	int j;
	cout<<endl;
	cout<<"x= ";
	for(j=1;j<=node;j++)
		cout<<x[a][j]<<"  ";		
	   
	cout<<endl;
	cout<<"y= ";
	for(j=1;j<=vehicle;j++)
		cout<<y[a][j]<<"  ";		
	    
	cout<<endl;
	cout<<"t= ";
	for(j=1;j<=vehicle;j++)
		cout<<int(t[a][j])<<":"<<int(  (t[a][j]-int(t[a][j]))*60   )<<"  ";
	cout<<endl; 
	cout<<"  dist  "<<distance(a)<<endl;
}


double min(double a, double b)
{
	if (a>b)
		return b;
	else return a;
}
double max(double a, double b)
{
	if (a<b)
		return b;
	else return a;
}

bool feasible(int m,double cr0)
{
	double cr(int );
	if(cr(m)>=cr0)
		return true;
	else 
		return false;
}



void TtimeCreate(int m){
 
	vector<int> Xveh (node+1,0);
	vector<double> distr(99);
	vector<double> temp( 99 );

   
	int i,j,k=0;
	for( int i=1;i<=vehicle;i++)
	{
		
		if(y[m][i]!=y[m][i-1])
		{
			for( j=k+1;j<=y[m][i];j++)
				Xveh[j]=i;
			k=j-1;
		}
		else ;
	}


	for( i=1; i<=node; i++)
	{		
		if(Xveh[i]!=Xveh[i-1])
		{
			 norm_distr(distr, tt[0][ x[m][i] ][0],  tt[0][ x[m][i] ][1] );	
			 vec_plus1(distr, t[m][ Xveh[i] ]);			
			 Ttime[m][i]=distr;			
		}


		else
		{
			norm_distr(distr, tt[ x[m][i]][ x[m][i-1] ][0],  tt[ x[m][i]][ x[m][i-1] ][1]  );
		
			temp =Ttime[m][i-1];
			vec_max(temp, tw[x[m][i-1] ][0]);
	
			vec_plus2(temp, distr);
			Ttime[m][i]=temp;			
		}
	
	}

	Xveh.clear();
	distr.clear();
	temp.clear();
}

		



void init(double cr0)
{
	int i,j,k,w,nm;
	double counter=0.;
	double cr_tmp;
	int ct, num=1;	
	double cr(int m);
	int distance(int m);
	double cro;
	int tt[node+1];
	int te;
	int judge;
	for(i=1;i<=node;i++)
		tt[i]=0;
		
	

	for(i=0;i<popsize;i++)	
	{
		nm=0;
		cro=0.;
		cr_tmp=0.;
		while(cro<(cr0))
		{

            //initializing X
			for(k=1;k<=node;k++)
				tt[k]=0;
			j=1;
			while(j<=node)
			{
				te=rdint(1,node);
				judge=1;
				for(k=1;k<=node;k++)
					if(tt[k]==te)
						judge=judge*0;
				if(judge==1)
				{
					x[i][j]=te;
					tt[j]=te;
					j++;
				}
				else
					continue;
			}

			
	      //initializing Y	
			for(k=1;k<vehicle;k++)
			{
				y[i][k]=rdint(0,node);
			}

			//sort
			for(j=1; j<vehicle; j++){
				for(k=j+1; k<vehicle; k++) {
					if(y[i][j] > y[i][k]) 
					{
						w = y[i][k];
						y[i][k] = y[i][j];
						y[i][j] = w;
					}
				}
			}
			y[i][0]=0;
			y[i][vehicle]=node;  


			//initializing T	
			for(k=1;k<=vehicle;k++){
				t[i][k]=double(rdint(0,24));
			}
			for(k=2; k<=vehicle; k++){
				if(y[i][k-1]==node){
					y[i][k]=0;
				}
			}
			cro=cr(i);
	
			
		
		}
		system("cls");
		counter=double(i)/popsize*100;
		ct=floor(counter);
		cout<<"Initializing..."<<ct<<"%"<<endl;
	}

	system("cls");
	cout<<"Initializing...100%"<<endl;

}


void crossover(int a, int b, double cr0)
{
	int i,tm;
	double c;
	
	double tp,cro1,cro2;
	int yp;
	
	bool feasible(int a,double cr0);
	void map(int a, int b);
	void bmap(int a, int b);
	double max(double a, double b);
	double cr(int m);

	tm=0;
	cro1=0.;
	cro2=0.;

	map(a,a);
	map(b,b);

	while((tm<5)&&(max(cro1,cro2)<cr0))
	{
		c=myu(0,1);
		for(i=1;i<=vehicle;i++)
		{
			tp =c*t[a][i]+(1-c)*t[b][i];
			t[b][i]=t[a][i]*(1-c)+t[b][i]*c;
			t[a][i]=tp;
		}//crossover  t[a], t[b]
		for(i=1;i<=vehicle;i++)
		{
			yp=y[a][i];
			y[a][i]=y[b][i];
			y[b][i]=yp;
		}//crossover y.

		cro1=cr(a);
		cro2=cr(b);

		if(cro1<cr0)
			bmap(a,a);
		
		if(cro2<cr0)
			bmap(b,b);
	
		tm++;
		
	}
		
	map(a,a);
	map(b,b);
	

}

void mutation(int a, double cr0)
{
	int i,j,k, n1, n2, nn, xx, w, tm;
	double M;

	void map(int a, int b);
	void bmap(int a, int b);

	
	double cr(int m);
	double dt[vehicle+1];
	tm=0;
	M=1.;
	
	map(a,a);

mutatet: n1=rdint(1,node);
	n2=rdint(1,node);
	if(n1>n2)
	{
		nn=n1;
		n1=n2;
		n2=nn;
	}
	
	for(i=n1;i<=n2;i++)
	{
		nn=rdint(i,n2);
		xx=x[a][i];
		x[a][i]=x[a][nn];
		x[a][nn]=xx;
	}
	//mutate x

	n1=rdint(1,vehicle-1);
	n2=rdint(1,vehicle-1);

	if(n1>n2)
	{
		nn=n1;
		n1=n2;
		n2=nn;
	}
	
	for(i=n1;i<=n2;i++)
	{
		y[a][i]=rdint(0,node);
	}

	for(j=1; j<vehicle; j++)
		for(k=j+1; k<vehicle; k++) 
			if(y[a][j] > y[a][k]) 
			{
				w = y[a][k];
				y[a][k] = y[a][j];
				y[a][j] = w;
			}	

	//mutate y
	 
	dt[0]=0.;
	for(i=1;i<=vehicle;i++)
	{
		dt[i]=myu(0.,60.) - 30.;
	}
	for(i=1;i<=vehicle;i++)
	{
		if(dt[i]>=0.)
			t[a][i]=min(460., t[a][i] + dt[i]);
		else
			t[a][i]=max(10., t[a][i] + dt[i]);
	}
	
	
	if(cr(a)>=cr0)
		;
	else if(tm<6)
	{	
		bmap(a,a);
		
		M=myu(0.,M);
		tm++;
		goto mutatet;
	}
	else bmap(a,a);
	
	//mutate t;

}


int distance(int m)//the m-th chrom's distance
{
	int dis=0;
	int i,j;



	for(i=0;i<node;i++)
		dis+=d[x[m][i]][x[m][i+1]];

	dis+=d[x[m][node]][0];

	
	for(j=1;j<vehicle;j++)
	{
		if(y[m][j]!=y[m][j-1])
		{
			dis+=d[x[m][y[m][j]]][0];
			dis+=d[x[m][y[m][j]+1]][0];
			dis-=d[x[m][y[m][j]]][x[m][y[m][j]+1]];
		}

	}	
	return dis;
	
}//distance tested.


void map(int a, int b)// original a-->temp b duplicate chrom.
{
	int i;
	for(i=0;i<=node;i++) 
		x_tmp[b][i]=x[a][i];
	for(i=0;i<=vehicle;i++)
	{
		y_tmp[b][i]=y[a][i];
		t_tmp[b][i]=t[a][i];
	}
}

void bmap(int a, int b)// temp a--> original b duplicate chrom.
{
	int i;
	for(i=0;i<=node;i++) 
		x[b][i]=x_tmp[a][i];
	for(i=0;i<=vehicle;i++)
	{
		y[b][i]=y_tmp[a][i];
		t[b][i]=t_tmp[a][i];
	}

}


double cr(int m)
{
	int i, key;
	double cre=1.0;
	vector<double> a(node+1,0);
	
	key=int(cr0/0.01);

	TtimeCreate(m);

	for( i=1; i<=node; i++){
		if (tw[x[m][i]][1]>=Ttime[m][i][98])
		{
			a[i]=0.99;
		}
		else{
			if(tw[x[m][i]][1]<Ttime[m][i][0])
			{
				a[i]=0;
			}
			else{
				for(int j=98; j>0; j--){
					if(Ttime[m][i][j-1]<=tw[x[m][i]][1]&&tw[x[m][i]][1]<Ttime[m][i][j])
					{
						a[i]=(double(j)/100);
					}
					else;
				}
			}
		}
	}

    for( i=1; i<=node; i++)
		cre= min(cre, a[i]);

	return (cre);
}



void GAfordis(double cr0)
{
	int r;
	double r1;
	int c=0;
	int i,j,k,times;
	int cross[2];
	double ob[popsize],obp[popsize];
	double fit[popsize];
	double sumob;
	int best;
	int bestob;

	best=0;
	
	
	
	int xb[node+1];
	int yb[vehicle+1];
	double tb[vehicle+1];  //save the best

	for(i=0;i<=node;i++)
		xb[i]=x[0][i];
	for(i=0;i<=vehicle;i++)
	{
		yb[i]=y[0][i];
		tb[i]=t[0][i];
	}
	bestob=distance(0);

	for(i=0;i<popsize;i++)
		fit[i]=0.;
	for(times=0;times<N;times++)
	{


		if(((times%1500)==0)&&(times>1))
			getchar();
	
		
		//crossover
		for(i=0;i<popsize;i++)
		{
			r=myu(0.,1.);
			//pk=1;
			if(r<pc)
			{
				cross[c]=i;
				c++;
			}
			if(c==2)
			{
				crossover(cross[0],cross[1],cr0);
				c=0;
			}
		}

		//mutation
		for(i=0;i<popsize;i++)
		{
			r=myu(0.,1.);
			
			if(r<pm)
			{
				map(i,i);
				mutation(i,cr0);
				if(cr(i)<cr0)
				{
					bmap(i,i);
				
				}			
			}
		}

		//objective value
		for(i=0;i<popsize;i++)
		{
			ob[i]=double(distance(i));
			obp[i]=100./(ob[i]-220.);
		}

		
				
		sumob=0.;
		for(i=0;i<popsize;i++)
		{
			sumob=sumob+obp[i];
		}

		//fitness
		fit[0]=obp[0]/sumob;
		for(i=1;i<popsize;i++)
		{
			fit[i]=fit[i-1]+(obp[i]/sumob);
		}
		
		//roulette wheel
		for(i=0;i<popsize;i++)
		{
			r1=myu(0., fit[popsize-1]);
			if(r1<=fit[0])
			{
				map(0,i);
			}
			else 
			{
				for(k=0;k<popsize-1;k++)
				{
					if((r1>fit[k])&&(r1<=fit[k+1]))
					{
						map(k+1,i);   //duplicate to temp
					}	
				}
			}
		}
		
		for(i=0;i<popsize;i++)
		{
			bmap(i,i);  //duplicate back
		}
				
		k=0;r=1000;		
		
		for(i=0;i<popsize;i++)
		{
			if(r>distance(i))
			{
				r=distance(i);
				k=i;
			}			
		}

		if(r<bestob)
		{
			best=k;
			for(j=0;j<=node;j++)
				xb[j]=x[k][j];
			for(j=0;j<=vehicle;j++)
			{
				yb[j]=y[k][j];
				tb[j]=t[k][j];
			}
			bestob=r;
			for(j=0;j<=node;j++)
			x[0][j]=xb[j];
			for(j=0;j<=vehicle;j++)
			{
				y[0][j]=yb[j];
				t[0][j]=tb[j];
			}
		}
		else
		{
			for(j=0;j<=node;j++)
				x[0][j]=xb[j];
			for(j=0;j<=vehicle;j++)
				{
					y[0][j]=yb[j];
					t[0][j]=tb[j];
				}
		}
	}
}