La cuerda doblada aguantó, la cuerda cortada se pasó un 9% — el final numérico de la disputa de 1747

La forma de representar la solución decide qué tipo de error aparece. Las ondas viajeras transportan la forma; la suma de modos oscila en cada esquina.

En 1747, tres hombres se dividieron por una sola cuerda#

d'Alembert escribió en 1747 la ecuación de la cuerda vibrante. Fue la primera ecuación en derivadas parciales de la historia. La pelea que vino después no era sobre la ecuación. Era sobre qué aspecto podían tener sus soluciones. Euler y Daniel Bernoulli presentaron respuestas distintas, y ninguno de los tres cedió en treinta años.

Los tres tenían razón. Las respuestas solo se separan al truncar la serie en un número finito de términos. Este artículo pone las dos representaciones en el mismo código y observa cómo la suavidad del dato inicial decide el tipo de error. El repique del 9% que aparece en métodos espectrales y esquemas de alto orden empieza aquí.

d'Alembert no desarrolló nada#

Al aplicar la segunda ley de Newton a un tramo pequeño de cuerda con tensión $T$ y densidad lineal $\rho$ se obtiene lo siguiente.

u_{tt} = c^2 u_{xx}, \qquad c = \sqrt{T/\rho}

$u$ es el desplazamiento transversal y $c$ la velocidad de onda. d'Alembert notó que el cambio $\xi = x - ct$ , $\eta = x + ct$ convierte la ecuación en $u_{\xi\eta} = 0$ . Dos integraciones y la solución cae sola.

u(x,t) = \tfrac{1}{2}\left[\,F(x - ct) + F(x + ct)\,\right]

$F$ es el desplazamiento inicial $f(x)$ extendido de forma impar respecto a ambos extremos, con periodo $2L$ . La condición de extremo fijo nunca se impone aparte. La propia extensión impar invierte el signo en cada pared. Sin desarrollo, sin coeficientes, sin frecuencias. Toda la solución consiste en partir la forma inicial por la mitad y llevarla hacia ambos lados.

Conviene manipular la simulación siguiente.

speed0.35

Turn halves off and you see only the string. Turn it on and the same picture is two rigid copies at half height sliding in opposite directions — no frequencies anywhere. Switch to jump: the corners stay razor sharp forever, because transport never smooths anything. The clock reads t·c/L = 0.00; at 2 the shape is exactly back.

Al apagar y encender halves, la única curva blanca resulta ser dos copias de media altura que se deslizan en sentidos opuestos. En el perfil jump, las esquinas siguen afiladas hasta el final. La advección no suaviza nada.

Bernoulli reescribió la misma respuesta con senos#

La objeción de Daniel Bernoulli vino de la música. Una cuerda suena con su fundamental y sus armónicos a la vez. Entonces la solución también debe ser una superposición de ondas estacionarias.

u(x,t) = \sum_{n=1}^{\infty} b_n \sin\!\frac{n\pi x}{L}\cos\!\frac{n\pi c t}{L}, \qquad b_n = \frac{2}{L}\int_0^L f(x)\,\sin\!\frac{n\pi x}{L}\,dx

$b_n$ es la amplitud del modo $n$ y $n\pi c/L$ su frecuencia angular. Que las frecuencias modales estén en la razón entera $1:2:3$ se sabía desde Pitágoras. En el París de entonces se discutía la teoría de la armonía de Rameau, y el propio d'Alembert estaba en el bando que la defendía con mecánica newtoniana. La base física de la escala musical salió de esta ecuación.

Lo que Euler cuestionó fue la palabra "función"#

Euler se puso del lado de d'Alembert, pero por otro motivo. Al pulsar una cuerda con el dedo, la forma inicial es un triángulo con un pico. En ese punto $u_{xx}$ no existe. d'Alembert ni siquiera quería admitir esa curva como solución, porque para él solo era función una curva expresable mediante una única fórmula analítica. Euler insistió en que cualquier curva dibujada a mano es un dato inicial legítimo.

