#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <armadillo>
#pragma warning(disable:4996)
using namespace std;
using namespace arma;

int rand_int(int n) {
	return randi<int>(distr_param(0, n - 1));
}
uvec random_sample(uvec x, int k) {
	int n = x.n_elem;
	uvec sample(k);
	if (k == 1)
		sample(0) = x(rand_int(n));
	else {
		uvec pom = shuffle(x);
		sample = pom(span(0, k - 1));
	}
	return sample;
}
uvec u_diff_v(uvec u, uvec v) {
	int n = u.n_elem;
	uvec r(n);
	int k = 0;
	for (int i = 0; i < n; i++)
	{
		if (!any(u(i) == v)) {
			r(k) = u(i);
			k++;
		}
	}
	if (k == 0)
		return uvec();
	else
		return r(span(0, k - 1));
}

uvec inicijalno_random_resenje(int n, int p) {
	uvec x(p);
	uvec v = linspace<uvec>(0, n - 1, n);
	x = random_sample(v, p);
	return x;
}
double f_cilja(uvec x, mat C) {
	int p = x.n_elem;
	int n = C.n_cols;
	double f = 0;
	mat Cx = C.cols(x);
	uvec all = linspace<uvec>(0, n - 1, n);
	uvec korisnici = u_diff_v(all, x);
	for (int i = 0; i < n - p; i++)
	{
		int korisnik = korisnici(i);
		double distance = min(Cx.row(korisnik));
		if (distance > f) {
			f = distance;
		}
	}
	return f;

}
//f-ja vraca nasumicnog suseda koji pripada okolini N_k(x)
//elementi N_k(x) se nalaze na Hammingovom rastojanju k od resenja x

uvec shaking(uvec x, int k, int n) {
	uvec x_shake = x;
	int p = x_shake.n_elem;
	uvec all = linspace<uvec>(0, n - 1, n);
	uvec korisnici = u_diff_v(all, x);

	uvec izlaze_iz_x = random_sample(linspace<uvec>(0, p - 1, p), k);//indeksi elemenata koji izlaze iz resenja
	uvec ulaze_u_x = random_sample(korisnici, k);
	for (int i = 0; i < k; i++)
	{
		x_shake(izlaze_iz_x(i)) = ulaze_u_x(i);
	}
	return x_shake;
}
//LS za BVNS i SVNS
//Localna pretraga okoline N_1
//koristi BI strategiju
uvec local_search_N1(uvec x, mat C, int n) {
	uvec x_best = x;
	double f_best = f_cilja(x_best, C);
	int p = x.n_elem;
	int poboljsanje = 1;
	uvec all = linspace<uvec>(0, n - 1, n);
	uvec korisnici = u_diff_v(all, x);
	uvec sused = x;
	double f_sused;

	while (poboljsanje) {
		poboljsanje = 0;
		//pretrazujemo okolinu N_1(x)

		for (int i = 0; i < p; i++)
		{
			for (int j = 0; j < n - p; j++) {
				sused(i) = korisnici(j);
			
				f_sused = f_cilja(sused, C);
				if (f_sused < f_best) {
					poboljsanje = 1;
					f_best = f_sused;
					x_best = sused; //azuriramo najbolje resenje u okolini N1 
				}
			}

			sused = x;
		}
	
		sused = x_best;//prelazimo u najbolje resenje i nastavljamo pretragu njegove okoline
		x = x_best;
		korisnici = u_diff_v(all, x);
	}
	cout << "f_best za LS: " << f_best << endl;
	return x_best;
}
//LS za BVNS i SVNS
//Localna pretraga okoline N_1
//koristi FI strategiju
uvec local_search_N1_FI(uvec x, mat C, int n) {
	uvec x_best = x;
	double f_best = f_cilja(x_best, C);
	int p = x.n_elem;
	int poboljsanje = 1;
	uvec all = linspace<uvec>(0, n - 1, n);
	uvec korisnici = u_diff_v(all, x);
	uvec sused = x;
	double f_sused;

	while (poboljsanje) {
		poboljsanje = 0;
		//pretrazujemo okolinu N_1(x)

		for (int i = 0; i < p; i++)
		{
			for (int j = 0; j < n - p; j++) {
				sused(i) = korisnici(j);

				f_sused = f_cilja(sused, C);
				if (f_sused < f_best) {
					poboljsanje = 1;
					f_best = f_sused;
					x_best = sused; //azuriramo najbolje resenje u okolini N1
					break;
				}
			}
			sused = x;
			if (poboljsanje)
				break;
		}
		sused = x_best;//prelazimo u najbolje resenje i nastavljamo pretragu njegove okoline
		x = x_best;
		korisnici = u_diff_v(all, x);
	}
	cout << "f_best za LS: " << f_best << endl;
	return x_best;
}
uvec BVNS(mat C, int n, int p,int iter_max, int k_max) {

	uvec x = inicijalno_random_resenje(n, p);
	double f = f_cilja(x, C);

	uvec x_best = x;
	double f_best = f;

	int iter = 0;
	uvec x_shake;
	uvec x_ls;
	double f_shake;
	double f_ls;
	int k;
	while (iter < iter_max) {
		k = 1;
		while (k < k_max) {
			x_shake = shaking(x_best, k, n);
			x_ls = local_search_N1(x_shake, C, n);
			f_ls = f_cilja(x_ls, C);
			if (f_ls < f_best) {
				x_best = x_ls;
				f_best = f_ls;
				k = 1;
			}
			else
				k++;
		}
		iter++;
	}
	cout << "f_best za BVNS: " << f_best << endl;
	return x_best;
}
uvec RVNS(mat C, int n, int p, int iter_max, int k_max) {

	uvec x = inicijalno_random_resenje(n, p);
	double f = f_cilja(x, C);

	uvec x_best = x;
	double f_best = f;

	int iter = 0;
	uvec x_shake;
	double f_shake;
		int k;
	while (iter < iter_max) {
		k = 1;
		while (k < k_max) {
			x_shake = shaking(x_best, k, n);
			f_shake = f_cilja(x_shake, C);
			if (f_shake < f_best) {
				x_best = x_shake;
				f_best = f_shake;
				k = 1;
			}
			else
				k++;
		}
		iter++;
	}
	cout << "f_best za RVNS: " << f_best << endl;
	return x_best;
}

