Cuando una microburbuja arde como una estrella bajo el agua — Keller–Miksis y la fuerza secundaria de Bjerknes#

¿Qué le ocurre a una burbuja de aire de 5 micrómetros de diámetro cuando la alcanza el ultrasonido? La burbuja se hincha y, al instante siguiente, queda aplastada a la décima parte de su tamaño. En ese breve colapso el gas interior se calienta hasta miles de grados. De ahí surge la sonoluminiscencia (emisión de luz de la burbuja), ese destello parecido a la luz de una estrella que salta dentro de un vaso de agua. El mismo colapso es el que limpia unas gafas en un baño ultrasónico y el que rompe cálculos con ultrasonido.

La ecuación que gobierna esta violenta oscilación radial es la ecuación de Keller–Miksis. Este texto sigue el artículo de Nagy y Hegedűs (2025): plantea la dinámica radial de una sola burbuja y la integra directamente en Python. Después reproduce incluso la fuerza secundaria de Bjerknes, con la que dos burbujas se empujan y se atraen entre sí.

“

Artículo: D. Nagy, F. Hegedűs, "Assessing the accuracy of the coupled-spherical-bubble approach for bubble pairs in an acoustic field", Ultrasonics Sonochemistry 123 (2025) 107651. DOI: 10.1016/j.ultsonch.2025.107651

La ecuación que gobierna una sola burbuja#

Si se considera la burbuja como una esfera perfecta, la única incógnita que queda es el radio $R(t)$. Tomando el líquido circundante como incompresible e integrando la cantidad de movimiento, se obtiene la ecuación de Rayleigh–Plesset (ecuación de movimiento radial de una burbuja esférica).

$\rho_L$ es la densidad del líquido, $R$ es el radio y la notación de punto indica derivada temporal. El lado izquierdo es la inercia del líquido que se desplaza radialmente. El lado derecho es la diferencia de presión que empuja la pared de la burbuja hacia adentro y hacia afuera.

La presión del líquido en la pared $p_L$ se separa en tres términos.

$p_G$ es la presión del gas dentro de la burbuja, $\sigma$ es la tensión superficial y $\mu_L$ es la viscosidad del líquido. El segundo término es la tensión superficial que contrae la burbuja; el tercero es la resistencia viscosa que se opone al movimiento de la pared.

Si se considera que el gas se comprime de forma casi adiabática, la relación politrópica cierra el sistema.

$R_E$ es el radio de equilibrio, $p_0$ es la presión ambiente y $n_G$ es el exponente politrópico del gas. Cuando el radio se reduce a la mitad, la presión del gas se dispara a $2^{3n_G}\approx 11$ veces. Este brusco rebote es la fuerza que expulsa el colapso.

El forzamiento se aplica a través de la presión lejana.

$p_A$ es la amplitud del ultrasonido y $f$ es la frecuencia. En el semiciclo en que la presión baja, la burbuja se hincha; en el semiciclo en que sube, queda aplastada.

El líquido también se comprime — la corrección de Keller–Miksis#

El punto débil de Rayleigh–Plesset es la suposición de que "el líquido es incompresible". En el instante del colapso, cuando la velocidad de la pared se acerca a la velocidad del sonido, esa suposición se rompe. La ecuación de Keller–Miksis incorpora hasta primer orden el hecho de que la velocidad del sonido en la pared $c_L$ es finita.

$c_L$ es la velocidad del sonido en el líquido. Cuando $\dot{R}/c_L$ es pequeño, todos los paréntesis convergen a 1 y se regresa a Rayleigh–Plesset. A medida que crece la velocidad de la pared, estos términos disipan la energía del colapso hacia ondas acústicas y atenúan la amplitud. Nagy y Hegedűs clasifican la precisión según el número de Mach de la pared. Rayleigh–Plesset es de orden 0, Keller–Miksis de orden 1 y el modelo de Gilmore de orden 2.

En la simulación de abajo es posible manipular todo directamente. Se alternan los dos modelos aplicando el mismo forzamiento.

Si la amplitud de forzamiento se eleva hasta 1,2 atmósferas, la burbuja se hincha mucho y luego colapsa de forma abrupta. Rayleigh–Plesset representa un rebote excesivo. Al cambiar a Keller–Miksis, el rebote tras el colapso se atenúa de manera notable. Esa es la energía que se escapó en forma de ondas acústicas.

Integración directa del colapso#

Se reduce la ecuación a un sistema de primer orden en el estado $[R,\dot{R}]$ y se integra con Runge–Kutta de cuarto orden. Como el colapso es muy agudo, el paso de tiempo requiere del orden de nanosegundos.

import numpy as np
 
P0, RHO, SIGMA, MU, C, KAPPA = 1.0e5, 998.0, 0.0725, 1.0e-3, 1481.0, 1.4
 
def bubble_rhs(y, t, RE, pA, f):
    R, V = y
    w = 2 * np.pi * f
    pG0 = P0 + 2 * SIGMA / RE
    pG = pG0 * (RE / R) ** (3 * KAPPA)
    pL = pG - 2 * SIGMA / R - 4 * MU * V / R
    pInf = P0 - pA * np.sin(w * t)
    # d/dt(pL - pInf), con el término R-ddot aislado del lado izquierdo
    dP = (-3 * KAPPA * pG * V / R + 2 * SIGMA * V / R**2
          + 4 * MU * V**2 / R**2 + pA * w * np.cos(w * t))
    num = (-1.5 * (1 - V / (3 * C)) * V**2
           + (1 + V / C) * (pL - pInf) / RHO
           + R * dP / (RHO * C))
    den = (1 - V / C) * R + 4 * MU / (RHO * C)
    return np.array([V, num / den])
 
