# ---
# jupyter:
#   jupytext:
#     formats: ipynb,md,py:percent
#     text_representation:
#       extension: .py
#       format_name: percent
#       format_version: '1.3'
#       jupytext_version: 1.16.1
#   kernelspec:
#     display_name: Python 3 (ipykernel)
#     language: python
#     name: python3
# ---

# %% [markdown]
"""
# Evolutionary Computation

**Master in Artificial Intelligence, UVigo, UdC, USC**

Academic year 2025/26
"""

# %% [markdown]
# Let us import all necessary modules:

# %%
import math
import matplotlib.pyplot as plt
import random
from sko.GA import GA_TSP


# %% [markdown]
# First we define some helper functions we will need in the on-going.

# %% [markdown]
# We want to sort in ascending order. Hence, a pair is sorted whenever the one at smaller index is less or equal to the one at larger index.

# %%
def PairIsSorted(v,i,j):
  assert i<j
  return v[i]<=v[j]


# %% [markdown]
# An entire sequence is sorted whenever all consecutive pairs are sorted.

# %%
def IsSorted(v):
  for i in range(len(v)-1):
    if not PairIsSorted(v,i,i+1):
      return False
  return True


# %% [markdown]
# We sort a pair by exchange whenever the pair was not sorted beforehand.

# %%
def PairSort(v,i,j):
  if not PairIsSorted(v,i,j):
    v[i],v[j]=v[j],v[i]


# %% [markdown]
# To select a random pair, we just draw two different indices in the corresponding range. Note that we always arrange the index pair such that the first one is the smallest one.

# %%
def RandomPair(v):
  i,j=random.sample(range(len(v)),2)
  i,j=(j,i) if j<i else (i,j)
  return i,j


# %% [markdown]
# Our objective function as given in the lecture: we multiply each value in the sequence with its index. We return the negative value as we want to minimize the objective function. Observe, this objective function works for all numeric values (not only for integers).

# %%
def SortObjective(v):
  f=0
  for i in range(len(v)):
    f+=v[i]*(i+1)
  return -f


# %% [markdown]
# So we can plot the evolution of the value of the objective function and store an image file.

# %%
def PlotFobj(fo,filename=None,show=False):
  fig,ax=plt.subplots(figsize=(12.8,9.6))
  ax.plot(fo,linewidth=3)
  for label in (ax.get_xticklabels()+ax.get_yticklabels()):
    label.set_fontsize(16)
  plt.xlabel('n pair-comparisons', fontsize=22)
  plt.ylabel('sort error', fontsize=22)
  if filename!=None:
    plt.savefig(filename,transparent=True,bbox_inches='tight')
  if show==True:
    plt.show()


# %% [markdown]
# Given a sequence of values, we compute the objective function and show whether the data is sorted or not.

# %%
def ShowSorted(w):
  print("objective",SortObjective(w),
    "sorted\n" if IsSorted(w) else "not sorted\n"
  )


# %% [markdown]
# Here we build our data set: just a vector of 100 entries as random integer values.

# %%
v = random.sample(range(100),100)
maxv=SortObjective(v)
print("initial sort objective",maxv)
rounds = 15*int(math.log2(len(v))+0.5)*len(v)
print("rounds",rounds)


# %% [markdown]
# We first implement the well-known bubble sort algorithm. It's a deterministic algorithm with worst case quadratic runtime (whenever the initial sequence is sorted in descending order) and best case linear runtime (whenever the initial sequence is already sorted). 

# %%
def BubbleSort(v):
  fo=[]
  while not IsSorted(v):
    for i in range(len(v)-1):
      PairSort(v,i,i+1)
    # we store fobj every round (n pair comparisons)
    fo.append(SortObjective(v))
  return fo


# %% [markdown]
# The test piece for bubble sort. We draw the convergence of the objetive function as relative error according to the achievable optimum.

# %%
print("BubbleSort")
w = v.copy()
fo=BubbleSort(w)
fobj_min=SortObjective(w)
fo_bs=[ (x-fobj_min)/-fobj_min for x in fo ]
PlotFobj(fo_bs,"BubbleSort.png")
ShowSorted(w)


# %% [markdown]
# A Monte Carlo sort just runs a certain number of rounds and sorts a random pair.

# %%
def MonteCarloSort(v,rounds):
  fo=[]
  m=0
  for j in range(0,rounds):
    i,j=RandomPair(v)
    PairSort(v,i,j); m+=1
    # we store fobj every n pair comparisons
    if m==len(v):
      fo.append(SortObjective(v)); m=0
  return fo


# %% [markdown]
# The test piece for Monte Carlo sort. We draw the convergence of the objective function as relative error according to the achievable optimum.

# %%
print("MonteCarlo sort")
w = v.copy()
fo=MonteCarloSort(w,rounds)
fo_mc=[ (x-fobj_min)/-fobj_min for x in fo ]
PlotFobj(fo_mc,"MonteCarloSort.png")
ShowSorted(w)


