I'm trying to figure out how my professor arrived at the following first order system for the Van der Pol equation $x''+ c(x^2-1)*x' + x = 0$.
It's supposed to be equivalent to the first-order system $$y'=-x/c$$ $$x'=c(y+x-x^3/3)$$ however, we are assumed to follow the reasoning without explanation and unfortunately, I do not. Is anyone able to illuminate this a bit for me? Thanks!
[Math] Second-order Van der Pol to First Order System
ordinary differential equations
Related Solutions
$$ \frac{d^2x}{dt^2} - a(1 - x^2)\frac{dx}{dt}+ x = 0 $$
let $\frac{dx}{dt} = y $ . which yields the two dimensional first order system:
$$\frac{dx}{dt} = y $$
$$\frac{dy}{dt} = -x + a(1 - x^2)y $$
The linearized system is easy to write down in this case:
$$\frac{dx}{dt} = y $$
$$\frac{dy}{dt} = -x + ay $$
clearly (0,0) is the equilibrium point
a plot of the equation near the origin with a as parameter . (You can play around with this quite a bit). The red solution curve is the Van der Pol Equation, the blue solution curve is the linearized system. It is spiraling out from the origin , but without a limit cycle.
$\newcommand\rx{\mathfrak{x}}\newcommand\ru{\mathfrak{u}}\newcommand\mA{\mathfrak{A}}$
) For the definition of the most common stability properties one refers to the properties of solutions of linear systems of ODE. Apply the RK4 method to some linear system $\ru'=\mA\ru$ (for instance the linearization of a non-linear ODE system around a fixed point). You should get exactly the first 5 terms (up to degree 4) of the exponential series.
The stability condition is that this propagator has to be non-expanding if the spectrum of $\mA$ is contained in the negative half-plane of the complex plane. On the level of singular eigenvalues $\lambda$ and $z=λh$ with step size $h$ this means that $z$ is in the stability region if $$|1+z+z^2/2+z^3/6+z^4/24|\le1.$$
The half-circle with $\Im(z)\le 0$, $|z|\le 2.5$ is contained in that set.
If $L$ is a Lipschitz constant (or $L=\|\mA\|$), then this gives a bound $h<2.5/L$. At the upper bound the numerical solution is merely not totally wrong, for a useful solution one needs a much smaller step size, $Lh\in[10^{-3},10^{-1}]$. Step sizes much smaller than that lead to an accumulation of floating point errors that dominate the method error.
See Maximum timestep for RK4 for visualizations of the stability region and the effect of too small or too large step sizes.
) There are two readily available strategies to explore the error picture of a numerical ODE solver.
One proceeds to compute a series of results over a series of doubled step sizes $Lh=(1,2,4,8)·10^{-2}$. Compare these to the expected formula $y_h=y^*+C·h^4+O(h^5)$, that is, apply Richardson extrapolation to compute values for $C$ and $y^*$ and check if they are stable over all the doublings.
Another strategy uses MMS (manufactured solutions). Here one would fix some nice function $p$ and consider the equation $$ y''+μ(y^2−1)y'+y=p''+μ(p^2−1)p'+p, ~~~y(0)=p(0),~~y'(0)=p'(0), $$ for instance $p(t)=2\sin(t)$ or with the next perturbation term added. Then $p$ is the exact solution that the numerical solution can be compared against.
If you want to actually prove the order, that is, verify the order conditions, read up on Butcher trees, for instance on Butcher's homepage https://www.math.auckland.ac.nz/~butcher/ODE-book-2008/Tutorials/. Or read the original derivation of the low order Runge-Kutta methods in W. Kutta (1901) https://archive.org/details/zeitschriftfrma12runggoog/page/n449
Best Answer
The first-order system : $$y'=-\frac{x}{c}$$ $$x'=c(y+x-\frac{x^3}{3})$$ leads to : $$x''=c(y'+x'-x^2x')$$ $$x''=c(-\frac{x}{c}+x'-x^2x')$$ $$x''+ c(x^2-1)x' + x = 0$$ This non-linear second order ODE of autonomous kind can be reduced to a first order non-linear ODE.