El salto de Felix Baumgartner en Python

Introducción

Supongo que todos estáis al tanto de la hazaña de Felix Baumgartner, el hombre que ha saltado desde una altura de más de 120000 pies desde un globo, convirtiéndose en el hombre que más alto ha saltado y el que ha alcanzado la mayor velocidad sin ayuda mecánica como parte de la misión Red Bull Stratos.

En Pybonacci somos tan frikis, que vamos a visualizar el salto supersónico de Baumgartner como mejor se nos da: con Python 😉

Nota: Esto es un artículo recreativo que he escrito en un par de horas y he hecho unas cuantas suposiciones que no tienen porqué coincidir con la realidad. Tómese esto en cuenta a la hora de valorar los resultados.

Formulación del problema

Nota: Esta es la parte aburrida. Las gráficas están más abajo 🙂
La ecuación diferencial que gobierna el movimiento en caída libre de Baumgartner es

\( \displaystyle m \frac{d^2 y}{d t^2} = -m g + D \)

donde el término \( -m g \) corresponde a la atracción gravitatoria y \( D \) es el rozamiento aerodinámico. La presencia de esta componente es fundamental, porque es la que determina la velocidad de equilibrio o velocidad terminal que se alcanza durante el salto. En los típicos problemas de escuela secundaria se desprecia, pero en Pybonacci hacemos las cosas con rigor 😉

El rozamiento aerodinámico tiene por expresión \( D = \frac{1}{2} \rho v^2 C_D A \), donde \( C_D \) y \( A \) son, respectivamente, el coeficiente de rozamiento y un área de referencia y son parámetros que dependen del cuerpo que estemos estudiando. Vemos también que el rozamiento es proporcial al cuadrado de la velocidad, de tal forma que cuanto más rápido vaya el cuerpo, mayor será esta fuerza: por eso se alcanza una velocidad terminal.

Hay una cosa más que hay que tener en cuenta, y es que al ser el salto desde una altura tan grande, la densidad del aire no puede considerarse constante.

¿Qué significa esto? Que, una vez alcanzada la velocidad terminal, el rozamiento irá aumentando conforme pase el tiempo, aquella irá disminuyendo… y el cuerpo caerá cada vez más despacio.

«Ehm… ¿cada vez más despacio? Chicos de Pybonacci, dedicaos a Python porque la física se os da fatal». Eso pensáis, ¿eh? ¡ya veremos!

Para completar la formulación del problema nos falta saber el valor del coeficiente de rozamiento \( C_D \), el área de referencia \( A \) y la ley de variación de la densidad \( \rho \) con la altitud. Para los dos primeros, con el permiso de Arturo Quirantes Sierra, vamos a suponer que \( C_D = 0.4 \) y \( A = 1.0~m^2 \).

En cuanto a la densidad, vamos a acudir utilizar el paquete AeroCalc, que nos permite, entre otras cosas, calcular propiedades de la atmósfera según el modelo de la Atmósfera Estándar Internacional de 1976 hasta 84.852 kilómetros.

¡Ya podemos resolver el problema! Para ello, vamos a hacer algo que no habíamos explicado todavía en el blog que es cómo resolver ecuaciones diferenciales ordinarias o EDOs en Python con SciPy, en concreto problemas de Cauchy o de valor inicial, como es este caso. Para ello utilizaremos la función odeint del paquete scipy.integrate. Si tenemos un sistema de EDOs del tipo

\( \displaystyle \frac{d y}{d t} = f(y, t) \)

la función odeint acepta, como mínimo, estos argumentos:

  • La función del sistema f(y, t0, …), que a su vez recibe como argumentos el vector y y el instante de tiempot0.
  • El array y0 de condiciones iniciales en el instante t[0].
  • El array t de valores temporales para los que resolver el sistema de ecuaciones diferenciales. El primer valor debe ser el que corresponde al vector de condiciones iniciales y0.