Bernoulli fue aún más lejos: toda curva, afirmaba, se puede escribir como suma infinita de senos. En su momento parecía optimismo sin fundamento. La cuestión no se zanjó hasta Fourier, en 1822. Igual que d'Alembert tuvo que aceptar un flujo potencial alrededor de un cilindro que da arrastre nulo, aquí también su intuición acertaba solo a medias.

Python, con las dos representaciones en el mismo instante#

Se preparan tres desplazamientos iniciales: una función caja con saltos, un triángulo con solo un pico y una campana suave. En el mismo instante $tc/L = 0.15$ se compara la solución de d'Alembert con la serie de senos de $N$ términos.

import numpy as np
 
L, c, T = 1.0, 1.0, 0.15
X = np.linspace(0.0, L, 2001)
 
def hat_profile(x, a=0.35, b=0.55):
    """Salto: vale 1 en [a,b] y 0 fuera - la banda que golpea el martillo"""
    return np.where((x >= a) & (x <= b), 1.0, 0.0)
 
def kink_profile(x, a=0.45):
    """Pico: continua, pero la pendiente salta en x=a - cuerda pulsada"""
    return np.where(x < a, x / a, (L - x) / (L - a))
 
def bell_profile(x, a=0.45, s=0.055):
    """Suave: campana infinitamente derivable"""
    return np.exp(-((x - a) ** 2) / (2 * s ** 2))
 
def odd_extend(f0, xq):
    """Interpola sobre la extension impar y 2L-periodica que exigen los extremos fijos"""
    xs = np.mod(xq, 2 * L)
    sign = np.where(xs > L, -1.0, 1.0)
    xs = np.where(xs > L, 2 * L - xs, xs)
    return sign * np.interp(xs, X, f0)
 
def dalembert_wave(f0, t):
    """Solucion de d'Alembert - dos ondas de media altura que viajan en sentidos opuestos"""
    return 0.5 * (odd_extend(f0, X - c * t) + odd_extend(f0, X + c * t))
 
def modal_wave(f0, t, n_modes):
    """Solucion de Bernoulli - superposicion de n_modes modos senoidales"""
    n = np.arange(1, n_modes + 1)[:, None]
    k = n * np.pi / L
    b = 2.0 / L * np.trapezoid(f0[None, :] * np.sin(k * X[None, :]), X, axis=1)
    return (b[:, None] * np.sin(k * X[None, :]) * np.cos(k * c * t)).sum(axis=0)
 
def overshoot_pct(u, exact):
    """Exceso maximo sobre el pico exacto, en % de la amplitud exacta"""
    return 100.0 * (u.max() - exact.max()) / (exact.max() - exact.min())
 
for name, f0 in (("jump  ", hat_profile(X)),
                 ("kink  ", kink_profile(X)),
                 ("smooth", bell_profile(X))):
    exact = dalembert_wave(f0, T)
    print(f"[{name}]  t*c/L = {T}")
    for n_modes in (8, 32, 128, 512):
        u = modal_wave(f0, T, n_modes)
        print(f"   N={n_modes:4d}   max|modal - dAlembert| = {np.abs(u - exact).max():.5f}"
              f"   overshoot = {overshoot_pct(u, exact):+6.2f} %")

[jump  ]  t*c/L = 0.15
   N=   8   max|modal - dAlembert| = 0.30480   overshoot = +21.01 %
   N=  32   max|modal - dAlembert| = 0.26454   overshoot = +12.07 %
   N= 128   max|modal - dAlembert| = 0.23847   overshoot =  +9.86 %
   N= 512   max|modal - dAlembert| = 0.18738   overshoot =  +8.94 %
[kink  ]  t*c/L = 0.15
   N=   8   max|modal - dAlembert| = 0.02073   overshoot =  -0.23 %
   N=  32   max|modal - dAlembert| = 0.00639   overshoot =  -0.19 %
   N= 128   max|modal - dAlembert| = 0.00157   overshoot =  -0.03 %
   N= 512   max|modal - dAlembert| = 0.00038   overshoot =  -0.01 %
