# Program Name: NSGA-II.py
# Description: This is a python implementation of Prof. Kalyanmoy Deb's popular NSGA-II algorithm
# Author: Haris Ali Khan 
# Supervisor: Prof. Manoj Kumar Tiwari

# Importing required modules
import math
import random
import matplotlib.pyplot as plt
from NYC.pylouvain import *
import NYC.variance as variance


# First function to optimize(区域数量方差)
def function1(areasize_weight, distance_weight):
    filename = "D:\\JetBrains\\PcharmWorkSpace\\Paper2\\NYC\\data\\relationMatrixNYCForCommunity.csv"
    pyl = PyLouvain.from_file(filename)
    partition, q = pyl.apply_method(areasize_weight, distance_weight, 1)
    result  = variance.get_areasize_variance(partition)
    print("areaSize方差为：", result)
    print("模块度1为：", q)
    return result


# Second function to optimize
def function2(areasize_weight, distance_weight):
    filename = "D:\\JetBrains\\PcharmWorkSpace\\Paper2\\NYC\\data\\relationMatrixNYCForCommunity.csv"
    pyl = PyLouvain.from_file(filename)
    partition, q = pyl.apply_method(areasize_weight, distance_weight, 1)
    distanceSum, distance_variance = variance.get_distance_variance(partition)
    print("区域调度距离方差为：", distance_variance)
    print("模块度2为：", q)
    return distance_variance


# Function to find index of list
def index_of(a, list):
    for i in range(0, len(list)):
        if list[i] == a:
            return i
    return -1


# Function to sort by values
def sort_by_values(list1, values):
    sorted_list = []
    while (len(sorted_list) != len(list1)):
        if index_of(min(values), values) in list1:
            sorted_list.append(index_of(min(values), values))
        values[index_of(min(values), values)] = math.inf
    return sorted_list


# Function to carry out NSGA-II's fast non dominated sort
def fast_non_dominated_sort(values1, values2):
    S = [[] for i in range(0, len(values1))]
    front = [[]]
    n = [0 for i in range(0, len(values1))]
    rank = [0 for i in range(0, len(values1))]

    for p in range(0, len(values1)):
        S[p] = []
        n[p] = 0
        for q in range(0, len(values1)):
            if (values1[p] < values1[q] and values2[p] < values2[q]) or (
                    values1[p] <= values1[q] and values2[p] < values2[q]) or (
                    values1[p] < values1[q] and values2[p] <= values2[q]):
                if q not in S[p]:
                    S[p].append(q)
            elif (values1[q] < values1[p] and values2[q] < values2[p]) or (
                    values1[q] <= values1[p] and values2[q] < values2[p]) or (
                    values1[q] < values1[p] and values2[q] <= values2[p]):
                n[p] = n[p] + 1
        if n[p] == 0:
            rank[p] = 0
            if p not in front[0]:
                front[0].append(p)

    i = 0
    while (front[i] != []):
        Q = []
        for p in front[i]:
            for q in S[p]:
                n[q] = n[q] - 1
                if (n[q] == 0):
                    rank[q] = i + 1
                    if q not in Q:
                        Q.append(q)
        i = i + 1
        front.append(Q)

    del front[len(front) - 1]
    return front


# Function to calculate crowding distance
def crowding_distance(values1, values2, front):
    distance = [0 for i in range(0, len(front))]
    sorted1 = sort_by_values(front, values1[:])
    sorted2 = sort_by_values(front, values2[:])
    distance[0] = 4444444444444444
    distance[len(front) - 1] = 4444444444444444
    for k in range(1, len(front) - 1):
        distance[k] = distance[k] + (values1[sorted1[k + 1]] - values2[sorted1[k - 1]]) / (max(values1) - min(values1))
    for k in range(1, len(front) - 1):
        distance[k] = distance[k] + (values1[sorted2[k + 1]] - values2[sorted2[k - 1]]) / (max(values2) - min(values2))
    return distance


# Function to carry out the crossover
def crossover1(a, b): #对AreaSize进行交叉变异
    r = random.random()
    if r > 0.5:
        return mutation1((a + b) / 2)
    else:
        return mutation1((a - b) / 2)


# Function to carry out the mutation operator
def mutation1(solution):
    mutation_prob = random.random()
    if mutation_prob < 1:
        solution = min_areaweight + (max_areaweight - min_areaweight) * random.random()
    return solution

