I'm working on simulating this equation (application is motor control, not that it matters):
$$\frac{d^2\theta}{dt^2}+b\frac{d\theta}{dt}=a \sin (x-\theta)$$
where $x = vt$ for $t > 0$, and I'm finding that for given initial conditions $\frac{d\theta}{dt}|_{t=0}$ and $\theta|_{t=0}$, there seems to be a critical value $v_{crit}$ such that:
- if $v < v_{crit}$, $x-\theta$ oscillates but settles down to its equilibrium value $\phi = \sin^{-1} \frac{bv}{a}$ ("capture")
- if $v > v_{crit}$, $x-\theta$ tends to a linearly increasing difference + a small oscillating term.
Predicting the value of $v_{crit}$ is important in my application, and I would like to understand what is going on.
My training in nonlinear differential equations is rather limited and rusty, and I never took advanced classes... I think there might be some insight using energy techniques (Hamiltonians? Lyapunov stability?) since there is similarity to a driven pendulum equation ($\frac{d^2\theta}{dt^2}+\frac{g}{l}\sin \theta = u(t)$), but I can't figure out what.
How can I figure out this critical value?
Can anyone point me towards some reference material (or even the right terms to look up) so I could learn a technique to solve my problem?
It seems like if I can show that $|x-\theta|$ reverses direction before it hits $\pi$, then capture is guaranteed.
okay, splitting into 1st-order systems:
$\begin{eqnarray} \dot{\omega} &=& -a \sin(\theta - x) - b\omega \cr \dot{\theta} &=& \omega \end{eqnarray}$
Change of variable $u = \theta - x$ so $\dot{u} = \omega - v$ and $\ddot{u} = \ddot{\theta} = \dot{\omega}$:
$\begin{eqnarray} \dot{\omega} &=& -a \sin u - b\omega \cr \dot{u} &=& \omega - v \end{eqnarray}$
If I try to write a Lyapunov equation $E = c\omega^2 + d \cos u$ I get
$$\begin{eqnarray} \dot{E} &=& 2c\dot{\omega}\omega - d \dot{u} \sin u \cr &=&2c\omega(-a\sin u - b\omega) - d (\omega - v) \sin u \cr &=&\omega \sin u (-2ac - d) -2bc\omega^2 + dv \sin u \end{eqnarray} $$
I can make the first term go away if I choose d=-2ac; the second term is negative if $c>0$ but I can't get rid of the third term.
Attempt #2: $E = (\omega-v)^2 - 2a \cos u$ I get
$$\begin{eqnarray} \dot{E} &=& 2\dot{\omega}(\omega-v) +2a \dot{u} \sin u \cr &=&2(\omega-v)(-a\sin u - b\omega) +2a (\omega - v) \sin u \cr &=&(\omega-v) \sin u (-2a + 2a) -2b\omega(\omega-v) \cr &=& -2b\omega(\omega-v) \cr &=& -2b\left(\left(\omega-\frac{v}{2}\right)^2 - \frac{v^2}{4}\right) \cr \end{eqnarray} $$
but that's not necessarily negative. Urk.
Some numerical sample points: I'm using a = 2.3086177e5, b = 1.78179103 (max v with equilibrium at $a/b \approx 129567$, but $v_{crit}$ tends to be much smaller in practice), and in my simulations I'm seeing:
$\omega|_{t=0} = 0, \theta|_{t=0} = 0 : v_{crit} \approx 958.929$
$\omega|_{t=0} = 0, \theta|_{t=0} = 0.5 : v_{crit} \approx 929.408$
$\omega|_{t=0} = 0, \theta|_{t=0} = 1.0 : v_{crit} \approx 842.336$
$\omega|_{t=0} = 0, \theta|_{t=0} = \pi/2 : v_{crit} \approx 680.156$
For the last case, here's a phase plot (it's a very underdamped system so the turns of the plotted curve come very close together; the cusp on the left is where things slow down for a moment):
and a timeseries plot:
Assuming what you're saying is correct, a wild guess is $v_{\text{crit}}=a/b$ (just because anything above that would have no arcsin, and your first case cannot work)?