import numpy as np
import matplotlib.pyplot as plt

#Exercício 1:  Movimento de um pêndulo simples

L = 1 #m
theta = 0 #rad
g = 9.8 #m/s^2

def acelera_pendulo(theta, L, g):
    # Aceleração do pêndulo
    a = -g / L * np.sin(theta)
    return a



# a) Simule o movimento do pêndulo usando o método de Euler-Cromer durante 10s.

theta0 = 0.1
v0 = 0

t0 = 0
tf = 10
dt = 0.01
n = int((tf - t0) / dt)
t = np.linspace(t0, tf, n)


def euler_cromer(theta, v, acelera, L, g, dt):
    # Calcula a aceleração
    a = acelera(theta, L, g)
    
    # Atualiza a velocidade
    v += a * dt
    
    # Atualiza o ângulo
    theta += v * dt
    
    return theta, v


# Inicializa os arrays para armazenar os resultados
theta_values = np.zeros(n)
v_values = np.zeros(n)
theta_values[0] = theta0
v_values[0] = v0

for i in range(1, n):
    theta_values[i], v_values[i] = euler_cromer(theta_values[i-1], v_values[i-1], acelera_pendulo, L, g, dt)


# Plota o movimento do pêndulo
plt.figure(figsize=(10, 5))
plt.plot(t, theta_values, label='Ângulo (rad)')
plt.title('Movimento do Pêndulo Simples')
plt.xlabel('Tempo (s)')
plt.ylabel('Ângulo (rad)')
plt.grid()
plt.legend()
plt.show()



#b)Para ângulos pequenos, a equação diferencial pode ser escrito -g/L * theta,
#cuja expressao analítica é theta(t) = A * cos(sqrt(g/L) * t + phi).
#compare esta solução com a solução numérica obtido em alinea a), escolhendo valores
#adequados para 𝐴 e 𝜙.

A = theta0  # Amplitude inicial
phi = 0  # Fase inicial
ω = np.sqrt(g / L)  # Frequência angular
theta_analitica = A * np.cos(ω * t + phi)

# Plota a comparação entre a solução analítica e numérica
plt.figure(figsize=(10, 5))
plt.plot(t, theta_values, label='Solução Numérica (Euler-Cromer)', color='blue')
plt.plot(t, theta_analitica, label='Solução Analítica', color='orange', linestyle='--')
plt.title('Comparação entre Solução Numérica e Analítica')
plt.xlabel('Tempo (s)')
plt.ylabel('Ângulo (rad)')
plt.grid()
plt.legend()
plt.show()


#c) Repita a simulação e a comparação com a teoria com ângulos iniciais de 0.3rad e 0.5rad
angles = [0.3, 0.5]

for angle in angles:
    theta0 = angle
    v0 = 0

    # Recalcula o movimento do pêndulo
    theta_values = np.zeros(n)
    v_values = np.zeros(n)
    theta_values[0] = theta0
    v_values[0] = v0

    for i in range(1, n):
        theta_values[i], v_values[i] = euler_cromer(theta_values[i-1], v_values[i-1], acelera_pendulo, L, g, dt)

    # Recalcula a solução analítica
    theta_analitica = A * np.cos(ω * t + phi)

    # Plota a comparação
    plt.figure(figsize=(10, 5))
    plt.plot(t, theta_values, label='Solução Numérica (Euler-Cromer)', color='blue')
    plt.plot(t, theta_analitica, label='Solução Analítica', color='orange', linestyle='--')
    plt.title(f'Comparação para Ângulo Inicial {angle} rad')
    plt.xlabel('Tempo (s)')
    plt.ylabel('Ângulo (rad)')
    plt.grid()
    plt.legend()
    plt.show()