Volví a mirar esto con más cuidado. Todavía no estoy convencido de que sea correcto, pero espero que esto sea al menos mejor que lo que tenía antes.
Dejemos que $h$ sea la altura de la columna de agua dentro de la paja. A medida que esta altura aumenta en una cantidad $\delta h$ el trabajo realizado en la columna es
$$\delta W = P A \delta h$$
donde $A$ es el área de la sección transversal. El cambio de energía potencial es
$$\delta U = \rho A \delta h g h$$
por la conservación de la energía,
$$\delta K = \delta W - \delta U = \left(P - \rho g h\right) A \delta h$$
Este exceso de energía cinética proviene de dos contribuciones: la masa añadida,
$$\delta K_1 = \frac{1}{2}mv^2 = \frac{1}{2}(\rho A \delta h)\dot{h}^2$$
y cualquier cambio en la velocidad de la columna de agua,
$$\delta K_2 = \frac{1}{2}m(2v\delta v) = (\rho A h)\dot{h} \delta\dot{h}$$
Si lo juntamos todo, obtenemos
$$P A \delta h - \rho A g h \delta h = \frac{1}{2}\rho A \dot{h}^2 \delta h + \rho A h\dot{h}\delta\dot{h}$$
Si asume (o demuestra) que $P$ depende de $h$ mediante el teorema de Bernoulli,
$$P + \frac{1}{2}\rho\dot{h}^2 = P_0$$
Sustituyendo (y cancelando el factor común de $A$ ), se obtiene
$$P_0 \delta h - \rho g h \delta h = \rho \dot{h}^2 \delta h + \rho h \dot{h}\delta\dot{h}$$
lo que al menos explica el misterioso factor de $\frac{1}{2}$ que aparecieron en versiones anteriores de mi respuesta.
Ahora, no creo que podamos simplemente asumir que $\delta h \neq 0$ y dividirlo, porque si hacemos eso, obtenemos un factor de $\frac{\delta\dot{h}}{\delta h}$ que es indefinido en las alturas inicial y máxima. (A grandes rasgos, la variación $\delta h$ es de segundo orden en esos puntos mientras que la variación $\delta\dot{h}$ sigue siendo de primer orden). En su lugar, dividiré por $\delta t$ que ciertamente no debería producir ninguna singularidad, para obtener
$$P_0 \dot{h} - \rho g h \dot{h} = \rho \dot{h}^3 + \rho h \dot{h}\ddot{h}$$
En las alturas inicial y máxima, $\dot{h} = 0$ por lo que la ecuación se satisface trivialmente allí. Pero consideremos la situación cuando se desplaza de la altura inicial o máxima una cantidad arbitrariamente pequeña, tal que $\dot{h}\neq 0$ . Aquí podemos anular $\dot{h}$ para conseguir
$$P_0 - \rho g h = \rho \dot{h}^2 + \rho h \ddot{h}$$
Desde $\dot{h}$ será infinitesimalmente pequeño alrededor de la altura máxima, podemos despreciar el primer término de la derecha, pero no el segundo. Así que nos queda
$$P_0 - \rho g h = \rho h \ddot{h}$$
Obsérvese que esto concuerda con un análisis simple utilizando la segunda ley de Newton. (Las fuerzas que actúan sobre la columna de agua en su altura máxima son la fuerza de presión $P_0 A$ actuando hacia arriba y la gravedad $\rho Ahg$ actuando hacia abajo, y la diferencia es igual a $ma = \rho Ah\ddot{h}$ .) Así que la ecuación diferencial pasa al menos una prueba de consistencia básica.
De todos modos, esta ecuación ya no admite la solución $h = \frac{P_0}{\rho g}$ . En cambio, tenemos
$$h = \frac{P_0}{\rho (g + \ddot{h})}$$
Desgraciadamente no se me ocurre una forma de determinar $\ddot{h}$ al máximo sin resolver la ecuación, así que por ahora estoy limitado a una solución numérica.
Para una estimación rápida, he introducido la ecuación diferencial completa de arriba en el programa de Mathematica NDSolve
función. Con las condiciones de contorno $h(0) = 0$ et $\dot{h}(0) = 0$ se quejaba de expresiones no definidas, por lo que utilicé condiciones de contorno en un tiempo no nulo,
$$h(\epsilon_t) = \epsilon_h$$
y
$$\dot{h}(\epsilon_t) = \epsilon_{\dot{h}}$$
para los valores de los distintos $\epsilon$ constantes que van desde $10^{-3}$ a $10^{-8}$ . En mis pruebas, obtengo este gráfico, aparentemente independiente de los valores de $\{\epsilon\}$ o las relaciones entre ellos:
Mathematica indica que el gráfico tiene un pico en $15.5\,\mathrm{m}$ Si este análisis es correcto, esa sería la altura máxima. (Por cierto, todavía estoy muy aunque sospecho de este cálculo)