#!/usr/bin/env python
# coding: utf-8

# <div class="alert alert-block alert-info">
# <center><font size=5><b>TP n°13 (algorithmique, complexité) : Éléments de correction</b></font></center></div>

# $\require{color}$
# $\require{xcolor}$
# $\newcommand{\myfbox}[1]{\fcolorbox{red}{white}{$\textrm{#1}$}}$
# $\require{stmaryrd}$
# $\newcommand{\ient}{[\![}$
# $\newcommand{\fient}{]\!]}$
# $\newcommand{\R}{\mathbb R}$
# $\newcommand{\K}{\mathbb K}$
# $\newcommand{\N}{\mathbb N}$
# $\newcommand{\C}{\mathbb C}$
# $\newcommand{\sh}{\operatorname{sh}}$
# $\newcommand{\ch}{\operatorname{ch}}$
# $\newcommand{\id}{\operatorname{Id}}$
# $\newcommand{\mat}{\operatorname{mat}}$
# $\newcommand{\sp}{\operatorname{Sp}}$
# $\newcommand{\In}{\operatorname{I}}$
# $\newcommand{\vect}{\operatorname{Vect}\ }$
# $\newcommand{\rg}{\operatorname{rg}}$
# $\newcommand{\tr}{\operatorname{Tr}}$
# $\newcommand{\dis}{\displaystyle}$
# $\renewcommand{\Im}{\operatorname{Im}}$
# $\newcommand{\im}{\operatorname{Im}}$
# $\newcommand{\bordermatrix}[3]{ \begin{matrix} \begin{matrix}#1\end{matrix} & \\  #3  & \hspace{-1em} \begin{matrix}#2\end{matrix} \end{matrix}}$
# $\newcommand{\fonction}[5]{#1\ \colon \left\{\begin{array}{ccl}#2 & \longrightarrow & #3\\#4 & \longmapsto & #5\end{array}\right.}$
# $\newcommand{\revdots}{\mathinner{\mkern1mu\raise1pt{\kern7pt\hbox{.}}\mkern2mu\raise4pt\hbox{.}\mkern2mu\raise7pt\hbox{.}\mkern1mu}}$
# $\newcommand{\tendvers}[2]{\xrightarrow[#1\to#2]{}}$
# $\newcommand{\q}[1]{\mathbf{Q#1}\  \triangleright}$

# In[1]:


import random as rd
import numpy as np
from time import perf_counter as toc
import matplotlib.pyplot as plt
from scipy import stats


# <div style="border-radius:15px;-webkit-border-radius:15px;-moz-border-radius:15px;width:37%;float:center;background-color:#CBD6E7;margin-left:0;padding:10px;text-align:center;">
#     <font size=3><b>Exercice 1 : tri par sélection</b></font>
# </div>

# $\q{1}$ Il n'y a que des boucles ```for``` dans cet algorithme. Cela assure sa $\myfbox{terminaison.}$

# $\q2$ On note $\mathscr P_{i} \colon$ "La liste ```L[0 : i]``` est triée et tous ses éléments sont plus petits ou égaux que les éléments de la liste ```L[i : n]```". Montrons qu'il s'agit d'un invariant de la boucle ```for``` dans l'algorithme de tri par sélection :
# 
# - au début de la première itération, quand ```i=0```, la liste ```L[0 : i]``` est vide, donc triée, et tous ses éléments sont plus petits ou égaux que les éléments de la liste ```L[0 : n]```, qui est en fait la liste ```L```. $\mathscr P_0$ est donc vérifiée.
# 
# - si $\mathscr P_{i}$ est vérifiée au début d'une itération, la liste ```L[0 : i]``` est triée et tous ses éléments sont plus petits ou égaux que les éléments de la liste ```L[i : n]```. Puisque ```j``` est l'indice du plus petit élément de la liste ```L[i : n]``` et que l'on vient le placer en position ```i``` grâce à l'échange ```L[i], L[j] = L[j], L[i]```, en fin d'itération, la liste ```L[0 : i+1]``` est triée et tous ses éléments sont plus petits ou égaux que les éléments de la liste ```L[i+1 : n]```.