def integrate_km(RE=5e-6, pA=1.2e5, f=60e3, dt=5e-10, steps=200000):
    y = np.array([RE, 0.0])
    t, Rlog = 0.0, []
    for _ in range(steps):
        k1 = bubble_rhs(y, t, RE, pA, f)
        k2 = bubble_rhs(y + 0.5 * dt * k1, t + 0.5 * dt, RE, pA, f)
        k3 = bubble_rhs(y + 0.5 * dt * k2, t + 0.5 * dt, RE, pA, f)
        k4 = bubble_rhs(y + dt * k3, t + dt, RE, pA, f)
        y = y + dt / 6 * (k1 + 2 * k2 + 2 * k3 + k4)
        y[0] = max(y[0], 0.1 * RE)   # protege contra el colapso singular
        t += dt
        Rlog.append(y[0])
    return np.array(Rlog)
 
R = integrate_km()
print(f"R_max / R_E = {R.max() / 5e-6:.2f}")
print(f"R_min / R_E = {R.min() / 5e-6:.2f}")
# R_max / R_E = 3.71
# R_min / R_E = 0.10

Se hincha hasta 3,7 veces el radio de equilibrio y luego queda aplastada a la décima parte. En el radio mínimo, la presión del gas alcanza miles de atmósferas y la temperatura miles de grados. Son las condiciones en las que surge la sonoluminiscencia.

Las burbujas se empujan y se atraen entre sí#

El verdadero tema del artículo no es una sola burbuja, sino el par de burbujas. Cuando una burbuja oscila, irradia ondas de presión hacia el líquido que la rodea. Bajo la hipótesis de incompresibilidad, la presión que la burbuja $j$ genera a una distancia $r$ es proporcional a la aceleración de volumen.

Esta presión se suma a la $p_\infty$ de la burbuja vecina y acopla las dos ecuaciones del radio. El resultado de ese acoplamiento es la fuerza secundaria de Bjerknes (fuerza promediada en el tiempo entre dos burbujas que oscilan). Su dirección la determina la diferencia de fase entre las dos burbujas.

$V_i$ es el volumen de cada burbuja y $D$ es la distancia entre centros. Si las dos burbujas laten en la misma fase, se atraen entre sí; si están en fase opuesta, se repelen. La fase se decide por la relación entre la frecuencia de forzamiento y la frecuencia de resonancia de cada burbuja. Aquí la resonancia es la frecuencia de Minnaert.

Cuanto más pequeña es la burbuja, más alta es $f_0$. Si la frecuencia de forzamiento cae entre las resonancias de las dos burbujas, una responde por debajo de su resonancia (en fase) y la otra por encima (en fase opuesta), de modo que quedan en fase contraria. En ese caso las dos burbujas se repelen.

Abajo es posible variar el radio de las dos burbujas y la frecuencia de forzamiento.

Si los dos radios son parecidos, responden en la misma fase a cualquier frecuencia y aparece la flecha verde (atracción). Al separar una burbuja grande de otra pequeña y luego empujar la frecuencia de forzamiento entre las dos resonancias, la diferencia de fase se abre y todo se invierte hacia la flecha roja (repulsión).

Dónde falla el modelo esférico#

Tanto Keller–Miksis como Gilmore parten de la premisa de "esfericidad". Nagy y Hegedűs señalan cuándo se rompe esa suposición comparándola con simulación numérica directa (DNS) del solver de flujo multifásico ALPACA. Surgen tres conclusiones.

Primera: cuando una burbuja aislada colapsa de forma suave, el modelo esférico es sorprendentemente preciso. Aunque no sea una esfera perfecta, acierta bien la presión de colapso. Segunda: en colapsos violentos, donde el número de Mach de la pared se acerca a 1, Gilmore queda más cerca de la DNS que Keller–Miksis. Esto se debe a que trata la compresibilidad del líquido de forma más fiel mediante la ecuación de estado. Tercera: si un par de burbujas colapsa con fuerza a corta distancia, se forma un chorro (jet: inyección de líquido que penetra por un lado). La esfericidad se rompe y el modelo esférico sobreestima la presión interna. En ese caso se necesita DNS.

Es decir, el modelo esférico es una aproximación barata y, en general, útil. Basta para predecir la fracción de vacío de una nube de burbujas o el rendimiento de un sonorreactor. Solo hay que tener claro su límite cuando una burbuja choca contra una pared o contra una vecina y genera un chorro.

Puntos para recordar#

• Rayleigh–Plesset es el modelo básico de una burbuja esférica, y Keller–Miksis incorpora la velocidad del sonido del líquido a primer orden para atenuar de forma fiel el rebote del colapso.
• El par de burbujas se acopla mediante ondas de presión radiadas, y la diferencia de fase determina el signo de la fuerza secundaria de Bjerknes. Si el forzamiento cae entre las dos resonancias, se repelen; si no, se atraen.
• El modelo esférico es preciso en colapsos suaves, pero cuando se forma un chorro en un colapso violento y cercano, la respuesta es la DNS.