Understanding solutions to the Mathieu Equation in the context of Classical Mechanics

Question

I have have two questions concerning the following system in Classical Mechanics. I will list them here and then give some context in which I will note where the questions occur.

The Questions

(1) Is it the case that the Mathieu Characteristic Functions have a "frequency" half that of the parametric driving term? If so, are there conditions that must be met or is it a general fact?

(2) Is there any way to facially distinguish between periodic and divergent solutions to this system? Why do small changes in the parameters of the Mathieu equation sometimes result in huge transitions in the solution space?

The Problem

Consider a simple pendulum whose hinge is attached to the end of a Hooke's Law spring hanging from the ceiling. Assume that the spring and the hinge of the pendulum are both massless and that the the hinge (and spring) is constrained to vertical motion only. The pendulum is free to swing in the plane.

The acceleration due to gravity is $g$, the length of the pendulum is $l$, the mass of the pendulum bob is $m$, and the spring constant is $k$. Taking $\theta$ to be the angle with which the pendulum swings from the vertical and $D$ to be the displacement of the spring from its equilibrium length, this system has the following Lagrangian: $$ \mathscr L=\frac{1}{2}m\dot D^2+\frac{1}{2}ml^2\dot\theta^2+ml\sin(\theta)\dot D \dot\theta - \frac{1}{2}kD^2-mgD+mlg\cos(\theta) $$

Applying the Euler-Lagrange Equation for each dynamical variable yields the following equations of motion: $$ \ddot\theta=-\frac{\sin(\theta)}{l}\ddot D-\frac{\sin(\theta)}{l}g\\ \ddot D+l\cos(\theta)\dot\theta^2+l\sin(\theta)\ddot\theta=-\frac{k}{m}D-g $$

Substituting the first equation into the second to eliminate $\ddot\theta$ and applying the small angle approximation (ignoring all $\theta$ terms of quadratic order or higher) yields $$ \ddot D+l\dot\theta^2=-\frac{k}{m}D-g $$ We can eliminate the non-linear term by requiring $\dot\theta\ll\sqrt{\frac{g}{l}}$. The remaining equation is the simple harmonic motion equation with a constant driving force. Solving this equation with appropriate initial conditions yields $$ D(t)=\frac{mg}{k}\cos(\sqrt{\frac{k}{m}}t)-\frac{mg}{k} $$

Taking two derivatives and substituting this solution into the first differential equation above (still under the small angle approximation, of course) yields $$ \ddot\theta+(\frac{g}{l}-\frac{g}{l}\cos(\sqrt{\frac{k}{m}}t))\theta=0 $$

This is a slightly modified Mathieu Equation which produces the following even solution given appropriate initial conditions (i.e. $\theta(0)=\theta_0$) $$ \theta (t)=\theta_0 C(\frac{4mg}{kl},\frac{2mg}{kl},\sqrt{\frac{k}{4m}}t) $$

I have been told that the pendulum ought to naturally oscillate with a frequency half that of the spring's. The reasoning being that the pendulum passes through its equilibrium position (where the tension in the rod is greatest) twice in every period. This position should correspond to the maximum downward displacement of the spring since the hinge is massless.

Assuming that what is meant by "frequency" is related to the pendulum's perceived, visual swing rather than a true mathematical repetition, then a little testing produces examples where this appears to be at least approximately true (example 1, example 2). However, I have heard that this halving of the frequency is closely related to Floquet's Theorem, but that doesn't seem right since it is easy to find examples that definitely don't have half the "frequency" if one does not restrict one's solutions too much (example 3).

Moreover, it is simple to find solutions which quickly diverge, and surprisingly, these solutions are not very different from bounded solutions (example 4, example 5).

The Questions Re-stated

Is it possible to identify bounded or divergent solutions simply by looking at their parameters? Is there a hard-and-fast relationship between the frequency of the spring and the "frequency" of the pendulum? How does this relationship relate to Floquet's Theorem?

Answer 1

Handbooks of mathematical functions are truly great tools. Abramowitz and Stegun is nice. More nice (for online viewing) is the NIST DLMF. I find myself using these most often when it's time to play with special functions.

In particular, chapter 20 in (the above linked version of) Abramowitz and Stegun, and DLMF chapter 28 are all about Mathieu functions. Everything I write here can be found in one of those references.

Now, regarding your questions.

Is it possible to identify bounded or divergent solutions simply by looking at their parameters?

The above picture shows where stable (bounded - colored) and unstable (unbounded - white) solutions exist in the 2 parameters of the Mathieu equation (separated by the characteristic curves). The model you're looking at has parameters that lie on the line $a=2q$ by this figure's notation. Picturing that line, the intervals of stability and instability are quite irregular, so it's no easy feat to determine which reregion you will be in for an arbitrary point on that line.

The characteristic curves can be computed by a Numerical continuation method by noting that these curves cross through the $q=0$ axis at square integer values of $a$ (also a fun fact you can find in said handbooks). This means that the stable and unstable regions on the line $a=2q$ can be computed reliably, if one were interested in doing so.

Is there a hard-and-fast relationship between the frequency of the spring and the "frequency" of the pendulum?

Definitely not hard and fast. As usual, it's complicated. Mathieu functions can be written in the form $F_\nu(t) = e^{i\nu t}P(t)$ where $P$ is a periodic function (with period $\pi$ in its standard form). So if we can get $\nu$, known as the characteristic exponent, then that really can give us all info about the solution.

However, how we write the Mathieu function in terms of exponentials boils down to the value of $\nu$. We can write $$ F_\nu(t) = \sum_{k=-\infty}^\infty c_ke^{i(\nu+2k)t} $$ due to the ambiguity in defining the periodic function $P$: $e^{i2kt}$ can be absorbed into the coefficients because they have period $\pi$. We can write the general solution to the ODE as $$ \theta(t) = AF_\nu(t) + BF_\nu(-t) $$ if $F_\nu(t)$ and $F_\nu(-t)$ are linearly independent (Note they are unless your parameters lie on one of the characteristic curves, which lead to $\nu$ being an integer).

Now we have a few different cases of what can happen depending on $\nu$. Because of the contributions of both $F_\nu(t)$ and $F_\nu(-t)$ in how we've written the general solution, if we have $AB\neq 0$, then solutions will blow up somewhere if $\nu$ is complex! This is precisely what is happening in the unstable regions of the above figure.

If $\nu$ is real (this happens in the stable regions shown above), then we have a few possibilities:

1. $\nu$ is an integer (won't go into detail, we can't write the solution as we assumed, but can construct a separate linearly independent solution)
2. $\nu$ is a fraction -> leads to periodic solutions (Exercise: What are their periods?)
3. $\nu$ is irrational -> bounded solution

That's all great, but how do we compute $\nu$?

How does this relationship relate to Floquet's Theorem?

In short, I think someone you've been talking to is trying to take too much from the statement

[Floquet's theorem] gives a coordinate change $y = Q^{−1} ( t ) x$ with $Q ( t + 2 T ) = Q ( t )$ that transforms the periodic system to a traditional linear system with constant, real coefficients.

The way we wrote the solution as a summation of exponentials is equivalent to what is used in Floquet analysis. Following from our linear combination for the general solution, it can be seen that the period of solutions can take on arbitrary values.

It's bedtime for me now, but perhaps I can elaborate on this in the days to come.