# À la fin de la dernière itération, la liste ```L[0 : n-1]``` est triée et tous ses éléments sont plus petits ou égaux que les éléments de la liste ```L[n-1 : n]```, qui est juste composée d'un unique élément ```L[n-1]```.
# 
# Ainsi, $\myfbox{la liste L que l'on renvoie est triée}$.

# $\q3$ Dans tous les cas (et donc en particulier dans le pire des cas), il y a $n-1$ itérations dans la boucle ```for``` lors de l'appel  ```select(L)``` ( pour $i\in\ient 0,n-2\fient$). Chacune de ces itérations donne lieu à un appel à ```minpos(L,i)```, qui effectue $n-i-1$ comparaisons (une par itération pour chacun des $k\in\ient i+1,n-1\fient$).
# 
# Il y a donc en tout :
# $$
# \sum_{i=0}^{n-2} (n-i-1) = \sum_{i=1}^{n-1} i = \frac{n(n-1)}{2} = \myfbox{$O\big(n^2\big)$} \ \text{comparaisons.}
# $$

# $\q4$ On reprend l'algorithme donné dans l'énoncé :

# In[2]:


def minpos(L, i):
    j = i
    m = L[i]
    for k in range(i + 1, len(L)):
        if L[k] < m:
            j = k
            m = L[k]
    return j


# In[3]:


def select(L):
    n = len(L)
    for i in range(n - 1):
        j = minpos(L, i)
        L[i], L[j] = L[j], L[i]
    return None


# On reprend également la fonction ```gen(n, N)``` vue dans le dernier TP :

# In[4]:


def gen(n, N):
    return [rd.randint(0, N) for i in range(n)]


# On peut maintenant écrire le script demandé :

# In[5]:


X = [np.log(i) for i in range(500, 1001, 5)]
Y = []
for i in range(500, 1001, 5):
    L = gen(i, i)
    t1 = toc()
    select(L)
    t2 = toc()
    Y.append(np.log(t2 - t1))

X = X[10:]
Y = Y[10:]
fit = stats.linregress(X, Y)
a, b, r = fit[:3]  # on ne garde que les 3 premiers elements de la liste fit
line = [a * X[i] + b for i in range(len(X))]
plt.plot(X, Y, 'o', X, line, 'r-')
plt.show()


# In[6]:


a, b


# On a donc en particulier $a\simeq 2$, d'où :
# $$ \ln(Y) \simeq 2\ln(X) + b$$
# donc :
# $$ Y \simeq e^{b} \times X^2$$
# 
# Si l'on considère que la complexité est proportionnelle au temps de calcul, on obtient une complexité (en moyenne) quadratique :
# $$
# \myfbox{$C(n) \simeq \text{cste}\times n^2$}
# $$

# <div style="border-radius:15px;-webkit-border-radius:15px;-moz-border-radius:15px;width:37%;float:center;background-color:#CBD6E7;margin-left:0;padding:10px;text-align:center;">
#     <font size=3><b>Exercice 2 : tri par insertion</b></font>
# </div>

# $\q1$ ```i``` est un variant de la boucle ```while``` : il s'agit en effet d'une quantité à valeurs entières, strictement décroissante et minorée par $0$. Ceci assure la terminaison de la boucle ```while```. Puisque l'autre boucle est une boucle ```for```, cela assure la terminaison de l'algorithme.

# $\q2$ On note $\mathscr P_k\colon$ "La liste ```L[0 : k]``` est triée". Montrons qu'il s'agit d'un invariant de la boucle ```for``` :
# - au début de la première itération, quand ```k=1```, la liste ```L[0 : k]``` ne possède qu'un seul élément, donc est triée.
# - si $\mathscr P_k$ est vérifiée au début de la $k^{\text{ème}}$ itération, alors la liste ```L[0 : k]``` est triée. La boucle ```while``` permet alors d'insérer ```L[k]``` dans sa bonne position dans la liste ```L[0 : k]```. La liste ```L[0 : k+1]``` est donc bien triée à la fin de l'itération.

