// Project Investment Optimization (CCP)
// Written by Microsoft Visual C++
// Copyright by UTLab @ Tsinghua University
// http://orsc.edu.cn/UTLab

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <conio.h>
#include "UTLab.h"
#include <time.h>

#define N 7    // number of decision variables
#define M 2     // number of objectives
#define TYPE -1 // 1=max; -1=min
#define GEN 3000
#define POP_SIZE 60
#define rate 0.006   // interest rate per month
#define P_MUTATION 0.1
#define P_CROSSOVER 0.2

#define Cycles 5000     // number of simulation cycles


double OBJECTIVE[POP_SIZE+1][M+1], q[POP_SIZE+1];
int CHROMOSOME[POP_SIZE+1][N+1];
static void   initialization(void);
static void   f(int x[N+1], double y[M+1]);
static void   evaluation(int gen);
static void   selection(void);
static void   crossover(void);
static void   mutation(void);
static void   Objectives(void);
static int    constraint_check(double x);
static void   rearrange(int x[]);

static void rearrange(int x[]) // to rearrange the schedule x!!
{
double temp[4];

//Path 1
temp[1]=x[1];
temp[2]=x[2];
temp[3]=x[5];
sortminmax(temp,1,3);
x[1]=temp[1];
x[2]=temp[2];
x[5]=temp[3];

//Path 2
temp[1]=x[1];
temp[2]=x[3];
temp[3]=x[5];
sortminmax(temp,1,3);
x[1]=temp[1];
x[3]=temp[2];
x[5]=temp[3];

//Path 3
temp[1]=x[1];
temp[2]=x[3];
temp[3]=x[6];
sortminmax(temp,1,3);
x[1]=temp[1];
x[3]=temp[2];
x[6]=temp[3];

//Path 4
temp[1]=x[1];
temp[2]=x[3];
temp[3]=x[7];
sortminmax(temp,1,3);
x[1]=temp[1];
x[3]=temp[2];
x[7]=temp[3];

//Path 5
temp[1]=x[1];
temp[2]=x[4];
temp[3]=x[7];
sortminmax(temp,1,3);
x[1]=temp[1];
x[4]=temp[2];
x[7]=temp[3];
}


static void f(int x[N+1], double y[M+1])
{
   int i,temp1;
   double sum[Cycles+1],P[6],xi[12];
   static double c[12]= {0,1500,1800,430,1600,800,500,2000,2100,550,530,630};
    
     temp1=0;
   for(i=1; i<=Cycles; i++) {
                        sum[i]=0;
        P[1]=x[1];
		P[2]=x[1];
		P[3]=x[1];
		P[4]=x[1];
		P[5]=x[1];

	        xi[1]=myn(9,3);
			
			xi[2]=myn(5,2);
						
			xi[3]=myn(10,3);
            
			xi[4]=myn(6,2);
			
			xi[5]=myn(8,2);
						
			xi[6]=myn(8,2);
			
			xi[7]=myn(9,3);
			
			xi[8]=myn(6,1);
			
			xi[9]=myn(10,2);
			
			xi[10]=myn(15,3);
			
			xi[11]=myn(11,2);
			
            P[1]+=xi[1]; 
			if(P[1]<x[2]) P[1]=x[2];
            P[1]+=xi[4]; 
			if(P[1]<x[5]) P[1]=x[5];
            
            P[2]+=xi[2]; 
			if(P[2]<x[3]) P[2]=x[3];
                        P[3]=P[2];
                        P[4]=P[2];
            P[2]+=xi[5]; 
			if(P[1]<P[2]) P[1]=P[2];

             P[1]+=xi[9]; 
			
			// P 1 represents P 1,2
			
	    			
            P[3]+=xi[6]; 
			if(P[3]<x[6]) P[3]=x[6];
            P[3]+=xi[10]; 

			// P 3 represents P 3

            P[4]+=xi[7]; 
			if(P[4]<x[7]) P[4]=x[7];
                        
                        
            P[5]+=xi[3];
			if(P[5]<x[4]) P[5]=x[4];
			P[5]+=xi[8];
			
            if(P[4]<P[5]) P[4]=P[5];
            P[4]+=xi[11];
			//P 4 represents P 4,5			
			
                     

                        if(P[1]<P[3]) P[1]=P[3];
						if(P[1]<P[4]) P[1]=P[4];
                        

                        if(P[1]-int(P[1])==0)            
			P[0]=P[1];
						else P[0]=int(P[1])+1;

			if(P[0]-x[1] <= 32) temp1=temp1+1;
				

            sum[i]+=c[1]*pow(1+rate,ceil(P[0]-x[1]));
			sum[i]+=c[2]*pow(1+rate,ceil(P[0]-x[1]));
			sum[i]+=c[3]*pow(1+rate,ceil(P[0]-x[1]));
			sum[i]+=c[4]*pow(1+rate,ceil(P[0]-x[2]));
			sum[i]+=c[5]*pow(1+rate,ceil(P[0]-x[3]));
			sum[i]+=c[6]*pow(1+rate,ceil(P[0]-x[3]));
			sum[i]+=c[7]*pow(1+rate,ceil(P[0]-x[3]));
			sum[i]+=c[8]*pow(1+rate,ceil(P[0]-x[4]));
			sum[i]+=c[9]*pow(1+rate,ceil(P[0]-x[5]));
			sum[i]+=c[10]*pow(1+rate,ceil(P[0]-x[6]));
			sum[i]+=c[11]*pow(1+rate,ceil(P[0]-x[7]));
               }
   
   y[1] = findmaxn(sum,1,Cycles,500);
   y[2] = temp1/5000.;  
  //printf("%f,  ",y[2]);
}