# %% [markdown]
# A Las Vegas sort sorts random pairs until the sequence is sorted.

# %%
def LasVegasSort(v):
  fo=[]
  m=0
  while not IsSorted(v):
    i,j=RandomPair(v)
    PairSort(v,i,j); m+=1
    # we store fobj every n pair comparisons
    if m==len(v):
      fo.append(SortObjective(v)); m=0
  return fo


# %% [markdown]
# The test piece for Las Vegas sort. We draw the convergence of the objective function as relative error according to the achievable optimum.

# %%
print("LasVegas sort")
w = v.copy()
fo=LasVegasSort(w)
fo_lv=[ (x-fobj_min)/-fobj_min for x in fo ]
PlotFobj(fo_lv,"LasVegasSort.png")
ShowSorted(w)


# %% [markdown]
# A Monte Carlo sort with Las Vegas condition just runs a certain number of rounds and sorts a random pair, but stops as well, whenever the sequence becomes sorted.

# %%
def MonteCarloWithLasVegasConditionSort(v,rounds):
  fo=[]
  m=0
  while not IsSorted(v) and rounds>0:
    rounds-=1
    i,j=RandomPair(v)
    PairSort(v,i,j); m+=1
    # we store fobj every n pair comparisons
    if m==len(v):
      fo.append(SortObjective(v)); m=0
  return fo


# %% [markdown]
# The test piece for Monte Carlo sort with Las Vegas condition. We draw the convergence of the objective function as relative error according to the achievable optimum.

# %%
print("MonteCarlo with LasVegas condition sort")
w = v.copy()
fo=MonteCarloWithLasVegasConditionSort(w,rounds)
fo_mclv=[ (x-fobj_min)/-fobj_min for x in fo ]
PlotFobj(fo_mclv,"MonteCarloLasVegasSort.png")
ShowSorted(w)


# %% [markdown]
# We implement the well-known merge sort algorithm in its iterative form. It's a deterministic algorithm with quasi-linear runtime independend from the input data.

# %%
def Merge(v,w,left,middle,right):
  i,j=left,middle
  for k in range(left,right):
    if j>=right or ( i<middle and PairIsSorted(v,i,j) ):
      w[k]=v[i]; i+=1
    else:
      w[k]=v[j]; j+=1

def MergeSort(w):
  fo=[]
  s=1
  n=len(w)
  while s<n:
    v=w[:]
    for left in range(0,n,2*s):
      Merge(v,w,left,left+s,min(left+2*s,n))
    s*=2
    # we store fobj every round (n pair comparisons)
    fo.append(SortObjective(w))
  return fo


# %% [markdown]
# The test piece for merge sort. We draw the convergence of the objective function as relative error according to the achievable optimum.

# %%
print("Merge sort")
w = v.copy()
fo=MergeSort(w)
fo_ms=[ (x-fobj_min)/-fobj_min for x in fo ]
PlotFobj(fo_ms,"MergeSort.png")
ShowSorted(w)


# %% [markdown]
# Now we use a genetic algorithm designed for TSP to sort. We only need to specify our special objective function. The genetic algorithm will intend to return a sequence of indices that minimizes that function, hence, the underlying data will be sorted!

# %%
def GASort(w):
  def fobj(routine):
    n, = routine.shape
    return -sum([w[int(routine[i%n])]*(i+1) for i in range(n)])
  ga=GA_TSP(
    func=fobj,n_dim=len(w),size_pop=100,max_iter=400,prob_mut=1
  )
  best_points,best_val=ga.run()
  v=w.copy()
  for i in range(len(best_points)):
    w[i]=v[int(best_points[i])]
  return ga.generation_best_Y


# %% [markdown]
# The test piece for the sorting with a genetic algorithm. We draw the convergence of the objective function as relative error according to the achievable optimum.

# %%
print("Genetic sort")
w = v.copy()
fo=GASort(w)
fo_ga=[ (x-fobj_min)/-fobj_min for x in fo ]
PlotFobj(fo_ga,"GeneticSort.png")
ShowSorted(w)

# %% [markdown]
# Finally, we draw all convergence data of the different algorithms into one plot in order to be able to compare the behavior.

# %%
fig,ax=plt.subplots(figsize=(12.8,9.6))
ax.plot(fo_ga,linewidth=3,label="genetic algorithm")
ax.plot(fo_bs,linewidth=3,label="bubble sort")
ax.plot(fo_lv,linewidth=3,label="las vegas sort")
ax.plot(fo_mc,linewidth=3,label="monte carlo sort")
ax.plot(fo_ms,linewidth=3,label="merge sort")
ax.legend(fontsize=16)
for label in (ax.get_xticklabels()+ax.get_yticklabels()):
  label.set_fontsize(16)
plt.xlabel('n pair-comparisons', fontsize=22)
plt.ylabel('sort error', fontsize=22)
plt.savefig("sort_fobj",transparent=True,bbox_inches='tight')
plt.show()

# %%
