//*************************************************************
//	UNIVERSIDADE FEDERAL DE MINAS GERAIS
//*************************************************************
// Implementacao CPLEX do algorimo de decomposicao de Benders
// para a solucao de Quadratic Assignment Problems
//*************************************************************
// Gilberto de Miranda Junior - 23/05/2002
//*************************************************************
//*************************************************************

# include <stdlib.h>
# include <sys/time.h>
# include <iostream.h>
# include <fstream.h>
# include <stdio.h>
# include <ilcplex/ilocplex.h>
# include "subproajcpxnew2.h"

struct timeval t1, t2;
int main(int argc, char *argv[]){

int h, i, j, k, l, n, n2, **a, *d;
double inf, sup, ub, ubb, *r, *u, *v, **coef, **xl, cv, cfcv, *alpha;
t1.tv_sec = t1.tv_usec = 0;
t2.tv_sec = t2.tv_usec = 0;


ifstream data(argv[1]);
ofstream result("resultADAMSJHONSON", ios::out|ios::app);
int m = atoi(argv[2]);
data>>n;
// alocando matrizes
xl = (double **)malloc((n+1)*sizeof(double*));
for(i = 0 ; i<= n; i++){
xl[i] = (double *)malloc((n+1)*sizeof(double));
}
a = (int **)malloc((n+1)*sizeof(int*));
for(i = 0 ; i<= n; i++){
a[i] = (int*)malloc((n+1)*sizeof(int));
}
/*c = (int **)malloc((n+1)*sizeof(int*));
for(i = 0 ; i<= n; i++){
c[i] = (int*)malloc((n+1)*sizeof(int));
}
coef = (double **)malloc((n+1)*sizeof(double*));
for(i = 0 ; i<= n; i++){
coef[i] = (double *)malloc((n+1)*sizeof(double));
}
b = (int **)malloc((n+1)*sizeof(int*));
for(i = 0 ; i<= n; i++){
b[i] = (int *)malloc((n+1)*sizeof(int));
}*/
n2=n*n;
n4 = n2*n2;
n3 = n2*n;
r = (double *)malloc((n*(n2)+1)*sizeof(double));
v = (double *)malloc((n*(n2)+1)*sizeof(double));
u = (double *)malloc((n*(n2)+1)*sizeof(double));
alpha = (double *)malloc((n2+1)*sizeof(double));
d = (int *)malloc((n4+1)*sizeof(int));

/*for(k = 1; k <= n; k++){
 for(i = 1; i <= n; i++){
        data>>b[k][i];
 }
}
for(k = 1; k <= n; k++){
 for(i = 1; i <= n; i++){
        data>>c[k][i];
 }
}*/

// lendo dados de entrada
for(k = 1; k <= n; k++){
 for(i = 1; i <= n; i++){
       data>>a[k][i];
	a[k][i] = m * a[k][i]/100;
 }
}

for(k=1; k<=n; k++){
 for(i=1; i<=n; i++){
  for(l=1; l<=n; l++){
   for(j=1; j<=n; j++){
    data>>d[(k-1)*n3+(i-1)*n2+(l-1)*n+j];
   }
  }
 }
}

for(k = 1; k <= n; k++){
 for(i = 1; i <= n; i++){
        coef[k][i]=0.0;
 }
}

gettimeofday(&t1,NULL);

ubb = 10000000000000;
sup = 1;
inf = 0;

// construindo ambiente de programacao matematica
IloEnv env;
IloModel mod(env);
IloNumVar t(env, 0, +IloInfinity, ILOFLOAT);
IloNumVarArray x(env, n2 + 1, 0.0, 1.0, ILOINT);

// Introduzindo funcao objetivo
IloExpr expr(env);
expr+= t;

for(k = 1; k <= n; k++){
 for(i = 1; i <= n; i++){
    expr+= a[k][i] * x[((k-1)*n + i)];
 }
}
mod.add(IloMinimize(env,expr));

// Introduzindo restricoes de Assignment

// somatorio em i
for(k=1;k<=n;k++){
 IloExpr expr(env);
 for(i=1;i<=n;i++) expr+= x[((k-1)*n + i)];
 mod.add(1<=expr<=1);
}

// somatorio em k
for(i=1;i<=n;i++){
 IloExpr expr(env);
 for(k=1;k<=n;k++) expr+= x[((k-1)*n + i)];
 mod.add(1<=expr<=1);
}

h = 1;
IloCplex cplex(mod);
while((ubb-inf)>0.0001){
	cplex.solve();
	for(k=1;k<=n;k++){
	 for(i=1;i<=n;i++){
	 	xl[k][i] = (double)cplex.getValue(x[((k-1)*n+i)]);
	 	if(xl[k][i]<=0.001) xl[k][i] = 0;
                if(xl[k][i]>=0.9) xl[k][i] = 1;
	 }
	}

	inf = cplex.getObjValue();
	ub =0;
        cout<<"Solucao primal (x): \n";
	for(k=1;k<=n;k++){
	 for(i=1;i<=n;i++){
	  ub = ub + a[k][i]*xl[k][i];
          cout<<" "<<xl[k][i]<<" ";
	 }
         cout<<"\n";
	}
	cout<<"Custo Linear:  "<<ub<<"\n";
        cout<<" Entrei em solve slave! \n";
	cv = solve_slave(n, b, c, xl, v, u, r, alpha);
        cout<<" Sai de solve slave! \n";
	cout <<"Custo Quadratico:  "<<cv<<"\n";
	cfcv = ub/cv;
	ub = ub + cv;
	if(ub<ubb) ubb = ub;
	sup = ub;

	// acrescentando corte de benders
	for(k=1;k<=n;k++){
	 for(i=1;i<=n;i++){
	  for(j=1;j<=n;j++){
	  if(i!=j)coef[k][i] = coef[k][i] - r[(k-1)*n2+(i-1)*n+j];
	  }
	 }
	}

	for(k=1;k<=n;k++){
	 for(l=(k+1);l<=n;l++){
	  for(i=1;i<=n;i++){
	  coef[k][i] = coef[k][i] - u[(k-1)*n2+(l-1)*n+i];
	  }
	 }
	}
	
	for(k=1;k<=n;k++){
	 for(l=(k+1);l<=n;l++){
	  for(j=1;j<=n;j++){
	  coef[l][j] = coef[l][j] - v[(k-1)*n2+(l-1)*n+j];
	  }
	 }
	}
	
	IloExpr expr(env);
	expr+= t;
	for(k = 1; k <= n; k++){
	 for(i = 1; i <= n; i++){
	    expr+= -1.0*coef[k][i] * x[(k-1)*n + i];
	 }
	}
	
	double rhs = 0.0;
	for(k = 1; k <= n; k++){
	 for(l = 1; l <= n; l++){
	    rhs+= alpha[(k-1)*n + l];
	 }
	}
	
	mod.add(rhs<=expr<=+IloInfinity);

 	for(k=1;k<=n;k++){
	 for(i=1;i<=n;i++){
	 coef[k][i] = 0.0;
	 }
	}

	cout<<" Contador h : "<<h<<endl;
	cout<<" Infimo     : "<<inf<<endl;
	cout<<" upperbound : "<<ubb<<endl;
	cout<<" supremo    : "<<sup<<endl;
	cout<<" cf/cv      : "<<cfcv<<endl;
	cout<<" gap(%)     : "<<(ubb-inf)/inf<<endl;
	h++;

}//end-while
for(k=1;k<=n;k++){
 for(i=1;i<=n;i++){
 cout<<xl[k][i]<<" ";
 }
 cout<<"\n";
}


gettimeofday(&t2, NULL);

long  extimes = (t2.tv_sec - t1.tv_sec);

	cout<<" contador h : "<<h<<endl;
	cout<<" Infimo     : "<<inf<<endl;
	cout<<" upperbound : "<<ubb<<endl;
	cout<<" supremo    : "<<sup<<endl;
	cout<<" cf/cv      : "<<cfcv<<endl;
	cout<<" gap(%)     : "<<(ubb-inf)/inf<<endl;
	cout<<" tempo de execucao [s] : "<< extimes<<"\n";

// impressao em arquivo
        result<<" Entrada: "<<argv[1];
	result<<" n : "<<n;
	result<<" Contador h : "<<h;
	result<<" Infimo     : "<<inf;
	result<<" Upperbound : "<<ubb;
	result<<" Supremo    : "<<sup;
	result<<" cf/cv      : "<<cfcv;
	result<<" Gap(%)     : "<<(ubb-inf)/inf;
	result<<" Tempo de execucao: [s] : "<< extimes<<"\n";

return 0;
}// end-of-file