static void Objectives(void)
{
	double y[M+1];
	int i,j,x[N+1];

	for(i=1; i<=POP_SIZE; i++) {
		for(j=1;j<=N;j++) { x[j] = CHROMOSOME[i][j];}
		f(x,y);
		OBJECTIVE[i][1] = y[1];
	    OBJECTIVE[i][2] = y[2];
	}
	for(i=1; i<=POP_SIZE; i++)
	  OBJECTIVE[i][0]= OBJECTIVE[i][1];
}


static int constraint_check(int x[])
{
	double y[M+1];
	f(x,y);
	if(x[1]<0) return 0;
	if(y[2]<0.9) return 0;
	return 1;
}


static void initialization(void) 
{
  int x[N+1];
  int i,j;

  for(i=1; i<=POP_SIZE; i++) {
	  mark:
	  x[7]=int(myu(0,20));
	  x[6]=int(myu(0,20));
	  x[5]=int(myu(0,20));
	  x[4]=int(myu(0,20));
	  x[3]=int(myu(0,20));
	  x[2]=int(myu(0,20));	  
	  x[1]=int(myu(0,20));

          rearrange(x); // rearrange the schedule x
	  
	  if(constraint_check(x) == 0) goto mark;
	  for(j=1; j<=N; j++) CHROMOSOME[i][j] = x[j]; 
	  printf("N=%d",i);
  }
}

int main( void )
{
   int i, j;
   double a;
   FILE *fp;

   srand((unsigned)time(NULL));
   
   q[0]=0.05; a = 0.05;
   for(i=1;i<=POP_SIZE;i++) {a=a*0.95; q[i]=q[i-1]+a;}
   
   fp=fopen("RESULT.dat","w");
   
   initialization();
   
   evaluation(0);
   for(i=1; i<=GEN; i++) {
	  selection();
	  crossover();
	  mutation();
	  evaluation(i);

      fprintf(fp,"No. %d,",i);
	  for(j=1; j<=N; j++) fprintf(fp,"%d,",CHROMOSOME[0][j]);
	  for(j=1; j<=M; j++) fprintf(fp," %8.15f,  ",OBJECTIVE[0][j]);
      fprintf(fp,"%8.15f,\n",OBJECTIVE[0][0]);

	  printf("\nGeneration NO.%d\n", i);
	  printf("x=(");
	  for(j=1; j<=N; j++) {
		  if(j<N) printf("%d,",CHROMOSOME[0][j]);
		  else printf("%d",CHROMOSOME[0][j]);
	  }
	  	  if(M==1) printf(")\n f=%8.15f\n", OBJECTIVE[0][1]);
	  else {
	      printf(")\nf=(");
	      for(j=1; j<=M; j++) {
		     if(j<M) printf("%8.15f,", OBJECTIVE[0][j]);
		     else printf("%8.15f", OBJECTIVE[0][j]);
		  }
          printf(")  Aggregating Value=%8.15f\n",OBJECTIVE[0][0]);
	  }
   }
   printf("\nThe final result has been written in RESULT.dat.\n\n");
   fclose(fp);
   return 1;
}


