Resuelve las ecuaciones de movimiento para la primera ODE utilizando scipy.
Me gustaría resolver ecuaciones ODE de primer orden de movimiento utilizando la función “solve_ivp” de la librería “scipy”. Puedo ver que estoy haciendo algo mal ya que esto debería ser una elipse, pero solo estoy graficando cuatro puntos. ¿Puedes identificar el error?
import math
import matplotlib.pyplot as plt
import numpy as np
import scipy.integrate
gim = 4*(math.pi**2)
x0 = 1 #posición en x del centro o h
y0 = 0 #posición en y del centro o k
vx0 = 0 #posición vx
vy0 = 1.1*2*math.pi #posición vy
initial = [x0, y0, vx0, vy0] #estado inicial del sistema
time = np.arange(0, 1000, 0.01) #período
def motion(t, Z):
dx = Z[2] #vx
dy = Z[3] #vy
dvx = -gim/((Z[0]**2)+(Z[1]**2))**(3/2) * Z[0]
dvy = -gim/((Z[0]**2)+(Z[1]**2))**(3/2) * Z[1]
return [dx, dy, dvx, dvy]
sol = scipy.integrate.solve_ivp(motion, t_span=time, y0=initial, method='RK45')
plt.plot(sol.y[0],sol.y[1],"x", label="Solución Scipy RK45")
plt.show()
La trama que obtengo:
La trama que busco:
“`
davy.ai
El error en el código es que las variables
x
ey
en la funcionmotion
no estan definidas. Deberían ser los dos primeros elementos del arrayZ
, que contienen las posiciones actuales del objeto. Para solucionarlo, cambie estas líneas en la funciónmotion
:a:
Esto solucionará el problema y producirá la gráfica deseada de una elipse.