int hamming_distance(uvec x1, uvec x2) {
	int p = x1.n_elem;
	int presek=0;
	for (int i = 0; i < p; i++)
	{
		if (any(x1(i) == x2))
			presek++;
	}
	return p - presek;
}

uvec SVNS(mat C, int n, int p, int iter_max, int k_max, double alpha) {

	uvec x = inicijalno_random_resenje(n, p);
	double f = f_cilja(x, C);

	uvec x_best = x;
	double f_best = f;

	int iter = 0;
	uvec x_shake;
	uvec x_ls;
	double f_shake;
	double f_ls;
	int k;
	int hd;
	while (iter < iter_max) {
		k = 1;
		while (k < k_max) {
			x_shake = shaking(x, k, n);
			x_ls = local_search_N1(x_shake, C, n);
			f_ls = f_cilja(x_ls, C);
			if (f_ls < f_best) {
				x_best = x_ls;
				f_best = f_ls;
			}
			hd = hamming_distance(x_ls, x);
			if (f_ls < f + alpha * hd) {
				x = x_ls;
				f = f_ls;
				k = 1;
			}
			else
				k++;
		}
		iter++;
	}
	cout << "f_best za SVNS: " << f_best << endl;
	return x_best;
}
// Pretrega okoline LK_l za VND
//x je tekuce resenje, l je redni broj okoline koja se pretrazuje

