import matplotlib.pyplot as plt
import numpy as np

# Paramètres en USI =====================
g = 9.81     # champ de pesanteur
m = 5        # masse du bloc
k = 20       # constante de raideur du ressort
L0 = 0.15    # longueur à vide du ressort

fs = 0.4     # coefficient de frottement statique bloc sur support
fd = 0.2     # coefficient de frottement dynamique bloc sur support

u = 0.5      # vitesse de l'extrémité A du ressort

x0 = -L0        # position initiale de l'extrémité B à compléter
v0 = 0        # vitesse initiale de l'extrémité B à compléter

# Nombre d’iterations, instant final et pas de temps en USI =====================
n = 10001
tf = 30  
dt = tf/(n-1)

# Listes =====================
# Listes des instants 
t = [i*dt for i in range(n)]
# Création d’une liste pour les positions et les vitesses
xB = [0 for i in range(n)]
vB = [0 for i in range(n)]
# Liste des modes glissement/non glissement ; initialement, on affecte tout à la chaîne de caractères '...'
mode=['...' for i in range(n)]

# Initialisation des listes =====================
xB[0] = x0
vB[0] = v0
mode[0]='non glissement'

# Remplissage des listes par récurrence =====================
for i in range(n-1):
    # pour remplir les listes positions et vitesses, il faut distinguer les cas non-glissement / glissement
    # En mode de non glissement :
    if mode[i]=='non glissement':
        # Remplissage des listes des positions et des vitesses dans le cas non glissement -> lignes 43-44 à compléter
        xB[i+1]=xB[i]
        vB[i+1]=0
        # on teste la condition de non-glissement pour remplir la liste des modes -> ligne 46 à compléter
        if abs(k*(xB[i+1]-u*t[i+1]+L0))>fs*m*g:
            mode[i+1]='glissement'
        else:
            mode[i+1]='non glissement'
    # En mode de  glissement selon +ux :
    else:
        # Remplissage des listes des positions et des vitesses dans le cas glissement à l'aide de la méthode d'Euler -> lignes 53-54 à compléter
        xB[i+1]=xB[i]+vB[i]*dt
        vB[i+1]=vB[i]+(-fd*g+k/m*(u*t[i]-L0-xB[i]))*dt
        # on teste la condition de glissement selon +ux pour remplir la liste des modes -> ligne 56 à compléter
        if vB[i+1]<=0:
            mode[i+1]='non glissement'
        else:
            mode[i+1]='glissement'

# Tracé de la courbe xB en fonction du temps -> ligne 62 à compléter =====================
plt.plot(t,xB)   
plt.xlabel("t (s)") 
plt.ylabel("xB(t) (m)") 
plt.grid()
plt.show()