Solve $y''(x)=\alpha y(x)^{-\gamma}-\beta$ when $\gamma>1$, $\alpha>0$, $\beta>0$, $y(0)>0$, $y'(0)=0$

Question

Consider the following ODE $y''=\alpha/y^{\gamma}-\beta$ where $\gamma>1$, $\alpha>0$, $\beta>0$. Initial conditions: $y(0)>0$, $y'(0)=0$.

Actually, $y(t)$ gives the position of a piston being pushed through a cylinder by the adiabatic expansion of some gas (and there is atmospheric pressure on one side of the piston).

The system is in equilibrium when $y''=0$ i.e $y=(\alpha/\beta)^{1/{\gamma}}$ so we expect the solution to be oscillatory about this equilibrium position. I don't think the equation is solvable but Wolfram Alpha gives some weird results which seem to suggest the piston periodically returns to its initial position. Is this true?

Also is it possible to calculate the maximum value of $y$?

Answer 1

For small perturbations $y=y^*+u$ around the stationary point $y^*=(α/β)^{1/γ}$ you get the power series expansion equation \begin{align} u'' = α (y^*+u)^{-γ}-β &= β(1+u/y^*)^{-γ}-β \\ &= -βγ(β/α)^{1/γ}u+O(u^2) \end{align} which in first order is an oscillation equation $$u''+ω^2u=0$$ for frequency $$ω^2=α^{-1/γ}β^{1+1/γ}γ.$$

---

Remember that the Taylor series for $(1+z)^{-γ}$ starts as $$ 1-γz+\frac{γ(γ+1)}2z^2-\frac{γ(γ+1)(γ+2)}6z^3+… $$

Answer 2

((Interesting silent downvote...))

Every solution solves $$\tfrac12(y')^2=y''y'=\alpha y^{-\gamma}y'-\beta y'=(-\alpha (\gamma-1)^{-1}y^{1-\gamma}-\beta y)'$$ hence, along every trajectory starting from $y(0)=y_0>0$ and $y'(0)=0$,

$$\tfrac12y'(x)^2+h(y(x))=h(y_0)\tag{$\ast$}$$

where $h$ denotes the function defined on $z>0$ by $$h(z)=\alpha(\gamma-1)^{-1}z^{1-\gamma}+\beta z$$ Introduce the equilibrium of the system, which is at $$y^*=\left(\frac{\alpha }\beta \right)^{1/\gamma}$$ Since $\alpha >0$, $\beta >0$ and $\gamma>1$, the function $h$ is decreasing from $h(0+)=+\infty$ to $h(y^*)$, then increasing from $h(y^*)$ to $h(+\infty)=+\infty$. For every $z>0$, $z\ne y^*$, let $\hat z$ denote the unique solution $\hat z\ne z$ of the equation $$h(z)=h(\hat z)$$ Then, for every initial condition $y_0$, the trajectory stays confined to the domain $h(\ )\leqslant h(y_0)$, that is, to the interval whose boundary is $$\{y_0,\hat y_0\}$$ For every $y_0<y^*$, the maximal value of $y(x)$ along the trajectory is

$$\hat y_0$$

There is no general expression of $\hat y_0$ in terms of $y_0$ except for some particular values of $\gamma>1$, for example, if $\gamma=2$ then $$\hat y_0=\frac{\alpha}{\beta y_0}$$ In full generality, when $z\to y^*$, $$h(z)-h(y^*)\sim\tfrac12h''(y^*)(z-y^*)^2$$ hence, when $y_0\to y^*$, $$\hat y_0-y^*\sim y^*-y_0$$ Thus, for every $y_0<y^*$ close to $y^*$, the maximal position $\hat y_0$ of the piston is approximately

$$2y^*-y_0$$

Again for every $y_0$ close to $y^*$, $$y''(x)=h'(y(x))\approx h''(y^*)\cdot(y(x)-y^*)$$ hence the frequency $\omega$ of the movement solves $$\omega^2\approx h''(y^*)=\alpha\gamma(y^*)^{-1-\gamma}=\frac{\beta \gamma}{y^*}$$ that is,

$$\omega^2\approx \gamma\frac{\beta^{1+1/\gamma}}{\alpha^{1/\gamma}}$$

Finally, for every $y_0<y^*$, integrating $(\ast)$ from $y_0$ to $\hat y_0$ yields one half of the period $\tau$, hence $$\tau=\sqrt2\int_{y_0}^{\hat y_0}\frac{dz}{\sqrt{h(z_0)-h(z)}}$$