import numpy as np
import matplotlib.pyplot as plt

# défnir les parmètres du circuit
f0 = 480
w0 = 2*np.pi*f0
T0 = 1/f0
sig = 0.1
A=3+sig
Usat= 15
Umax= Usat/A

tmax = 25*T0 # pour afficher environ 25 pseudo-oscillations

# discrétiser le temps et définir le pas 
N = 10000
t=np.linspace(0,tmax,N+1)
dt = tmax/N  # pas de calcul

# créer le tableau des tensions solution de l'éq diff et définir les CI
u=np.empty(N+1)
u[0]=0.05  # tension initiale de 50 mV
u[1]=u[0] # dérivée de la tension initiale égale à 0

# calcul des valeurs de u(t) par schéma d'Euler explicite
for k in range (1,N):
    if (abs(u[k])<= Umax):
        u[k+1]=1/(1-w0*sig*dt)*(u[k]*(2-sig*w0*dt-w0**2*dt**2)-u[k-1])
    else:
        u[k+1]=1/(1+3*w0*dt)*(u[k]*(2+3*w0*dt-w0**2*dt**2)-u[k-1])
    
plt.plot(t,u)
plt.xlabel('Temps (s)')
plt.ylabel('Ue (V)')
plt.title('Croissance des oscillations de Ue(t)')
plt.show()