#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Sep 20 16:02:50 2022

@author: julie
"""

import matplotlib.pyplot as plt
import numpy as np

"Fonction un saut"
def un_saut(x):
    posi = x + np.random.choice([-1,1])
    return posi


#nombre de sauts
K = 50

#Tableau des positions successives et du numéro du saut
pos_part = np.zeros(K+1)
num_saut = np.zeros(K+1)

#Initialisation des sommes pour les calculs de la position moyenne et quadratique sur l'ensemble des sauts
xmoy = 0
r2moy = 0

# Détermination des positions succesives
for i in range(K) :
    pos_part[i+1] = un_saut(pos_part[i])
    num_saut[i+1] = num_saut[i] + 1
    xmoy = xmoy + pos_part[i+1]/K
    r2moy = r2moy + pos_part[i+1]**2/K
    
# Affichage de xmoy et r2moy
print('position moyenne =', xmoy)
print('poition moyenne quadratique =', r2moy)


plt.plot(num_saut,pos_part,color='red',ls='--',marker='+',label="une particule")
plt.xlabel('temps réduit ou numéro du saut')
plt.ylabel('position')
plt.ylim(-K,K)
plt.legend( loc="best")
plt.title('Évolution d une particule')
plt.show()


'Passage à N particules'
N = 1000
k=100

# Tableau des positions des particules à t
positions = np.zeros(N , dtype=int)
dates = np.zeros(k+1 , dtype=int)
# Tableau des valeurs quadratique moyenne des N paricules en focntion du temps
r2moyen = np.zeros(k+1 , dtype=float)

for j in range(k) :
    dates[j+1] = dates[j] + 1
    sommer2 = 0

    for i in range (N-1) :
        positions[i] = un_saut(positions[i])
        sommer2 = sommer2 + positions[i]**2
# Autres rédaction spossibles :
 #   positions[:] = un_saut(positions[:])
 #   r2moyensc[j+1] = np.dot(positions,positions) / N
    r2moyen[j+1] = sommer2/N
#    print(positions)
#print('positions finales :' , positions)
#print(dates)

#print(r2moyen)
#print(r2moyensc)


plt.plot(dates,r2moyen,color='red',ls='--',marker='+')
plt.plot(dates,dates,color='blue')
plt.xlabel('temps réduit')
plt.ylabel('x**2 moyen')
plt.title('Écart quadratique moyen en fonction du temps')
plt.show()
#plt.close('all')


#Dans le plan
posx = np.zeros(N,dtype=float)
posy = np.zeros(N,dtype=float)
tps = np.zeros(k+1,dtype = int)
r2 = np.zeros(k+1 , dtype=float)


for j in range(k) :
    tps[j+1] = tps[j] + 1
    rmoy2 = 0

    for i in range (N-1) :
        angle = np.random.uniform(0,2*np.pi)
        saut = np.random.choice([-1,1])
        posx[i] = posx[i] + saut*np.cos(angle)
        posy[i] = posy[i] + saut*np.sin(angle)
        rmoy2 = rmoy2 + posx[i]**2 + posy[i]**2
#    plt.scatter(posx,posy)
    r2[j+1] = rmoy2/N

    if (j==1) or (4*j)%k == 0  :
        plt.scatter(posx,posy,s=10)

plt.show()
plt.plot(tps,r2,color='red',ls='--',marker='+')
plt.plot(tps,tps,color='blue')
plt.xlabel('temps réduit')
plt.ylabel('r**2 moyen')
plt.show()