[smooth]  t*c/L = 0.15
   N=   8   max|modal - dAlembert| = 0.04435   overshoot =  -6.25 %
   N=  32   max|modal - dAlembert| = 0.00000   overshoot =  -0.00 %
   N= 128   max|modal - dAlembert| = 0.00000   overshoot =  +0.00 %
   N= 512   max|modal - dAlembert| = 0.00000   overshoot =  -0.00 %

Los tres bloques cuentan tres historias distintas. La campana suave toca el cero de doble precisión con 32 términos. El triángulo con pico divide su error por cuatro cada vez que $N$ se multiplica por cuatro: primer orden. Solo el salto conserva un sobrepico, y este se asienta cerca del 9% en lugar de reducirse.

El 9% no baja aunque se añadan términos#

Ese número es el fenómeno de Gibbs. Wilbraham lo encontró en 1848, Gibbs lo redescubrió en 1899 y el nombre quedó. El valor teórico es el 8.95% del salto. El 8.94% que imprime la corrida con 512 términos es esa constante.

Lo importante es que el sobrepico se estrecha sin bajar nunca de altura. Al añadir términos, la zona de oscilación se encoge como $1/N$ . Por eso una integral, o una norma $L^2$ , sí converge. Medido en norma del máximo, no. Esa brecha entre normas es lo que en la práctica aparece como densidad negativa.

modes N24speed0.25

Drag modes N to the right on jump: the red ringing narrows but its height parks near 9% (now 0.00%), and the right-hand curve flattens out. Switch to kink and then smooth with N untouched — the same slider that bought nothing now buys everything, and the error curve turns into a cliff (max error 0.0000).

Conviene llevar el deslizador modes N hasta el extremo derecho. En jump el repique rojo se afina y mantiene su altura, y la curva de convergencia de la derecha se aplana por abajo. El mismo deslizador sobre smooth hace que la curva caiga por un precipicio. Lo único que cambió fue el dato inicial.

La tasa de decaimiento de los coeficientes lo explica todo. Con un salto, $b_n \sim n^{-1}$ ; con solo un pico, $b_n \sim n^{-2}$ ; con un perfil suave, más rápido que cualquier potencia. La cola truncada es el error, así que cuanto más gruesa la cola, mayor la cicatriz que deja el corte.

1755, cuando el mismo hombre se topó con lo no lineal#

Euler escribió el movimiento de un fluido como ecuación en derivadas parciales por primera vez en 1755: las ecuaciones de Euler. A diferencia de la ecuación de ondas, aquí la pendiente de una característica depende de la propia solución. Por suave que sea el dato inicial, el cruce de características genera una discontinuidad en tiempo finito. Es la misma estructura que explica por qué un flujo supersónico no sabe nada de lo que hay aguas arriba.

Por eso, en cálculo compresible, elegir un dato inicial suave no sirve de nada. La onda de choque fabrica su propio salto, y desde ese instante un esquema de alto orden vuelve al problema de 1747. Que von Neumann emborronara los choques a propósito en 1950 apuntaba exactamente a este repique. Los limitadores TVD y los pesos WENO bajan el orden solo cerca del choque, porque recuperar suavidad de forma local es la única manera de borrar ese 9%.

Revisar primero la suavidad del dato inicial#

Cuando un esquema nuevo repica, conviene mirar el dato inicial y el de contorno antes de sospechar del esquema. ¿Se sembró el campo inicial como constantes por celda? ¿Se metió la interfaz como un salto? ¿El perfil de entrada es solo $C^0$ en el tiempo? Si alguna respuesta es afirmativa, esa oscilación no es un error de código. Es el precio de la representación.

El mismo criterio elige los casos de verificación. Una solución suave devuelve el orden de diseño sin sorpresas. En cuanto entra un salto, la convergencia en norma del máximo desaparece y solo queda $L^1$ . Los tres hombres de 1747 no tenían vocabulario para esa distinción. Nosotros la llamamos elección de norma.