# Function to carry out the crossover
def crossover2(a, b): #对distance进行交叉变异
    r = random.random()
    if r > 0.5:
        return mutation2((a + b) / 2)
    else:
        return mutation2((a - b) / 2)


# Function to carry out the mutation operator
def mutation2(solution):
    mutation_prob = random.random()
    if mutation_prob < 1:
        solution = min_distanceweight + (max_distanceweight - min_distanceweight) * random.random()
    return solution


# Main program starts here
pop_size = 100
max_gen = 200 #最大迭代次數

# Initialization
min_areaweight = 1
max_areaweight = 10
min_distanceweight = 1
max_distanceweight = 10

areaweight = [min_areaweight + (max_areaweight - min_areaweight) * random.random() for i in range(0, pop_size)]
distanceweight = [min_distanceweight + (max_distanceweight - max_distanceweight) * random.random() for i in
                  range(0, pop_size)]
gen_no = 0
while (gen_no < max_gen):
    function1_values = [function1(areaweight[i], distanceweight[i]) for i in range(0, pop_size)]
    function2_values = [function2(areaweight[i], distanceweight[i]) for i in range(0, pop_size)]
    non_dominated_sorted_solution = fast_non_dominated_sort(function1_values[:], function2_values[:])
    # print("The best front for Generation number ", gen_no, " is")
    # for valuez in non_dominated_sorted_solution[0]:
    #     print(round(solution[valuez], 3), end=" ")
    # print("\n")
    crowding_distance_values = []
    for i in range(0, len(non_dominated_sorted_solution)):
        crowding_distance_values.append(
            crowding_distance(function1_values[:], function2_values[:], non_dominated_sorted_solution[i][:]))
    areaweight2 = areaweight[:]
    distanceweight2 = distanceweight[:]
    # Generating offsprings
    while (len(areaweight2) != 2 * pop_size):
        a1 = random.randint(0, pop_size - 1)
        b1 = random.randint(0, pop_size - 1)
        areaweight2.append(crossover1(areaweight[a1], areaweight[b1]))
    while (len(distanceweight2) != 2 * pop_size):
        a1 = random.randint(0, pop_size - 1)
        b1 = random.randint(0, pop_size - 1)
        distanceweight2.append(crossover2(distanceweight[a1], distanceweight[b1]))
    # 精英选择策略
    function1_values2 = [function1(areaweight2[i], distanceweight2[i]) for i in range(0, 2 * pop_size)]
    function2_values2 = [function1(areaweight2[i], distanceweight2[i]) for i in range(0, 2 * pop_size)]
    non_dominated_sorted_solution2 = fast_non_dominated_sort(function1_values2[:], function2_values2[:])
    crowding_distance_values2 = []
    for i in range(0, len(non_dominated_sorted_solution2)):
        crowding_distance_values2.append(
            crowding_distance(function1_values2[:], function2_values2[:], non_dominated_sorted_solution2[i][:]))
    new_solution = []
    for i in range(0, len(non_dominated_sorted_solution2)):
        non_dominated_sorted_solution2_1 = [
            index_of(non_dominated_sorted_solution2[i][j], non_dominated_sorted_solution2[i]) for j in
            range(0, len(non_dominated_sorted_solution2[i]))]
        front22 = sort_by_values(non_dominated_sorted_solution2_1[:], crowding_distance_values2[i][:])
        front = [non_dominated_sorted_solution2[i][front22[j]] for j in
                 range(0, len(non_dominated_sorted_solution2[i]))]
        front.reverse()
        for value in front:
            new_solution.append(value)
            if (len(new_solution) == pop_size):
                break
        if (len(new_solution) == pop_size):
            break
    areaweight = [areaweight2[i] for i in new_solution]
    distanceweight = [distanceweight2[i] for i in new_solution]
    gen_no = gen_no + 1

# Lets plot the final front now
function1 = [i for i in function1_values]
function2 = [j for j in function2_values]
print("最优areaSize的权重为:", areaweight[0])
print("最优distance的权重为:", distanceweight[0])
plt.xlabel('Function 1', fontsize=15)
plt.ylabel('Function 2', fontsize=15)
plt.scatter(function1, function2)
plt.show()