import numpy as np
import math
# creando la clase planeta
class Planeta:
def __init__ (self,posicion,velocidad,masa):
self.posiciones=[[posicion[0]],[posicion[1]],[posicion[2]]]
self.posicion_actual= posicion
self.velocidad = velocidad
self.momento = masa*velocidad
self.masa = masa
self.aceleracion = np.array([0,0,0])
## importando librerias
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
# creando los ejes
fig = plt.figure(figsize=(12,5))
ax = plt.subplot(1,2,1)
ax.set_xlim((-5,5))
ax.set_ylim((-5,5))
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_title('simulacion 2d')
# creando objetos que van a cambia en la animacion . inicialmente estan vacios
# y actualizan sus coordenadas en cada iteracion
punto1, = ax.plot([], [], 'g.', ms=20)
line1, = ax.plot([], [], 'y', lw=2)
punto2,= ax.plot([], [], 'b.', ms=20)
line2, = ax.plot([], [], 'r', lw=2)
ax.legend(['pt1','','pt2','']);
# esconder el plano vacio
plt.close()
## inicialmente la simulacion sera en 2d por lo que la coordenada z sera 0 en velocidad y posicion
planeta1 = Planeta(np.array([4.5,0,0]),np.array([0,2,0]),1)
planeta2 = Planeta(np.array([-4.5,0,0]),np.array([0,-2,0]),1)
G = 100 # Usaremos una constante de gravedad exagerada para las primeras simulaciones de prueba
dt = 0.01 ## segundos
def drawframe(n):
r12 = planeta2.posicion_actual - planeta1.posicion_actual
r12_unitario = r12/(math.sqrt(sum(i**2 for i in r12)))
F12 = ((-1*G*planeta1.masa*planeta2.masa)/(math.sqrt(sum(i**2 for i in r12)))**2)*r12_unitario
F21 = -1*F12 ## la fuerza que el planeta 2 ejerce sobre el planeta 1 sera la misma en magnitud pero en sentido contrario que la furerza que el planeta 1 ejerce en el planta 2
planeta1.momento = planeta1.momento + F21*dt
planeta2.momento = planeta2.momento + F12*dt
planeta1.posicion_actual = planeta1.posicion_actual+ (planeta1.momento*dt)/planeta1.masa
planeta2.posicion_actual = planeta2.posicion_actual+ (planeta2.momento*dt)/planeta2.masa
for k in range(0,2):
planeta1.posiciones[k].append(planeta1.posicion_actual[k])
planeta2.posiciones[k].append(planeta2.posicion_actual[k])
line1.set_data(planeta1.posiciones[0],planeta1.posiciones[1])
punto1.set_data(planeta1.posicion_actual[0] , planeta1.posicion_actual[1])
punto2.set_data(planeta2.posicion_actual[0] , planeta2.posicion_actual[1])
line2.set_data(planeta2.posiciones[0],planeta2.posiciones[1])
return (punto1,line1,punto2,line2)
from matplotlib import animation
# blit=True solo re dibuja las partes que cambian
anim = animation.FuncAnimation(fig, drawframe, frames=1000, interval=20, blit=True)
from IPython.display import HTML
HTML(anim.to_html5_video())
## funcion para la simulación de n cuerpos
## librerias necesarias
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from matplotlib import animation
from IPython.display import HTML
def animation2d(coleccion_planetas,constante_G,delta_t,tiempo_simulado,largo,vectores = False,intervalos= 20):
# crear la figura y los ejes
fig = plt.figure(figsize=(12,5))
ax = plt.subplot(1,2,1)
#tabla = plt.subplot(1,3,3)
#leyenda
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_title('simulacion 2d')
txt_title = ax.set_title('')
# limites de dimensiones
ax.set_xlim((-largo,largo))
ax.set_ylim((-largo,largo))
# Esconde el fondo
plt.close()
#configuracion de puntos y trayectorias
colors = sns.color_palette("husl",len(coleccion_planetas))
puntos = sum([ax.plot([],[],'o',color=colors[k],ms=15,label=f"cuerpo {k+1}") for k in range(0,len(coleccion_planetas))],[])
trayectorias = sum([ax.plot([],[],color=colors[k],lw=2) for k in range(0,len(coleccion_planetas))],[])
if vectores == True:
vectores_vel = []
vectores_acc = []
for i in coleccion_planetas:
tempvec_v = ax.quiver(i.posicion_actual[0],i.posicion_actual[1],i.velocidad[0],i.velocidad[1],color='lightseagreen',scale_units='xy', scale=0.6)
tempvec_a = ax.quiver(i.posicion_actual[0],i.posicion_actual[1],0,0,color='darkred',minlength=3,scale_units='xy', scale=0.5)
vectores_vel.append(tempvec_v)
vectores_acc.append(tempvec_a)
ax.legend(handles=puntos,loc='center left',bbox_to_anchor=(1, 0.5),labelspacing=1,borderpad=2)
#data_tabla = [[np.linalg.norm(i.velocidad),np.linalg.norm(i.aceleracion),f"x:{i.posicion_actual[0]:.2f},y:{i.posicion_actual[1]:.2f}"] for i in coleccion_planetas]
#etiquetas_fila = [f"cuerpo {i+1}" for i in range(0,len(coleccion_planetas))]
#tabla.axis('tight')
#tabla.axis('off')
#tabla.table(cellText=data_tabla,rowLabels=etiquetas_fila ,colLabels=['velocidad','aceleracion','posicion'],colLoc='center',loc='center')
#simulacion de cada iteracion
def draw(n):
for i in coleccion_planetas:
acum_acc = np.array([0,0,0])
for j in coleccion_planetas:
if i != j:
##Calculo del vector que va de i a j
r_ij = j.posicion_actual -i.posicion_actual
r_ij_magnitud = math.sqrt(sum(k**2 for k in r_ij))
r_ij_unitario = r_ij/r_ij_magnitud
##calculando aceleracion con un solo cuerpo
acc_neta_i = ((constante_G*j.masa)/(r_ij_magnitud**2))*r_ij_unitario
#subiendolo al acumulador para sumar la celeracion que produce cada cuerpo
acum_acc =acum_acc + acc_neta_i
#asignando aceleracion calculada
i.aceleracion = acum_acc
##creo que esta bien
#calcular la nueva posicion de cada planeta
for i,index in zip(coleccion_planetas,range(0,len(coleccion_planetas))):
i.velocidad = i.aceleracion*delta_t + i.velocidad
i.posicion_actual = i.posicion_actual + i.velocidad*delta_t
for k in range(0,2):
i.posiciones[k].append(i.posicion_actual[k])
puntos[index].set_data(i.posicion_actual[0],i.posicion_actual[1])
trayectorias[index].set_data(i.posiciones[0],i.posiciones[1])
if vectores == True:
vectores_vel[index].set_UVC(i.velocidad[0],i.velocidad[1])
vectores_vel[index].set_offsets([i.posicion_actual[0],i.posicion_actual[1]])
vectores_acc[index].set_UVC(i.aceleracion[0],i.aceleracion[1])
vectores_acc[index].set_offsets([i.posicion_actual[0],i.posicion_actual[1]])
#datos_actualizados= [[f"{np.linalg.norm(i.velocidad):.2f} m/s",f"{np.linalg.norm(i.aceleracion):.2f}m/s^2",f"x:{i.posicion_actual[0]:.2f},y:{i.posicion_actual[1]:.2f}"] for i in coleccion_planetas]
#tabla.table(cellText=datos_actualizados,rowLabels=etiquetas_fila,colLabels=['velocidad','aceleracion','posicion'],colLoc='center',loc='center')
tiempo = n*delta_t
txt_title.set_text(f"{tiempo:.2f} segundos , Δt = {delta_t} seg.")
return puntos + trayectorias
iteraciones = math.ceil(tiempo_simulado/delta_t)
##los intervalos afectan la duracion de la animacion bajos intervalos son una animacion mas rapida y altos intervalos son una animacion mas lenta
anim = animation.FuncAnimation(fig, draw, frames=iteraciones, interval=intervalos, blit=True)
return HTML(anim.to_html5_video())
estrella1 = Planeta(np.array([4.5,0,0]),np.array([0,2,0]),1)
estrella2 = Planeta(np.array([-4.5,0,0]),np.array([0,-(2+(3/1000)),0]),1)
planeta = Planeta(np.array([25,0,0]),np.array([0,3,0]),1/1000)
cuerpos = [estrella1,estrella2,planeta]
G = 100#constante de gravedad
dt = 0.1#intervalo de tiempo en segundos
t = 90 #segundos de simulacion
animation2d(cuerpos,G,dt,t,40)
sol = Planeta(np.array([0,0,0]),np.array([0,0,0]),1.989e30)
tierra = Planeta(np.array([1.5e11,0,0]),np.array([0,29800,0]),5.97e24)
luna = Planeta(np.array([(1.5e11+3.48e8),0,0]),np.array([0,29790,0]),7.34767e22)
cuerpos = [sol,tierra,luna]
G = 6.67e-11#constante de gravedad
dt = 360#intervalo de tiempo en segundos
t = 3600*24*10 #segundos de simulacion
animation2d(cuerpos,G,dt,t,160000000000)
planeta1 = Planeta(np.array([4.5,0,0]),np.array([0,2,0]),1)
planeta2 = Planeta(np.array([-4.5,0,0]),np.array([0,-2,0]),1)
G = 100
def drawframe(n):
r12 = planeta2.posicion_actual - planeta1.posicion_actual
r12_unitario = r12/(math.sqrt(sum(i**2 for i in r12)))
accNeta1 = ((G*planeta2.masa)/(math.sqrt(sum(i**2 for i in r12)))**2)*r12_unitario
r21 = planeta1.posicion_actual - planeta2.posicion_actual
r21_unitario = r21/(math.sqrt(sum(i**2 for i in r12)))
accNeta2 = ((G*planeta1.masa)/(math.sqrt(sum(i**2 for i in r21)))**2)*r21_unitario
##F21 = -1*F12
planeta1.aceleracion = accNeta1
planeta2.aceleracion = accNeta2
planeta1.velocidad = planeta1.aceleracion*dt + planeta1.velocidad
planeta2.velocidad = planeta2.aceleracion*dt + planeta2.velocidad
planeta1.posicion_actual = planeta1.posicion_actual+ planeta1.velocidad*dt
planeta2.posicion_actual = planeta2.posicion_actual+ planeta2.velocidad*dt
for k in range(0,2):
planeta1.posiciones[k].append(planeta1.posicion_actual[k])
planeta2.posiciones[k].append(planeta2.posicion_actual[k])
#actualizando trayectorias
line1.set_data(planeta1.posiciones[0],planeta1.posiciones[1])
line2.set_data(planeta2.posiciones[0],planeta2.posiciones[1])
punto1.set_data(planeta1.posicion_actual[0] , planeta1.posicion_actual[1])
punto2.set_data(planeta2.posicion_actual[0] , planeta2.posicion_actual[1])
##print(planeta2.posicion_actual)
return (punto1,line1,punto2,line2)
#from matplotlib import animation
# blit=True re-draws only the parts that have changed.
anim = animation.FuncAnimation(fig, drawframe, frames=1000, interval=20, blit=True)
#from IPython.display import HTML
HTML(anim.to_html5_video())
planeta_1 = Planeta(np.array([-10,0,0]),np.array([0,-5,0]),10)
planeta_2 = Planeta(np.array([10,0,0]),np.array([0,5,0]),10)
estrella = Planeta(np.array([0,0,0]),np.array([0,0,0]),1)
sistemad = [planeta_1,planeta_2,estrella]
G = 100#constante de gravedad
dt = 0.05#intervalo de tiempo en segundos
t = 90 #segundos de simulacion
animation2d(sistemad,G,dt,t,15,vectores=True)
planeta_1 = Planeta(np.array([-10,0,0]),np.array([0,-6,0]),10)
planeta_2 = Planeta(np.array([10,0,0]),np.array([0,6,0]),10)
estrella = Planeta(np.array([0,0,0]),np.array([0,0,0]),1)
sistemad = [planeta_1,planeta_2,estrella]
G = 100#constante de gravedad
dt = 0.05#intervalo de tiempo en segundos
t = 90 #segundos de simulacion
animation2d(sistemad,G,dt,t,15)
v = 1.5
cuerpo1 = Planeta(np.array([10,0,0]),np.array([v*math.cos(1.0472),v*math.sin(1.0472),0]),1)
cuerpo2 = Planeta(np.array([-10,0,0]),np.array([v*math.cos(1.0472),-v*math.sin(1.0472),0]),1)
cuerpo3 = Planeta(np.array([0,17.3205,0]),np.array([-v,0,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 100#constante de gravedad
dt = 0.01#intervalo de tiempo en segundos
t = 90 #segundos de simulacion
animation2d(cuerpos,G,dt,t,20)
cuerpo1 = Planeta(np.array([-0.9892620043,0,0]),np.array([0,1.9169244185,0]),1)
cuerpo2 = Planeta(np.array([2.2096177241,0,0]),np.array([0,0.1910268738,0]),1)
cuerpo3 = Planeta(np.array([-1.2203557197,0,0]),np.array([0,-2.1079512924,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.005#intervalo de tiempo en segundos
t = 10 #segundos de simulacion
animation2d(cuerpos,G,dt,t,3)
cuerpo1 = Planeta(np.array([0.97000436,-0.24308753,0]),np.array([0.93240737/2,0.86473146/2,0]),1)
cuerpo2 = Planeta(np.array([-0.97000436,0.24308753,0]),np.array([0.93240737/2,0.86473146/2,0]),1)
cuerpo3 = Planeta(np.array([0,0,0]),np.array([-0.93240737,-0.86473146,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.01#intervalo de tiempo en segundos
t = 5 #segundos de simulacion
animation2d(cuerpos,G,dt,t,1.5)
cuerpo1 = Planeta(np.array([-1,0,0]),np.array([0.392955,0.097579,0]),1)
cuerpo2 = Planeta(np.array([1,0,0]),np.array([0.392955,0.392955,0]),1)
cuerpo3 = Planeta(np.array([0,0,0]),np.array([-0.78591,-0.195158,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.0001#intervalo de tiempo en segundos
t = 10 #segundos de simulacion
animation2d(cuerpos,G,dt,t,5)
cuerpo1 = Planeta(np.array([-1,0,0]),np.array([0.080584,0.588836,0]),1)
cuerpo2 = Planeta(np.array([1,0,0]),np.array([0.080584,0.588836,0]),1)
cuerpo3 = Planeta(np.array([0,0,0]),np.array([-0.161168,-1.177672,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.0001#intervalo de tiempo en segundos
t = 5 #segundos de simulacion
animation2d(cuerpos,G,dt,t,3,intervalos=3)
cuerpo1 = Planeta(np.array([-1,0,0]),np.array([0.464445,0.39606,0]),1)
cuerpo2 = Planeta(np.array([1,0,0]),np.array([0.464445,0.39606,0]),1)
cuerpo3 = Planeta(np.array([0,0,0]),np.array([-0.92889,-0.79212,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.001#intervalo de tiempo en segundos
t = 30 #segundos de simulacion
animation2d(cuerpos,G,dt,t,3,intervalos=3)
cuerpo1 = Planeta(np.array([0.419698802831 ,1.190466261252,0]),np.array([0.1022945660031,0.687248445943,0]),1)
cuerpo2 = Planeta(np.array([0.076399621771, 0.296331688995,0]),np.array([0.148950262064,0.240179781043,0]),1)
cuerpo3 = Planeta(np.array([0.100310663856,-0.729358656127,0]),np.array([-0.251244828060,-0.9274282269779,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.0001#intervalo de tiempo en segundos
t = 6 #segundos de simulacion
animation2d(cuerpos,G,dt,t,2)
cuerpo1 = Planeta(np.array([-1 ,0,0]),np.array([0.2869236336,0.0791847624,0]),1)
cuerpo2 = Planeta(np.array([1, 0,0]),np.array([0.2869236336,0.0791847624,0]),1)
cuerpo3 = Planeta(np.array([0,0,0]),np.array([-1.1476945344,-0.3167390496,0]),1/2)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 1#constante de gravedad
dt = 0.001#intervalo de tiempo en segundos
t = 6 #segundos de simulacion
animation2d(cuerpos,G,dt,t,2,intervalos=10)
cuerpo1 = Planeta(np.array([-2 ,0,0]),np.array([0,0,0]),1)
cuerpo2 = Planeta(np.array([2, 0,0]),np.array([0,0,0]),1)
cuerpo3 = Planeta(np.array([0,3.46,0]),np.array([0,0,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 0.5#constante de gravedad
dt = 0.01#intervalo de tiempo en segundos
t = 30 #segundos de simulacion
animation2d(cuerpos,G,dt,t,10,intervalos=3)
estrella = Planeta(np.array([0,0,0]),np.array([0,0,0]),200)
planeta1 = Planeta(np.array([60,0,0]),np.array([0,20,0]),1/1000)
planeta2 = Planeta(np.array([100,0,0]),np.array([0,15,0]),1/1000)
cuerpos = [estrella,planeta1,planeta2]
G = 100#constante de gravedad
dt = 0.1#intervalo de tiempo en segundos
t = 100 #segundos de simulacion
animation2d(cuerpos,G,dt,t,150)
estrella1 = Planeta(np.array([15,0,0]),np.array([0,10,0]),100)
estrella2 = Planeta(np.array([-15,0,0]),np.array([0,-10,0]),100)
planeta1 = Planeta(np.array([60,0,0]),np.array([0,20,0]),1/1000)
planeta2 = Planeta(np.array([100,0,0]),np.array([0,15,0]),1/1000)
cuerpos = [estrella1,estrella2,planeta1,planeta2]
G = 100#constante de gravedad
dt = 0.01#intervalo de tiempo en segundos
t = 100 #segundos de simulacion
animation2d(cuerpos,G,dt,t,150)
# sistema 1
vel = 1.5
cuerpo1 = Planeta(np.array([10,0,0]),np.array([0,vel,0]),1)
cuerpo2 = Planeta(np.array([20*math.cos(2.888903884)-10,20*math.sin(2.888903884),0]),np.array([2*vel,0,0]),1)
cuerpo3 = Planeta(np.array([20*math.cos(2.888903884)-10,-20*math.sin(2.888903884),0]),np.array([-2*vel,0,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 100#constante de gravedad
dt = 0.1#intervalo de tiempo en segundos
t = 150 #segundos de simulacion
animation2d(cuerpos,G,dt,t,100,vectores=True)
# sistema 2
vel = 1
cuerpo1 = Planeta(np.array([10.5,0,0]),np.array([0,vel,0]),1)
cuerpo2 = Planeta(np.array([20*math.cos(2.888903884)-10,20*math.sin(2.888903884),0]),np.array([2*vel,0,0]),1)
cuerpo3 = Planeta(np.array([20*math.cos(2.888903884)-10,-20*math.sin(2.888903884),0]),np.array([-2*vel,0,0]),1)
cuerpos = [cuerpo1,cuerpo2,cuerpo3]
G = 100#constante de gravedad
dt = 0.1#intervalo de tiempo en segundos
t = 150 #segundos de simulacion
animation2d(cuerpos,G,dt,t,100,vectores=True)
# aqui hubo un problema por el intervalo de tiempo
estrella1 = Planeta(np.array([10,0,0]),np.array([0,10,0]),100)
estrella2 = Planeta(np.array([-10,0,0]),np.array([0,-10,0]),100)
planeta1 = Planeta(np.array([25,0,0]),np.array([0,3,0]),1/1000)
cuerpos = [estrella1,estrella2,]
G = 100#constante de gravedad
dt = 0.1#intervalo de tiempo en segundos
t = 90 #segundos de simulacion
animation2d(cuerpos,G,dt,t,40)