nota: si vota a favor (o incluso si no lo hace), no olvide desplazarse hacia abajo y ver también la excelente respuesta: ¡es hermosa!
El problema pitagórico de los tres cuerpos, también conocido como problema de Burrau, es un caso especial del problema general de los tres cuerpos, donde los tres cuerpos tienen masas de 3, 4 y 5, y las condiciones iniciales son tales que comienzan en reposo, en el vértices de un triángulo rectángulo 3-4-5.
He pegado algunas capturas de pantalla de los documentos vinculados aquí .
Puedes ver y leer más en este post
Y mire este video: parece que el tiempo que se muestra en la trama del video es tiempo en el papel.
Originalmente, la idea era que podría tener algún significado especial, pero no parece tenerlo. Sin embargo, plantea un gran desafío para los integradores numéricos porque da como resultado varios muy cercanos (~ ) pasa entre pares, y muchos integradores comunes no responderán lo suficientemente rápido con la reducción del tamaño de paso para mantener la precisión numérica.
Esto es lo que me sucedió al usar el integrador ODE predeterminado estándar en SciPy.
Hay algunos trucos para probar dentro de SciPy y, por supuesto, otros integradores disponibles en python, y en realidad puedo implementar algunos métodos de Runge-Kutta de orden superior y escribir mi propio controlador de tamaño de paso hipervigilante . No tiene que ser rápido porque muy pronto, uno de los tres es expulsado y los otros dos se establecen en rotación de dos cuerpos. Esto es bastante común en situaciones de tres cuerpos, en computadoras y en sistemas estelares ternarios que no son lo suficientemente jerárquicos.
Lo que necesito ahora es comparar los resultados con la solución numérica correcta: una tabla con una selección de algunas coordenadas precisas frente al tiempo. Comparar con YouTube no es tan preciso, ¡y tampoco hay garantías de que sean correctos!
¿Alguien sabe dónde puedo encontrar esos números ?
nota: El comentario señala que debo tener cuidado con la palabra "correcto". Estoy buscando resultados usando un solucionador de ODE que funcione bien con ecuaciones rígidas (ver aquí también ) que pueden ser numéricamente inestables , y en este caso se espera que tengan una precisión de, digamos, seis dígitos de precisión por .
Aquí hay una salida de muestra y un script. Esto está mal. Puede encontrar buenas soluciones que se muestran en YouTube y otros lugares, pero no puedo encontrar los resultados numéricos para ayudar a mi depuración.
Si desea sugerir una mejora de Python, puede dejar una respuesta o comentario en mi pregunta en stackoverflow
def deriv(X, t):
Y[:6] = X[6:]
r34, r35, r45 = X[2:4]-X[0:2], X[4:6]-X[0:2], X[4:6]-X[2:4]
thing34 = ((r34**2).sum())**-1.5
thing35 = ((r35**2).sum())**-1.5
thing45 = ((r45**2).sum())**-1.5
Y[6:8] = r34*thing34*m4 + r35*thing35*m5
Y[8:10] = r45*thing45*m5 - r34*thing34*m3
Y[10:12] = -r35*thing35*m3 - r45*thing45*m4
return Y
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint as ODEint
# Pythagorean Three Body Problem
# This script WILL NOT solve it yet, just for illustration of the problem
m3, m4, m5 = 3.0, 4.0, 5.0
x0 = [1.0, 3.0] + [-2.0, -1.0] + [1.0, -1.0]
v0 = [0.0, 0.0] + [ 0.0, 0.0] + [0.0, 0.0]
X0 = np.array(x0 + v0)
t = np.linspace(0, 60, 50001)
Y = np.zeros_like(X0)
tol = 1E-9 # with default method higher precision causes failure
hmax = 1E-04
answer, info = ODEint(deriv, X0, t, rtol=tol, atol=tol,
hmax=hmax, full_output=True)
xy3, xy4, xy5 = answer.T[:6].reshape(3,2,-1)
paths = [xy3, xy4, xy5]
plt.figure()
plt.subplot(2, 1, 1)
for x, y in paths:
plt.plot(x, y)
for x, y in paths:
plt.plot(x[:1], y[:1], 'ok')
plt.xlim(-6, 6)
plt.ylim(-4, 4)
plt.title("This result is WRONG!", fontsize=16)
plt.subplot(4,1,3)
for x, y in paths:
plt.plot(t, x)
plt.ylim(-6, 4)
plt.subplot(4,1,4)
for x, y in paths:
plt.plot(t, y)
plt.ylim(-6, 4)
plt.show()
Lo acabo de ejecutar, y el mío se parece bastante a los del periódico.
Ver algunas coordenadas en la parte inferior.
Aquí hay algunas coordenadas {x, y} en los tiempos en la columna izquierda:
0. {1.,3.} {-2.,-1.} {1.,-1.}
5. {2.46917,-1.22782} {-2.2782,-0.20545} {0.34106,0.901049}
10. {0.77848,0.141392} {-2.02509,0.0972194} {1.15299,-0.162611}
15. {1.41845,0.686214} {-2.00654,0.0599408} {0.754159,-0.459681}
20. {3.00429,0.511925} {-1.38863,-0.470476} {-0.691674,0.0692257}
25. {2.2699,-0.0832} {-2.63692,-0.426417} {0.747596,0.391054}
30. {0.85634,2.28709} {-0.877984,-0.865964} {0.188583,-0.679485}
35. {0.0273748,0.895529} {0.942553,-1.60223} {-0.770468,0.744467}
40. {-0.622004,1.85832} {0.173545,-2.36841} {0.234367,0.779737}
45. {-0.657058,2.53557} {1.61355,-1.23947} {-0.896608,-0.529771}
50. {-2.70146,-3.79723} {1.50595,0.960811} {0.416122,1.50969}
55. {-2.75171,-4.29907} {1.72673,0.97731} {0.269648,1.7976}
60. {0.743681,1.93961} {0.263967,-0.731477} {-0.657382,-0.578586}
65. {4.05348,11.7131} {-1.0722,-3.92197} {-1.57432,-3.8903}
70. {6.93108,20.2566} {-1.99418,-6.87252} {-2.5633,-6.65594}
Eso fue todo con una precisión de trabajo de 30 dígitos. Comprobando la energía total final y el momento angular total contra las condiciones iniciales, con 30 dígitos de trabajo, los resultados son buenos para 10 dígitos. Con 50 dígitos de trabajo, los resultados son buenos para 20 dígitos. Con la precisión de la máquina (alrededor de 15 dígitos de trabajo), los resultados son buenos para cinco o seis dígitos, lo que sigue siendo bastante bueno teniendo en cuenta los acercamientos cercanos.
NDSolve
con InterpolationOrder -> All, WorkingPrecision -> 30, MaxSteps -> 10^5
.SymplecticPartitionedRungeKutta
opción disponible, pero no la usé. Usé los métodos predeterminados, que eligen un predictor-corrector y un método de diferenciación hacia atrás, según la rigidez. Entonces, la energía total final es realmente una buena medida de la calidad del resultado, ya que no hay nada explícito en el método de integración, aparte de las ecuaciones de movimiento, que asegure su conservación.
russell borogove
UH oh
russell borogove
UH oh
tol
puesto grandesUH oh
russell borogove