#include <iostream>
#include <vector>
#include <cmath>
#include <ctime>
#include <random>
#include "timing.hpp"
#include "graph.h"
#include "DEBUG.h"
using namespace std;
#define APPLY_ELITISM YES
typedef unsigned long long ull;
#define DECAY_RATE 0.9
enum META_HEURISTIC_INFO{
    DELTA = 4,
    BETA = 6,
    POP_SIZE = 100,
    NUM_LOOP = 256
};
#define Q 1.7
#define PUNISHMENT 0.3
#define N 55
int k = 34;
#define S 0
#define T (N - 1)
Graph G;
string CSVFILE(){
    return "undirected_weighted_graph" + to_string(N)+".csv";
 }
int get_tour_length(Graph& G, Tour& t){
    int cost = 0;
    int n = t.size();
    for(int i = 1; i < n; i++){
        vertex u = t[i-1], v = t[i];
        cost += G[u][v];
        if(v == T) break;
    }
    return cost;
}
vector<vector<double>> pheromone, heuristic;
random_device rd;
static default_random_engine eng(rd());
static uniform_real_distribution<> getReal(0, 1);
class Ant{
    public:
    int location;
    vector<bool> node_visited;
    Tour trace;
    bool tour_complete;
    Ant(){
        node_visited = vector<bool>(N , false);
        location = S;
        node_visited[S] = true;
        trace = {S};
        tour_complete = false;
    }
    double f(double tao, double neta){
        return pow(tao, DELTA) + pow(neta, BETA);
    }
    vertex pick_next_node(){
        vector<vertex> go;
        for(int i = 0; i < N; i++)
            if(i != this->location && !node_visited[i]) 
                go.push_back(i);
        int sz = go.size();
        if(sz == 0 ){
            // cout<<"warning\n";
            return -1;
        }
        //an implementation of roulette wheel.
        vector<double> possiblities;
        double sum_possi = 0.0;
        for(vertex& v : go){
            double d1 = pheromone[location][v];
            double d2 = heuristic[location][v];
            possiblities.push_back(f(d1, d2));
        }

        for(double& p : possiblities) 
            sum_possi += p;
        for(double& p : possiblities) 
            p /= sum_possi;

        double rnd = getReal(eng);
        for(int i = 0; i < (N - 1); i++){
            rnd -= possiblities[i];
            if (rnd <= 0)
                return go[i];
        }
    }
};
vector<Ant> population;
Tour best_tour = {};
int optimal_cost = INF;
void set_init_info(){
    //ant population init.
    population = vector<Ant>(POP_SIZE);
    for(Ant& ind : population)
        ind = Ant();
    //pheromone array init.
    pheromone = vector<vector<double>>(N, vector<double>(N, 0.0));
    heuristic = pheromone;
    for (int i = 0; i < N; i++){
        for (int j = 0; j < N; j++){//边长倒数{
                if(i == j) continue;
                else 
                  pheromone[i][j] = (G[i][j] == 0) ? 0.0 : (100.0/(double)G[i][j]);
            }
    }
    // heuristic info array init. begins with compuation of MEAN
    double H0 = 0.0;
    int edgeWsum = 0;
    int n_edge = 0;
    for(int i = 0; i < N; i++){
        for (int j = 0; j < N; j++){
            if(G[i][j] != 0) {
                edgeWsum += G[i][j];
                n_edge ++;
            }
        }
    }
    double mean = ((double)edgeWsum) / n_edge;
    H0 = 1.0 /(N * mean);
    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++)
            heuristic[i][j] *= H0;

}
void update_pheromone(){
    //the pheromone decays over time.
    // cout<<"updating\n";
    for(int i = 0; i < N; i++)
        for(int j = 0; j < N; j++)
            pheromone[i][j] *= DECAY_RATE;

    for(Ant& a : population){
        int cost = get_tour_length(G, a.trace);
        double inv = 1.0 /(double)(cost);
        int save = optimal_cost;
        #ifdef APPLY_ELITISM
        if(cost >= save) 
            inv *= PUNISHMENT;
        else{
            inv *= Q;  
            optimal_cost = min(cost, optimal_cost);
        }
        #endif 
        int n = a.trace.size();
        for(int i = 1; i < n; i++){
            vertex u = a.trace[i-1], v = a.trace[i];
            pheromone[u][v] += inv;
        }
    }
    // #ifdef APPLY_ELITISM
    //     int L_best = optimal_cost;
    //     double invv = 100.0 / (double)L_best;
    //     invv *= Q;
    //     for(int i = 1; i < best_tour.size(); i++){
    //         vertex u = best_tour[i-1], v = best_tour[i];
    //         pheromone[u][v] += (invv);
    //     }
    // }
    // #endif
}

void construct_path(){
    for(Ant& ind : population){
        for(int step = 0; step < N; step++){
            if(ind.tour_complete) continue;
            vertex node =  ind.pick_next_node(); 
            if (node == -1) {
                continue;
            }
            ind.location = node;
            ind.node_visited[node] = true;
            ind.trace.push_back(node);
            //if an ant accidentally bulids a path.
            if (node == T){
            //comprate with existing solution.
                if(ind.trace.size() >= k){ // has at least k nodes.
                    // int cost = get_tour_length(G, ind.trace);
                    ind.tour_complete = true;
                    if(optimal_cost == INF)
                    // if(optimal_cost > cost){
                        optimal_cost = get_tour_length(G, ind.trace);
                        best_tour = ind.trace;
                    // }
                }else{
                    --step;
                    ind.trace.pop_back();
                    ind.node_visited[T] = false;
                    ind.location = ind.trace.back();
                }
            }
        }
    }
}
void boot(){
    for(int x = 0; x < NUM_LOOP; x++){
        set_init_info();
        construct_path();
        update_pheromone();
    }
}

int main(){
    // cout<<"Loading Graph from csv\n";
    G = read_from_csv(N, CSVFILE());
    // cout<<"Load complete.\n";
    TIMER_START(x);
    boot();
    TIMER_STOP(x);
    cout<<"time elpased :"<<TIMER_SEC(x)<<endl;
    cout<<"soultion :"<<optimal_cost<<endl;
    print_tour(best_tour);
    return 0;
}