static void evaluation(int gen)
{
  double a;
  int   b, i, j, k, label;

  Objectives();

  if(gen == 0){
	 for(k=0; k<=M; k++) OBJECTIVE[0][k] = OBJECTIVE[1][k];
	 for(j = 1; j <= N; j++) CHROMOSOME[0][j] = CHROMOSOME[1][j];
  }

  for(i=0; i<POP_SIZE; i++){
	  label = 0;  a = OBJECTIVE[i][0];
	  for(j=i+1; j<=POP_SIZE; j++)
		 if((TYPE*a)<(TYPE*OBJECTIVE[j][0])) {
			 a = OBJECTIVE[j][0];
			 label = j;
		 }
	  if(label != 0) {
		 for(k=0; k<=M;k++) {
			 a = OBJECTIVE[i][k];
			 OBJECTIVE[i][k] = OBJECTIVE[label][k];
			 OBJECTIVE[label][k] = a;
		 }
		 for(j = 1; j <= N; j++) {
			 b = CHROMOSOME[i][j];
			 CHROMOSOME[i][j] = CHROMOSOME[label][j];
			 CHROMOSOME[label][j] = b;
		 }
	  }
  }
}

static void selection()
{
  double rr;
  int   i, j, k, temp[POP_SIZE+1][N+1];

  for(i=1; i<=POP_SIZE; i++) {
	  rr = myu(0, q[POP_SIZE]);
	  for(j=0; j<=POP_SIZE; j++) {
		  if(rr <= q[j]) { 
			  for(k=1; k<=N; k++) temp[i][k] = CHROMOSOME[j][k];
			  break;
		  }
	  }
  }

  for(i=1; i<=POP_SIZE; i++)
	 for(k=1; k<=N; k++)
		 CHROMOSOME[i][k] = temp[i][k];
}

static void crossover()
{
  int   i, j, jj, k, x[N+1], y[N+1];
  double rr;

  for(i = 1; i <= POP_SIZE/2; i++) {
	 if(myu(0,1)>P_CROSSOVER) continue;
	 j = (int)myu(1,POP_SIZE);
	 jj = (int)myu(1,POP_SIZE);
	 rr = myu(0,1);
	 for(k = 1; k <= N; k++) {
		 x[k] = int(rr*CHROMOSOME[j][k]+(1-rr)*CHROMOSOME[jj][k]);
		 y[k] = int(rr*CHROMOSOME[jj][k]+(1-rr)*CHROMOSOME[j][k]);
	 }
	 rearrange(x);
	 if(constraint_check(x) == 1) 
		 for(k=1; k<=N; k++) CHROMOSOME[j][k] = x[k];
		 rearrange(y);
     if(constraint_check(y) == 1)
		 for(k=1; k<=N; k++) CHROMOSOME[jj][k] = y[k];
  }
}



static void mutation(void)
{
  int i, j, y[N+1];
  double infty, direction[N+1];
  double precision = 0.0001;
  double INFTY = 30;

  for(i=1; i<=POP_SIZE; i++) {
	  if(myu(0,1)>P_MUTATION) continue;
	  for(j=1; j<=N; j++) direction[j] = myu(-1,1);
	  infty = myu(0,INFTY);
	  while(infty>precision) {
		  for(j=1; j<=N; j++) y[j] = int(CHROMOSOME[i][j]+infty*direction[j]);
		  rearrange(y);
		  if(constraint_check(y) == 1) { 
			 for(j=1; j<=N; j++) CHROMOSOME[i][j] = y[j];
			 break;
		  }
		  infty = myu(0,infty);
	  }
  }
}