uvec pretraga_LK_l(uvec x, int l, mat C) {
	int n = C.n_rows;
	int p = x.n_elem;
	
	uvec all = linspace<uvec>(0, n - 1, n);
	uvec korisnici = u_diff_v(all, x);
	int broj_korisnika = korisnici.n_elem;
	
	umat P = zeros<umat>(p, l);
	//l resenje smestamo u matricu P po kolonama
	//vrednosti f-je cilja za ta resenja pamtimo u vektoru F
	vec F = zeros<vec>(l);

	int i_in;
	int i_out;
	uvec sused = x;
	uvec sused_novi;
	double f = datum::inf;
	for (int k = 0; k < l; k++)
	{
		for (int i = 0; i < p; i++)
		{	
			for (int j = 0; j < broj_korisnika; j++)
			{
				sused(i) = korisnici(j);
				if (f_cilja(sused, C) < f) {
					f = f_cilja(sused, C);
					sused_novi = sused;
					i_in = i;
					i_out = j;
				}
				sused = x;
			}
		}
		//pamtimo generisanog suseda
		P.col(k) = sused_novi;
		F(k) = f;
		//azuriramo za sledecu iteraciju
		f = datum::inf;
		sused = sused_novi;
		korisnici(i_out) = 1000;
		korisnici = korisnici(find(korisnici < 1000));
		broj_korisnika--;


	}

	int ind_min = F.index_min();
	return P.col(ind_min);


}
//VND
void VND(uvec& x_best, double& f_best, int n, int p, mat C, int l_max) {
	//inicijalno resenje
	uvec x = inicijalno_random_resenje(n, p);
	double f = f_cilja(x, C);

	x_best = x;
	f_best = f;

	uvec x_new;
	double f_new;
	int l = 1;

	while (l < l_max) {
		x_new = pretraga_LK_l(x_best, l, C);
		f_new = f_cilja(x_new, C);
		if (f_new < f_best) {
			x_best = x_new;
			f_best = f_new;
			l = 1;
		}
		else
			l++;

	}
	
	cout << "VND f_best: " << f_best<<endl;
}

uvec GVNS(mat C, int n, int p, int iter_max, int k_max, int l_max) {

	uvec x = inicijalno_random_resenje(n, p);
	double f = f_cilja(x, C);

	uvec x_best = x;
	double f_best = f;

	int iter = 0;
	uvec x_shake;
	uvec x_vnd;
	double f_shake=f;
	double f_vnd;
	int k;
	uvec x_new;
	double f_new;
	int l = 1;
	while (iter < iter_max) {
		cout << iter << endl;
		k = 1;
		while (k < k_max) {
			x_shake = shaking(x_best, k, n);
			f_shake = f_cilja(x_shake, C);
			//VND
			l = 1;

			while (l < l_max) {
				x_new = pretraga_LK_l(x_shake, l, C);
				f_new = f_cilja(x_new, C);
				if (f_new < f_shake) {
					x_shake = x_new;
					f_shake = f_new;
					l = 1;
				}
				else {
					l++;
				}
					
			}
			f_vnd = f_shake;
			x_vnd = x_shake;

			if (f_vnd < f_best) {
				x_best = x_vnd;
				f_best = f_vnd;
				k = 1;
			}
			else
				k++;
		}
		iter++;
	}
	cout << "f_best za GVNS: " << f_best << endl;
	return x_best;
}




int main(int argc, char** argv) {
	ifstream data("berlin52.txt");
	int n, p;
	data >> n;
	data >> p;
	mat C = zeros<mat>(n, n);

	for (int i = 0; i < n; i++)
	{
		for (int j = 0; j < n; j++) {
			data >> C(i, j);
		}
	}
	uvec x = inicijalno_random_resenje(n, p);
	//cout << "Inicijalno: " << x << endl;
	//cout << "Shaking: " << shaking(x, 3,n) << endl;
	//cout << local_search_N1(x, C, n) << endl;
	//cout << BVNS(C,n,p,100,round(p/2))<<endl;
	//cout << RVNS(C, n, p, 1500, round(p / 2)) << endl;
	//cout << SVNS(C, n, p, 50, round(p / 2),0.3) << endl;
	uvec x_best;
	double f_best;
	VND(x_best, f_best, n, p, C, n-p);
	cout << GVNS(C, n, p, 500, round(p / 2),2) << endl;

}