How to use method of multiple scale in coupled differential equation?

Lets consider a coupled differential equation-

$\ddot{z} + \omega_{o}^{2}z + 2\beta\dot{z} + \frac{\epsilon}{2m}\theta = 0$

and,

$\ddot{\theta} + \omega_{o}^{2}\theta + 2\gamma\dot{\theta} + \frac{\epsilon}{2I}z = 0$

where, $\omega_{o}^{2}$, $I$, $\beta$, $\gamma$ and, $m$ are constants.

My approach:

Very similar to the undamped version, substituting one of the coordinates into the other, but equations become more messy, like (even when assuming $\gamma=0$) -

$\ddddot{\theta}+2\beta\dddot{\theta}+2\omega_{o}^{2}\ddot{\theta}+2\beta\omega_{o}^{2}\dot{\theta}+(\omega_{o}^{4}-\frac{\epsilon^{2}}{4mI})\theta=0$

Assuming solution, $\theta(t)=Ae^{\alpha t}e^{i\omega t}$, because $\theta$ should decay as there is damping in the system, I get following -

$g(\omega,\alpha)=(\alpha+i\omega)^{4}+2\beta\cdot(\alpha+i\omega)^{3}+2\omega_{o}^{2}\cdot(\alpha+i\omega)^{2}+2\beta\omega_{o}^{2}\cdot(\alpha+i\omega)+(\omega_{o}^{4}-\frac{\epsilon^{2}}{4mI})=0$

this would be satisfied only when,

$Re(g(\omega,\alpha))=0$ & $Im(g(\omega,\alpha))=0$,

I get -

$Re(g(\omega,\alpha))=0\Rightarrow\omega^{4}-(6\alpha^{2}+3\alpha \cdot 2\beta-2\omega_{o}^{2})\cdot\omega^{2}+(\alpha^{4}+2\beta\alpha^{3}-2\omega_{o}^{2}\alpha^{2}+2\beta\omega_{o}^{2}\alpha+\omega_{o}^{4}-\frac{\epsilon^{2}}{4mI})=0$

and,

$Im(g(\omega,\alpha))=0\Rightarrow4\alpha^{3}\omega+6\alpha^{2}\beta\omega-\alpha\cdot(4\omega^{3}+4\omega\omega_{o}^{2})+2\beta\omega_{o}^{2}\omega-2\beta\omega^{3}=0$

Equations are too complicated to solve, although $Re(g(\omega,\alpha))= 0$ can be solved with quadratic formula and I get,

$\omega_{1,2}^{2}=\frac{1}{2}[6\alpha^{2}+6\beta\alpha+2\omega_{o}^{2}\pm\sqrt{32\alpha^{4}+64\alpha^{3}\beta+16\alpha^{3}\omega_{o}^{2}+36\alpha^{2}\beta^{2}+16\alpha \beta\omega_{o}^{2}+\frac{\epsilon^{2}}{4mI}}]$

Now I am stuck, how should I proceed with $Im(g(\omega,\alpha))=0$ to get $\alpha$ in terms of other constants?

Also, how can I apply the method of multiple scales (or Poincare-Lindstedt method) in these equations?

Can $(2\beta\dot{z} + \frac{\epsilon}{2m}\theta)$ and $(2\gamma\dot{\theta} + \frac{\epsilon}{2I}z)$ be considered as perturbations if $\beta$, $\gamma$, $\epsilon$ are small?

• Welcome to MSE. It will be more likely that you will get an answer if you show us that you made an effort. – José Carlos Santos Mar 13 '18 at 6:58

$$\begin{cases}\ddot{z} + \omega^{2}z + \frac{\epsilon}{2m}\theta = 0\\ \ddot{\theta} + \omega^{2}\theta + \frac{\epsilon}{2I}z = 0\end{cases}$$ Without approximation (whatever $\epsilon$ is), the system can be solved classically :

$\epsilon\theta =-2m(\ddot{z} + \omega^{2}z)$ that we put into the second equation :

$$-2m(\ddddot{z} + \omega^{2}\ddot z)-2m\omega^2(\ddot{z} + \omega^{2}z) +\frac{\epsilon^2}{2I}z=0$$

$$\ddddot{z} +2\omega^{2}\ddot z +\omega^4z =\frac{\epsilon^2}{4mI}z$$

This is a fourth order linear ODE which can be classically solved.

Nevertheless, if $\epsilon$ is small, you can consider $\frac{\epsilon^2}{4mI}z$ as a perturbation term.

Since it is a single ODE, I suppose that you known how to proceed.

You will obtain the result on the form $\quad z(t)\simeq z_0(t)+\epsilon z_1(t)$.

Then solve $\quad\ddot{\theta} + \omega^{2}\theta + \frac{\epsilon}{2I}z_0 = 0\quad$ on the same manner, in order to obtain $\quad \theta(t)\simeq \theta_0(t)+\epsilon \theta_1(t)$

• You could also solve the characteristic equation $r^4+2ω^2r^2+ω^4=ϵ^2a^2ω^4$ to $r^2+ω^2=\pm ϵaω^2$ so that $r=\pm i\sqrt{1\pm a}ω$ and then find the combinations of the trigonometric basis solutions that solve the original ODE. – LutzL Mar 13 '18 at 10:19
• Thanks for the answer @JJacquelin, Lutzl and sorry. Actually, by mistake I asked wrong doubt. I had already solved this using similar approach and found the analytical solution. My doubt is regarding the damped case (see the edited doubt). I have tried to solve, but unfortunately could not succeed. – luv_phy Mar 14 '18 at 9:03

To put it in more general terms, you have a system $$p(D)x = ϵay\\ q(D)y = ϵbx$$ where $p,q$ are polynomials with constant real coefficients. Eliminating $y$ you get $$q(D)p(D)x=ϵ^2abx$$ If $r$ is a simple root of $p(r)=0$ but not a root of $q$, $q(r)\ne 0$, then there is a root $r+ϵ^2h+O(ϵ^4)$ of the combined system close-by where $h=\frac{ab}{p'(r)q(r)}$. The basis solution to that root is $$x=Ae^{(r+ϵ^2h)t}+O(ϵ^4),\\ y=Be^{(r+ϵ^2h)t}+O(ϵ^4)$$ where $B=ϵA\frac{p'(r)h}{a}=ϵC\frac{b}{q(r)}$. Similar for simple roots of $q$ that are not roots of $p$.

This should cover your situation as long as the friction constants are different. If they are equal, then a simple root $r$ of $p$ is also a simple root of $q$ and the combined polynomial has a root $r+ϵh+O(ϵ^2)$ where $$p'(r)q'(r)h^2=ab\implies h=\pm\sqrt{\frac{ab}{p'(r)q'(r)}}$$ and the basis solution $$x=Ae^{(r+ϵh)t}+O(ϵ^2),\\ y=Be^{(r+ϵh)t}+O(ϵ^2)$$ has its coefficients connected via $p'(r)h\,A=a\,B$, $B=A\frac{p'(r)h}a$.

• Thanks @LutzL. So I should assume solution to be of exponential form and proceed with that like usual ODE, also D is different in q(D) and coefficient D? – luv_phy Mar 14 '18 at 16:24
• If you want to connect the characteristic roots to the roots of the unperturbed problem then something similar has to be done. But in the main yes, the usual exponential approach. And yes, first $D$ is the time derivative operator, in the end the letter is used for a constant factor. – LutzL Mar 14 '18 at 17:13