Como vemos la forma del sistema dada más arriba no corresponde con la forma en la que tenemos nosotros expresada la ecuación: al ser una ecuación de segundo orden, hay que transformarla en un sistema de dos ecuaciones de primer orden. Para eso, definimos la variable

\( v = \frac{d y}{d t} \)

y el sistema queda de la siguiente manera:

\( \displaystyle \frac{d}{d t} \begin{pmatrix} y \\ v \end{pmatrix} = \begin{pmatrix} v \\ -g – \frac{1}{2 m} \rho(y) v^2 C_D A \end{pmatrix} \)

¡Ya podemos empezar a calcular!

Caída libre

Este es el código del programa. Se explica por sí solo:

Altitud y velocidad de caída en función del tiempo

Ajá, ¿no os creíais que cada vez caía más despacio? 😛 Esto evidentemente tenía cierto truco, pero ¡ya se ven los resultados en las gráficas!

La barrera del sonido

Y ahora la gran pregunta, ¿superó Baumgartner la barrera del sonido? Para saberlo, necesitamos saber, aparte de su velocidad real, la distribución de la velocidad del sonido en el aire, que no es constante al no serlo la temperatura. Con esto podemos calcular el número de Mach:

\( \displaystyle M = \frac{u}{c} \)

siendo \( u \) la velocidad, \( c = \sqrt{\gamma R T} \) la velocidad del sonido en el aire, \( \gamma \) y \( R \) dos parámetros de valores conocidos (la razón de calores específicos y la constante específica del aire) y \( T \) la temperatura, que conocemos gracias al modelo de la Atmósfera Estándar Internacional. Veamos:

Número de Mach en función del tiempo

Así que sí, ¡Felix Baumgartner superó la barrera del sonido! Según Pybonacci, por supuesto 😉

19 pensamientos sobre “El salto de Felix Baumgartner en Python”

  1. jusjusjus… Me has ‘quitado’ el artículo de los dedos. Muy ilustrativo y bien explicado.
    Ahora me toca exprimirme el cerebro para idear algo nuevo :-0

  2. Excelente, me ha gustado mucho el análisis!
    Por cierto, ¿para qué sirve la línea ax1.plot(t, np.zeros_like(t), 'k--')? No he ejecutado el código aún, así que si lo hubiera hecho igual lo veía evidente…
    Un saludo

  3. Más allá de la hazaña de este hombre! que me parece mucha perdida de tiempo y plata (cada uno tiene sus ideales ¿no?), me pareció un artículo muy interesante.
    Nunca resolví ODE con python y me viene bien para cuando lo necesite.
    La herramienta que usaste para conocer la densidad del aire (o atmósfera) esta interesante, también estaría interesante un gráfico de densidad vs altura para tener entender mejor lo que vos insistís en el artículo.
    ¡Buen artículo!
    ha… algo que quería comentar era sobre la herramienta que usas para renderizar las formulas (ecuaciones latex) ¿cual es? ¿han evaluado cambiarla por mathjax?

    1. Gracias por el comentario! A mí no me parece una pérdida ni de tiempo ni de plata, pero creo que no es el lugar para discutir eso. Gracias por el comentario! 🙂
      Y oye, lo de la gráfica de la densidad frente a la altura es una buena idea. Mañana lo haré un poco más despejado. Un saludo!

  4. Muy buen artículo. Todavía no llego creerme del todo que se pueda “romper” la barrera del sonido en un medio que no parece suficientemente denso como para transmitir el sonido (al menos todo el rango de frecuencias audibles), ni mucho menos para transmitir el estruendo caracterítico que dan por hecho en los telediarios. Pero la matemáticas son indiscutibles. Gran trabajo.

    1. Bueno, creo que entiendo un poco a qué te refieres, pero si quieres puedes olvidarte del sonido como percepción auditiva por un momento y considerarlo simplemente como ondas transversales de presión. Y sí, las matemáticas son indiscutibles 😛 Gracias!

  5. Pingback: Bitacoras.com

  6. Pingback: Ecuaciones de Lotka-Volterra: modelo presa depredador | Pybonacci

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *

+ 21 = twenty seven