# Uniform perturbative solutions to the Mathieu equation

The Mathieu equation is a second-order linear differential equation given by $$y''(t) + [a - 2q\cos(2t)]y(t) = 0$$ There are two special functions defined as linearly independent solutions to Mathieu's equation:

1. The even Mathieu cosine function which I will denote as $MC(a,q,t)$, satisfying the initial conditions $MC(a,q,0)=1$ and $MC'(a,q,0)=0$ (where prime denotes $t$ derivative).

2. The odd Mathieu sine function, $MS(a,q,t)$, satisfying $MS(a,q,0)=0$ and $MS'(a,q,0)=1$.

For $q=0$, Mathieu's equation reduces down to the familiar differential equation $$y''(t) + ay(t) = 0 \tag{*}$$ for which we have even and odd solutions $\cos(\sqrt{a}t)$ and $\sin(\sqrt{a}t)$. Let $\epsilon = \frac{2q}{a}$ and suppose that $2q \ll a$ so that $\epsilon \ll 1$. Therefore Mathieu's equation can be written as $$\frac{1}{a}y''(t) + [1-\epsilon \cos(2t)]y(t) = 0$$ where we view it as a perturbed form of equation $(*)$. It follows that for the case of small $\epsilon$, the functions $MS$ and $MC$ should be approximated by $\sin$ and $\cos$ plus small perturbative terms.

The problem with such a perturbative solution is the appearance of secular terms at $\mathcal{O}(\epsilon^2)$, i.e. terms of the form $t\cos(\sqrt{a} t)$ or $t\sin(\sqrt{a} t)$ which diverge with large $t$, whereas the functions $MC$ and $MS$ appear to be uniformly bounded.

I have tried several standard methods for obtaining uniform perturbative solutions, i.e. the Poincaré - Lindstedt method and multiple scale analysis but neither method successfully gets rid of the secular terms. Does anyone know how I can (or if it's even possible) to obtain an approximation for $MC$ and $MS$ which is reasonably accurate for all $t$?

• One approach is to expand $y(t)$ in a Fourier expansion in some number of modes, then solve the resulting system of equations. That's very good as a numerical method (since it turns it into an eigenvalue problem) but it also works as an analytic approach if you only do a few modes (more than that and it becomes unwieldy.) – Semiclassical Aug 8 '14 at 11:58

It is possible to get a leading-order approximation to the solution via multiple scales relatively easily, as follows. First, for convenience, define $$\tau = a^{1/2}t$$, such that the equation reads $$y''(\tau)+(1-\epsilon \cos(2a^{-1/2}\tau))y(\tau)=0.$$
Since the unperturbed problem has a constant frequency, we might suspect a multiple scales approach to work. Define a slow time scale $$T=\epsilon \tau$$ and consider $$y\equiv y(\tau,T)$$, such that $$\frac{\mathrm{d}}{\mathrm{d}\tau} = \partial_\tau + \epsilon \partial_T$$. The ODE becomes $$(\partial_\tau+\epsilon\partial_T)^2 y(\tau,T) + (1-\epsilon\cos(2a^{-1/2}\tau))y(\tau,T)=0.$$ Expand the unknown function $$y(\tau,T)=y_o(\tau,T)+\epsilon y_1(\tau,T)+o(\epsilon)$$. At order one, we find the unperturbed problem $$\partial_\tau^2y_o+y_o=0,$$ whose solution is $$y_o(\tau,T)=A(T)e^{i\tau}+c.c.$$, where $$c.c.$$ means "complex conjugate". At order $$\epsilon$$, $$\partial_\tau^2y_1+y_1 = -2 \partial_\tau \partial_T y_o + \cos(2a ^{-1/2}\tau)y_o.$$ The RHS is $$-2\partial_\tau \partial_T y_o + \cos(2a ^{-1/2}\tau)y_o = -2(iA_T e^{i\tau} +c.c.)+\frac{1}{2} \left(e^{2a^{-1/2}i\tau}+e^{-2a^{-1/2}i\tau} \right) \left( A e^{i\tau} +c.c.\right)$$ There are two cases to be distinguished at this point. If $$|2a^{-1/2}\pm 1| \neq 1$$, then there are no resonances stemming from the product term. In this case, the solvability condition gives simply $$A_T =0,$$ s.t. the leading-order result is $$y_o(\tau,T)\equiv y_o(\tau)= A e^{i\tau} +c.c.=A e^{i\sqrt{a}t} +c.c.$$, with $$A\in \mathbb{C}$$ a constant. If $$a=1$$, however, then there's a resonant contribution from the product term and the solvability condition reads $$\partial_T A = - \frac{1}{4}i A^*$$ which, upon differentiating again w.r.t. $$T$$, yields $$\partial_T^2 A = \frac{1}{16} A$$, which is easily solved by $$A(T) =(\alpha-i\beta)e^{T/4}+(\gamma-i\delta)e^{-T/4}$$ for real constants $$\alpha,\beta,\gamma,\delta$$. Plugging this into the zeroth order solution gives $$y_o(\tau,T)= \left(\alpha e^{T/4}+\gamma e^{-T/4}\right) \cos(\tau) + \left(\beta e^{T/4} + \delta e^{-T/4}\right)\sin(\tau).$$ Using that in this case $$t=\tau$$, by definition $$T=\epsilon \tau$$ back in, $$y_o(t,T)= \left(\alpha e^{\epsilon t/4}+\gamma e^{-\epsilon t/4}\right) \cos(t) + \left(\beta e^{\epsilon t/4} + \delta e^{-\epsilon t/4}\right)\sin(t).$$ Both approximations at $$a=1$$ or $$a\neq 1$$ are valid until $$t=O(1/\epsilon)$$. If higher order accuracy is needed, one can introduce a third time-scale $$\tilde{T}=\epsilon^2 t$$ and go one order further.
• Under the re-definition of $t$, you would need to change it in all places, thus $y''+(1-ϵ\cos(2a^{1/2}t))y=0$. Thus you get cancellation in your method if $(a^{1/2}-1)t$ is on the scale of $T=ϵt$, and by $A_T=0$ in the non-resonant cases. – Lutz Lehmann Jan 6 at 13:48