Fórmula triginométrica para calcular los resultados de una colisión elástica en 2D

Question

Estoy intentando escribir un código que gestione las colisiones elásticas entre dos círculos en el espacio 2D pero no consigo localizar una fórmula que funcione.

Para la colisión en sí las variables conocidas son:

• masa de ambas partículas
• velocidades iniciales de ambas partículas
• ángulos de movimiento iniciales para ambas partículas
• ángulo de contacto

Lo que intento establecer son las velocidades finales y los ángulos de movimiento de las partículas que colisionan.

Según mis investigaciones, casi todas las fuentes que intentan resolver este problema hacen referencia a la siguiente fórmula disponible en Wikipedia ( enlace )

\begin{align} v'_{1x}&=\frac{v_{1}\cos(\theta_1-\varphi)(m_1-m_2)+2m_2v_{2}\cos(\theta_2-\varphi)}{m_1+m_2}\cos(\varphi) \\[0.2em] &\quad+v_{1}\sin(\theta_1-\varphi)\cos(\varphi+\frac{\pi}{2}) \\[0.8em] v'_{1y}&=\frac{v_{1}\cos(\theta_1-\varphi)(m_1-m_2)+2m_2v_{2}\cos(\theta_2-\varphi)}{m_1+m_2}\sin(\varphi) \\[0.2em] &\quad+v_{1}\sin(\theta_1-\varphi)\sin(\varphi+\frac{\pi}{2}) \end{align}

La fórmula en sí está tomada de este sitio web . Y al parecer ya ha causado problemas a los usuarios de este intercambio de pilas anteriormente ( un enlace a una pregunta anterior utilizando esta fórmula )

El problema que tengo es que al replicar esta fórmula en mi código no siempre proporciona lo que parece ser un resultado correcto, pero funciona como se espera en algunos de los casos.

Estaba verificando los resultados esperados con la simulación disponible aquí . Los resultados producidos por esta simulación difieren de los producidos por la fórmula, a pesar de que la misma fórmula también se muestra en el sitio web que alberga la simulación.

Mi pregunta es: ¿es realmente correcta esta fórmula y, si no lo es, cómo sería la fórmula correcta? Además, sé que es posible resolver esto con vectores, pero estoy buscando una solución trigonométrica.

Gracias

---

UPD: algunos ejemplos de problemas

Ejemplo 1 (parece funcionar como se espera):

m1 = 100; m2 = 10; v1 = 6; v2 = 1; ángulo1 = 0; ángulo2 = 0; ángulo_de_contacto = 0

resultados dados:

• v1x = 5,090909
• v1y = 0
• v2x = 10,09091
• v2y = 0

Ejemplo 2 (parece funcionar como se espera):

m1 = 30; m2 = 20; v1 = 40; v2 = 50; ángulo1 = 315; ángulo2 = 135; ángulo_de_contacto = -45

resultados dados:

• v1x = -22,62742
• v1y = 22,62742
• v2x = 41,01219
• v2y = -41,01219

Ejemplo 3 (parece que está mal):

m1 = 10; m2 = 10; v1 = 6; v2 = 1; angle1 = 0; angle2 = 180; contact_angle = -45

resultados dados:

• v1x = 2,581109
• v1y = -2,416636
• v2x = 2,486482
• v2y = -2,513894

¿Seguro que no pueden caer los dos después de la colisión?

---

UPD2: El ejemplo 3 fue producido por un código incorrecto. Mi respuesta a continuación incluye una versión que produce resultados correctos

Answer 1

La fórmula parece ser correcta o al menos produce los mismos resultados que la simulación en "scincalculators" ( enlace ) que estaba usando como referencia.

El problema en mi ejemplo 3 se debe a un error de codificación.

Aquí hay un código Python que produce los resultados correctos ( enlace al fragmento de código ), espero que a alguien le sea útil:

import math

mass1 = 10
mass2 = 10
velocity1 = 6
velocity2 = 1
angle1 = 0
angle2 = 180
contact_angle = -45
pi = math.pi

def final_vx(m1, m2, v1, v2, O1, O2, phi):
  numerator_p1 = v1 * math.cos((O1 - phi) * pi / 180) * (m1 - m2)
  numerator_p2 = 2 * m2 * v2 * math.cos((O2 - phi) * pi / 180)
  numerator = numerator_p1 + numerator_p2
  denominator = m1 + m2
  addition = v1 * math.sin((O1 - phi) * pi / 180) * math.cos((phi * pi / 180) + pi / 2)
  result = (numerator / denominator) * math.cos(phi * pi / 180) + addition

  return(result)

def final_vy(m1, m2, v1, v2, O1, O2, phi):
  numerator_p1 = v1 * math.cos((O1 - phi) * pi / 180) * (m1 - m2)
  numerator_p2 = 2 * m2 * v2 * math.cos((O2 - phi) * pi / 180)
  numerator = numerator_p1 + numerator_p2
  denominator = m1 + m2
  addition = v1 * math.sin((O1 - phi) * pi / 180) * math.sin((phi * pi / 180) + pi / 2)
  result = (numerator / denominator) * math.sin(phi * pi / 180) + addition

  return(result)

v1x = final_vx(mass1, mass2, velocity1, velocity2, angle1, angle2, contact_angle)
v1y = final_vy(mass1, mass2, velocity1, velocity2, angle1, angle2, contact_angle)
v2x = final_vx(mass2, mass1, velocity2, velocity1, angle2, angle1, contact_angle)
v2y = final_vy(mass2, mass1, velocity2, velocity1, angle2, angle1, contact_angle)
print("v1x = " + str(v1x))
print("v1y = " + str(v1y))
print("v2x = " + str(v2x))
print("v2y = " + str(v2y))