import numpy as np # llama a la libreria numpy (operaciones numericas) import matplotlib.pyplot as plt # llama a la libreria grafica from scipy.integrate import odeint # llama al integrador numerico 'odeint' # vamos a definir una funcion maricial que realiza la siguente cuenta # # [0 1].[ y(t) ] # [-w^2 0] [dy(t)/dt] # def sistema(Y,t,omega): y,ydot=Y return ydot,-omega**2*y t=np.linspace(0,2*np.pi,101) # puntos donde se evalua la funcion y0=[1,0] # posicion inicial in 1 y velocidad nula omega=2.0 # frecuencia angular del pendulo w^2=g/l f=odeint(sistema,y0,t,args=(omega,)) # integrador (de Euler) y,ydot=f[:,0],f[:,1] plt.figure() # inicio nuevo grafico #plt.plot(y,ydot) # grafico de diagrama de fases plt.plot(t,y,t,ydot) # grafico y vs. t e ydot vs. t plt.title("Solucion para omega= %5g" % omega) # titulo plt.show() # mostrar en pantalla