How to estimate the oscillation range of this differential equation Consider the system of differential equation
$$\left\{\begin{aligned}\frac{\mathrm{d}x}{\mathrm{d}t} &= y\\\frac{\mathrm{d}y}{\mathrm{d}t} &= x - x^3 - y^3\end{aligned}\right.$$
with initial value $x(0) = y(0) = 1$.
Numerical simulation shows that $(x, y)$ oscillates around $(1, 0)$, with no sign of convergence even at $t = 1000$. I would like to estimate the range of this oscillation.
My attempt: Let $t_0$ be the smallest positive $t$ at which $y(t) = 0$. Then when $0 \leq t \leq t_0$ we have $y > 0$, hence $x$ is increasing, so $x \geq 1$. Then $x - x^3 < 0$, so $y$ is decreasing, $y \leq 1$.
Note that
$$\frac{\mathrm{d}x}{\mathrm{d}t} + \frac{\mathrm{d}y}{\mathrm{d}t} = x + y - x^3 - y^3 = (x + y)(1 - x^2 - xy - y^2)$$
Since $x \geq 1$, we have $1 - x^2 \leq 0$, so $x + y$ is decreasing, so $x + y \leq 2$. We have $\mathrm{d}y/\mathrm{d}t \geq x - x^3 - (2-x)^3 =
 -6x^2+13x-8$.
Let $t_1$ be the smallest positive $t$ such that $x(t_1) = 7/6$. When $0 \leq t \leq t_1$, we have $\mathrm{d}y/\mathrm{d}t \geq -1$, so $y \geq 1-t, x \geq -t^2/2 + t + 1$. Hence $t_1 \leq 1 - \sqrt{6}/3$.
On the other hand $t_0 \geq 1 - \sqrt{6}/3$, because whenever $t \leq t_0$ and $t \leq 1/2$, we have $x \leq 1 + t \leq 3/2$, so $\mathrm{d}y/\mathrm{d}t \geq -2$, hence $y \geq 1 - 2t$. If $t_0 < 1/2$, then $y(t_0) \geq 1 - 2t > 0$, contradicting $y(t_0) = 0$. Hence $t_1 \leq 1 - \sqrt{6}/3 < 1/2 \leq t_0$.
I'm not sure how to proceed. Numerical experiment seems to show that $\mathrm{d}y/\mathrm{d}t \leq -1$ when $t_1 \leq t \leq t_0$, but I don't know how to prove that. If that's the case, then since $y(t_1) \leq 2 - x(t_1) = 5/6$, we have $t_0 \leq t_1 + 5/6 \leq 11/6 - \sqrt{6}/3$. Since $y(t_1 + t) \leq 5/6 - t$, we have $x(t_1 + t) \leq 1 + t_1 + 5/6t - 1/2t^2$. Plug in $t_1 \leq 1 - \sqrt{6}/3$ and $t \leq 5/6$ to get $x(t_0) \leq x_0$ where $x_0 \approx 1.531$. That's pretty close to what I saw in the numerical experiment.
How to estimate the next $x$ value where $\mathrm{d}x/\mathrm{d}t = y = 0$? The long-term behavior of this equation seems to be that $x$ oscillates between these two extreme values. Also, can we say anything about the behavior when $t$ goes to negative? Mathematica thinks this is a stiff system when $t \leq -0.55$
 A: Mechanical interpretation
As per the comment, you can see that as a second order equation of the form
$$
x''+(x'^2)x'+(x^3-x)=0,
$$
which can be interpreted/visualized as a mechanical system with a conservative force $F(x)=-V'(x)$ where $V=\frac14(x^2-1)^2$ and a friction/dissipation term that has a coefficient $c(x,x')=x'^2\ge 0$. This means that the system is continuously losing energy as long as it moves, crossing the level curves of
$$
E(x,x')=\frac12x'^2+V(x)
$$
towards ever lower energies to end up in one of the minima $x=\pm 1$ of the potential energy.
$$
\frac{d}{dt}E(x,x')=x'(x''+V'(x))=-x'^4.
$$
However, the closer to the minimum the solution gets, the slower it gets and thus the smaller the energy loss term. Thus it may look as if the solution becomes periodic. In the other time direction the dissipation term now "collects" energy from the environment, making the solution behave more and more lively. Now the 4th power of $x'$ gives a positive feed-back which in general leads to a divergence to infinity in finite time. Then also all derivatives take large values, so the system may appear as "stiff", but more correctly it would be "run-away".
A first approach to the energy decay
To be a little more precise, close to the equilibrium one can try $x'^2\approx 2E$ for the qualitatively close equation
$$
E'=-4E^2\implies E=\frac1{c+4t}.
$$
This is likely quantitatively deviating from numerical results, what should stay true is that $E^{-1}$ is approximately linear.
A second approach using averaging over a period
As the distance to the equilibrium decays slowly, one can average over periods. The linearization close to $x=1$ gives $u''+2u=0$ for $x=1+u$. This has solutions $x=1+r\sin(\sqrt2 t+\phi)$ and $x'=\sqrt2r\cos(\sqrt2 t+\phi)$ so that $E=r^2+O(r^3)$. Then the average of
$$
x'^4=4r^4\cos^4(\sqrt2 t+\phi)=\frac{r^4}2(\cos(4\sqrt2 t+4\phi)+4\cos(2\sqrt2 t+2\phi)+3)
$$
over the period $\frac{2\pi}{\sqrt2}$ is the constant term, $\frac32r^4$. This implies that $E^{-1}$ or $r^{-2}$ grows like $c+\frac32t$, which is slower than the first guess above.
The plot of $E^{-1}-\frac32t$ about confirms this behavior.