# $\q3$ En termes de nombres de comparaisons, le cas le pire est celui où la liste ```L``` est déjà triée par ordre décroissant. Il y a alors :
# - $n-1$ itérations sur $k$ dans la boucle ```for``` (pour $k\in\ient 1,n-1\fient$)
# - à chaque itération sur $k$, il y a $k$ itérations dans la boucle ```while``` (la condition ```L[i - 1] > elt``` étant toujours vérifiée), chacune des itérations donnant lieu à une comparaison (pour le ```if```).
# 
# Le nombre total de comparaisons est donc :
# $$
# \sum_{k=1}^{n-1} k = \frac{n(n-1)}{2} = \myfbox{$O\big(n^2\big)$}
# $$

# In[7]:


def insertion(L):
    n = len(L)
    for k in range(1, n):
        i = k
        elt = L[k]
        while i > 0 and L[i - 1] > elt:
            L[i] = L[i - 1]
            i -= 1
        L[i] = elt
    return L


# In[8]:


X = [i for i in range(500,3001,50)]
Y1 = []
Y2 = []
for i in range(500,3001,50):
    L = list(range(i))
    t1 = toc()
    insertion(L)
    t2 = toc()
    Y1.append((t2-t1))
    L.reverse() # on inverse la liste 
    t1 = toc()
    insertion(L)
    t2 = toc()
    Y2.append((t2-t1))
    
    
X = X[10:]
Y1 = Y1[10:]
Y2 = Y2[10:]
plt.plot(X, Y1, 'bo', label="liste triée par ordre croissant")
plt.legend()
plt.show()
plt.plot(X, Y2, 'ro',label="liste triée par ordre décroissant")
plt.legend()
plt.show()


# **Explication :** On a déjà vu que dans le cas le pire (liste triée par ordre décroissant); la complexité est en $O(n^2)$. Dans le meilleur des cas, celui d'une liste déjà triée par ordre croissant, le nombre de comparaisons effectuées est en $O(n)$ : il n'y a qu'une itération dans chaque boucle ```while```. 

# <div style="border-radius:15px;-webkit-border-radius:15px;-moz-border-radius:15px;width:37%;float:center;background-color:#CBD6E7;margin-left:0;padding:10px;text-align:center;">
#     <font size=3><b>Exercice 3 : exponentiation rapide</b></font>
# </div>

# $\q1$ On peut proposer l'algorithme naïf suivant :

# In[9]:


def naivepow(a, n):
    r = 1
    for k in range(n):
        r *= a
    return r


# In[10]:


naivepow(3, 4)


# $\q2$ Il y a $n$ itérations dans la boucle ```for```, et chaque itération donne lieu à $1$ multiplication.
# 
# Le nombre de multiplications effectuées est donc $n$.

# $\q3$ On implémente en Python l'algorithme proposé :

# In[11]:


def fastpow(a, n):
    r = 1
    while n > 0:
        if n % 2 == 0:
            a *= a
            n = n // 2
        else:
            r *= a
            n -= 1
    return r


# In[12]:


fastpow(3,4)


# $\q4$ Soit $k = \left\lfloor \log_{2}(n)\right\rfloor$. Alors $k$ vérifie :
# $$ 2^k \leqslant n < 2^{k+1}$$
# 
# Au cours de 2 itérations successives de la boucle ```while```, l'entier $n$ est au moins divisé par $2$.
# 
# Ainsi, si $p$ est le nombre d'itérations dans la boucle ```while```, alors $\frac{p}{2}$ est inférieur au nombre de fois où l'on peut diviser $n$ par $2$ avant d'obtenir un résultat nul, *i.e.* :
# $$ p \leqslant 2(k+1)$$
# Finalement :
# $$ p = O(k)$$
# Chaque itération donnant lieu à une multiplication, le nombre total de multiplications effectuées est en $O(k) = \myfbox{$O\big(\log_2(n)\big)$.}$

# $\q5$ On l'adapte en une fonction récursive. Le principe est le suivant :
# $$\texttt{recfastpow(a, n)}= \begin{cases}
# \texttt{recfastpow(a*a, n//2)}\quad\text{si $n$ est pair}\\
# a\times\texttt{recfastpow(a, n-1)}\quad\text{sinon}
# \end{cases}$$
# avec le cas de base :
# $$
# \texttt{recfastpow(a,0)} = 1
# $$

# In[13]:


def recfastpow(a, n):
    if n == 0:
        return 1
    else:
        if n % 2 == 0:
            return recfastpow(a*a, n // 2)
        else:
            return a * recfastpow(a, n - 1)


# In[14]:


recfastpow(3,4)


# $\q6$ On compare les temps d'exécution des 3 versions pour le calcul de $3^n$ (j'ai même rajouté la comparaison avec la fonction ```pow()``` de Python) :

# In[15]:


N = 500
Yn = [] # implémentation naive
Yf = [] # exponentiation rapide itérative
Yr = [] # exponentiation rapide récursive
Yp = [] # fonction native en Python
X = list(range(0,N,50))
for n in X:
    t1 = toc()
    naivepow(3,n)
    t2 = toc()
    Yn.append((t2-t1))
    t1 = toc()
    fastpow(3,n)
    t2 = toc()
    Yf.append((t2-t1))
    t1 = toc()
    recfastpow(3,n)
    t2 = toc()
    Yr.append((t2-t1))
    t1 = toc()
    pow(3,n)
    t2 = toc()
    Yp.append((t2-t1))
    
plt.plot(X, Yn, 'o', label="méthode naïve")
plt.plot(X, Yf, 'o',label="exponentiation rapide")
plt.plot(X, Yr, 'o',label="exponentiation rapide récursive")
plt.plot(X, Yp, 'o',label="fonction pow() de Python")
plt.legend()
plt.show()


# <div style="border-radius:15px;-webkit-border-radius:15px;-moz-border-radius:15px;width:37%;float:center;background-color:#CBD6E7;margin-left:0;padding:10px;text-align:center;">
#     <font size=3><b>Exercice 4 : méthode de Hörner</b></font>
# </div>

# $\q1$ On peut calculer la somme plus ou moins naïvement :

# In[16]:


def eval(A, x):
    s = 0
    for k in range(len(A)):
        s += A[k] * x**k
    return s


# In[17]:


eval([1, -2, 1], 4)


# $\q2$ On recalcule précédemment toutes les puissances de $x$ à partir de zéro. Pour éviter cela, on peut les stocker progressivement :

# In[18]:


def eval(A, x):
    s = A[0]
    p = x
    for k in range(1,len(A)):
        s += A[k] * p
        p *= x
    return s


# In[19]:


eval([1, -2, 1], 4)


# $\q3$ Il y a $d$ itérations de la boucle ```for```, chacune donnant lieu à deux multiplications.
# 
# Le nombre total de multiplications est donc $2d$.

# $\q4$ On implémente la méthode de Hörner de manière itérative :

# In[20]:


def horner(A, x):
    d = len(A)-1
    s = A[d]
    for k in range(1, d + 1):
        s = s * x + A[d - k]
    return s


# In[21]:


horner([1, -2, 1], 4)


# $\q5$ Il y a $d$ itérations dans la boucle ```for```. Chacune donne lieu à une seule multiplication.
# 
# Le nombre total de multiplications est donc exactement $d$.

# <div class="alert alert-block alert-info">
# <center><font size=3><b>Bonus</b></font></center></div>

# In[22]:


N = 1000
Yn = []  # implémentation naive
Yh = []  # methode de Horner
X = list(range(N, 10*N, 200))

x = 3
for n in X:
    A = gen(n, n)
    t1 = toc()
    eval(A, x)
    t2 = toc()
    Yn.append((t2 - t1))
    t1 = toc()
    horner(A, x)
    t2 = toc()
    Yh.append((t2 - t1))

X = X[20:]
Yn = Yn[20:]
Yh = Yh[20:]

plt.plot(X, Yn, 'o', label="méthode naïve")

plt.plot(X, Yh, 'o', label="méthode de Hörner")


plt.legend()
plt.show()