Question about driven harmonic oscillator

[kelly0303]

Hello! I wanted to make sure I am doing this right. In general for a driven damped harmonic oscillator:

$$\ddot{x}+\gamma\dot{x}+\omega_0^2x=\frac{F}{m}\sin{\omega_dt}$$
where $\gamma$ is the damping, $\omega_0$ is the resonant frequency, $m$ is the mass of the particle, $F$ is the applied force and $\omega_d$ is the frequency of the applying force. Assuming $F=0$ (and $\gamma<<\omega_0$, which is my case) the homogeneous solution is:

$$x(t)=Ae^{-\gamma t/2}\sin{\omega_0t+\phi}$$
If we include the force, the final solution is:

$$x(t) = \frac{F}{m}\left[\frac{\omega_0^2-\omega_d^2}{(\omega_0^2-\omega_d^2)^2+(\gamma\omega_d)^2}\sin{\omega_dt} + \frac{\gamma\omega_d}{(\omega_0^2-\omega_d^2)^2+(\gamma\omega_d)^2}\cos{\omega_dt} \right]$$
So the final solution is the sum of the 2 expressions above. In my situation, in the ideal case, I want to start with a particle at the very center and with zero velocity and by applying this force I want to take it to a fixed amplitude, $A_d$ then let it evolve freely. In my case I also have that $\omega_0 = \omega_d$ and $\gamma$ is basically zero. Thus I can assume that the homogeneous solution is identically zero, while for the driven part (for the right value of $F_d$) I can ignore the second term, given that $\gamma \sim 0$ and end up with:

$$x(t) = A_d\sin{\omega_dt} =A_d\sin{\omega_0t}$$

Thus I will need to simply apply this force for $T_0 = \frac{\pi}{2\omega_0}$.

(Is this right?)

Now, let's assume that instead of starting with zero velocity and position, the particle has a random initial position and velocity. I know that the possible amplitude of this motion, $A$, is much smaller than $A_d$ i.e. $A<<A_d$. What I want to know is the uncertainty in the final amplitude in this case. I assume that in this case, the homogeneous solution can be in general parameterized as:

$$x(t)=A\sin{(\omega_0t+\phi)}$$
and under the same assumption as above, the general solution is given by:

$$x(t)=A\sin{(\omega_0t+\phi)}+A_d\sin{\omega_0t}$$

(Is this right?)

The applied force will be for the same amount as before (as I don't know the initial parameters), so I end up with:
$$x(T_0)=A\sin{(\pi/2+\phi)}+A_d$$
$$\dot{x}(T_0)=A\omega_0\cos{(\pi/2+\phi)}$$

For $\phi=0$, the uncertainty in the amplitude would simply be $x(T_0)-A_d = A$. In general, the new amplitude would be defined as the point where the initial velocity is zero, which means $\pi/2+\phi+\omega_0t = 0$ or $\pi/2+\phi+\omega_0t = 3\pi/2$ (depending on the sign of $\phi$). Which gives an initial position of $\pm A + A_d$, so even in the general case, the uncertainty in the amplitude is given by $A$.

(Is this right?)

Thank you!

 

[vanhees71]

No, it's not right. The final general solution is the superposition of the solution of the homogeneous equation (i.e., the equation with the right-hands side set to 0) and an arbitrary solution of the inhomogeneous equation. The latter can be found by the ansatz that the system oscillates with the frequency of the force, $\omega_d$, which leads to the solution you quote.

Since all homogeneous solutions go to $0$ for $t \rightarrow 0$ (due to frictional energy loss), this particular solution of the inhomogeneous equations is the "stationary state" at long times, i.e., $x$ oscillates with the frequency of the driving force, $\omega_d$ but with a phase shift relative to the driving force.

The entire analysis is tremendously simplified by working with complex exponential functions, i.e., you write
$$\ddot{z} + \gamma \dot{z} + \omega_0^2 z=\frac{F}{m} \mathrm{i} \exp(-\mathrm{i} \omega_d t).$$
Then, at the end of the calculation you get as the solution of your equation $x=\mathrm{Re} z$.

 

[mpresic3]

It sounds to me like you want to start the system with zero initial conditions, ie velocity and position = 0. Int his case, you would have only the driven solution, but you cannot take damping gamma out of the systen as you did. Just because the damping at your t = 0 is zero, it does not mean it will remain zero, when the force sets the mass in motion. You need to look for a driven contribution that still may involve gamma

 

[vanhees71]

With 0 initial conditions you do not get only the stationary solution. To fullfill the boundary conditions of the full solution you need to add the free solutions with the appropriate integration constants chosen accordingly, i.e., initially you have a transient state. The free solutions decay with a typical lifetime $1/\gamma$.

 

[kelly0303]

No, it's not right. The final general solution is the superposition of the solution of the homogeneous equation (i.e., the equation with the right-hands side set to 0) and an arbitrary solution of the inhomogeneous equation. The latter can be found by the ansatz that the system oscillates with the frequency of the force, $\omega_d$, which leads to the solution you quote.

Since all homogeneous solutions go to $0$ for $t \rightarrow 0$ (due to frictional energy loss), this particular solution of the inhomogeneous equations is the "stationary state" at long times, i.e., $x$ oscillates with the frequency of the driving force, $\omega_d$ but with a phase shift relative to the driving force.

The entire analysis is tremendously simplified by working with complex exponential functions, i.e., you write
$$\ddot{z} + \gamma \dot{z} + \omega_0^2 z=\frac{F}{m} \mathrm{i} \exp(-\mathrm{i} \omega_d t).$$
Then, at the end of the calculation you get as the solution of your equation $x=\mathrm{Re} z$.

Sorry I am not sure I understand. In my case I can assume that $\gamma = 0$ (for context my system is an ion oscillating in vacuum and during the measurement time the amplitude doesn't change at all, at least not to an amount that can be measured experimentally). If I assume that the ion has initially zero velocity and starts at zero position then apply a driving force to take it to $A_d$ then stop the force, shouldn't the ion just oscillate as $x(t) = A_d\cos(\omega_0t)$? What else do I need to add to this? Thank you!

 

[kelly0303]

It sounds to me like you want to start the system with zero initial conditions, ie velocity and position = 0. Int his case, you would have only the driven solution, but you cannot take damping gamma out of the systen as you did. Just because the damping at your t = 0 is zero, it does not mean it will remain zero, when the force sets the mass in motion. You need to look for a driven contribution that still may involve gamma

But in my case $\gamma$ is zero (for context my system is an ion oscillating in vacuum and during the measurement time the amplitude doesn't change at all, at least not to an amount that can be measured experimentally).

 

[mpresic3]

If your gamma is zero, why have it in the problem to begin with. Cant you just deal with the driven undampled harmonic oscillator. The question is is there damping in the system or not. The context does not neccessarily tell is there is no damping. Are there energy losses due to radiation from the oscillating ion. I expect you want to neglect these losses, but it is not clear whether you have no damping in the system in which case, you can go back to your first equation and solve it without the gamma term, otherwise, you must consider gamma, even if you set it in oscillation with zero initial conditions.

 

[vanhees71]

Sorry I am not sure I understand. In my case I can assume that $\gamma = 0$ (for context my system is an ion oscillating in vacuum and during the measurement time the amplitude doesn't change at all, at least not to an amount that can be measured experimentally). If I assume that the ion has initially zero velocity and starts at zero position then apply a driving force to take it to $A_d$ then stop the force, shouldn't the ion just oscillate as $x(t) = A_d\cos(\omega_0t)$? What else do I need to add to this? Thank you!

For $\gamma=0$ the transient motion never dies out. If the driving force is oscillating with the eigenfrequency $\omega_0$ you even run into the resonance catastrophe. For your problem, where you have an arbitrary external force $m A(t)$ you can use the retarded Green's function to analyze the problem.

The Green's function is defined by the equation
$$\ddot{g}(t)+\omega_0^2 g(t)=\delta(t)$$
with the "causality condition" $g(t)=0$ for $t<0$. To solve this problem you have
$$g_<(t)=0, \quad g_>(t)=C_1 \cos(\omega_0 t) + C_2 \sin(\omega_0 t)],$$
where $g_<$ is the function for $t<0$ and $g_>$ for $t>0$.

To get $C_1$ and $C_2$ you need conditions at $t \rightarrow 0^+$. This we get from the fact that the 2nd derivative must lead to the $\delta$ distribution and thus $\dot{g}$ must have a jump at $t=0$ and $g$ must be continous at $t=0$. The latter condition leads to $g_>(0^+)=0$, i.e., $C_1=0$.

To determine $C_2$ integrate the equation of motion over an intervale $(-\epsilon,\epsilon)$ and then take $\epsilon \rightarrow 0^+$, giving
$$\omega_0 \dot{g}_{>}(0^+)=0 \; \Rightarrow\; C_2 \omega_0=1\; \Rightarrow C_2=\frac{1}{\omega_0},$$
and thus
$$g(t)=\frac{\Theta(t)}{\omega_0} \sin(\omega_0 t).$$
Then for $t>0$ the solution of your problem is given by
$$x(t)=x_0(t)+\int_0^{\infty} \mathrm{d} t' g(t-t') A(t')=x_0(t) + \int_0^t \mathrm{d} t' \frac{1}{\omega_0} \sin[\omega_0(t-t')] A(t').$$
$x_0$ is an arbitrary solution of the homogeneous equation.

Now
$$\dot{x}(t)=\dot{x}_0(t) + \frac{1}{\omega_0} \int_0^t \mathrm{d} t' \cos[\omega(t-t') A(t')].$$
With your boundary conditions $x(0)=0$ and $\dot{x}(0)=0$ you get $x_0(t) \equiv 0$, i.e.,
$$x(t)=\int_0^t \mathrm{d} t' \frac{1}{\omega_0} \sin[\omega_0(t-t')] A(t').$$
If $A(t)=0$ for $t>T$ of course you get what you expect, because then for $t>T$ you have
$$x(t)=\int_0^T \mathrm{d} t' \sin[\omega_0(t-t')] A(t') = \int_0^T \mathrm{d} t' \frac{1}{\omega_0} \left \{\sin(\omega_0 t) \cos(\omega_0 t') -\cos(\omega_0 t) \sin(\omega' t)]A(t') \right \}=C_1 \cos(\omega_0 t) +C_2 \sin(\omega_0 t)$$
with
$$C_1=-\frac{1}{\omega_0} \int_0^T \mathrm{d}t' \sin(\omega' t) A(t'), \quad C_2=\frac{1}{\omega_0} \int_0^T \mathrm{d}t' \cos(\omega' t) A(t').$$

 

[sophiecentaur]

You need to look for a driven contribution that still may involve gamma

This was bothering me too. If you have a lossless oscillator and driving source, it will have infinite Q. Wouldn't that imply a zero response to any drive frequency except its natural resonance? So the oscillator would remain stationary.

 

[vanhees71]

No, as my solution shows you always have the force kicking the oscillator. If it's undamped, i.e., $\gamma=0$ and you use a harmonic force with the eigenfrequency, $\omega_0$, of the oscillator the amplitude will grow linearly with time. That's the socalled "resonance catastrophe". The damped case is better behaved. It always leads for $t \rightarrow \infty$ to harmonic motion with the imposed frequency of the harmonic driving force, $\omega_d$. The eigenoscillations are damped out exponentially.

 

[kelly0303]

For $\gamma=0$ the transient motion never dies out. If the driving force is oscillating with the eigenfrequency $\omega_0$ you even run into the resonance catastrophe. For your problem, where you have an arbitrary external force $m A(t)$ you can use the retarded Green's function to analyze the problem.

The Green's function is defined by the equation
$$\ddot{g}(t)+\omega_0^2 g(t)=\delta(t)$$
with the "causality condition" $g(t)=0$ for $t<0$. To solve this problem you have
$$g_<(t)=0, \quad g_>(t)=C_1 \cos(\omega_0 t) + C_2 \sin(\omega_0 t)],$$
where $g_<$ is the function for $t<0$ and $g_>$ for $t>0$.

To get $C_1$ and $C_2$ you need conditions at $t \rightarrow 0^+$. This we get from the fact that the 2nd derivative must lead to the $\delta$ distribution and thus $\dot{g}$ must have a jump at $t=0$ and $g$ must be continous at $t=0$. The latter condition leads to $g_>(0^+)=0$, i.e., $C_1=0$.

To determine $C_2$ integrate the equation of motion over an intervale $(-\epsilon,\epsilon)$ and then take $\epsilon \rightarrow 0^+$, giving
$$\omega_0 \dot{g}_{>}(0^+)=0 \; \Rightarrow\; C_2 \omega_0=1\; \Rightarrow C_2=\frac{1}{\omega_0},$$
and thus
$$g(t)=\frac{\Theta(t)}{\omega_0} \sin(\omega_0 t).$$
Then for $t>0$ the solution of your problem is given by
$$x(t)=x_0(t)+\int_0^{\infty} \mathrm{d} t' g(t-t') A(t')=x_0(t) + \int_0^t \mathrm{d} t' \frac{1}{\omega_0} \sin[\omega_0(t-t')] A(t').$$
$x_0$ is an arbitrary solution of the homogeneous equation.

Now
$$\dot{x}(t)=\dot{x}_0(t) + \frac{1}{\omega_0} \int_0^t \mathrm{d} t' \cos[\omega(t-t') A(t')].$$
With your boundary conditions $x(0)=0$ and $\dot{x}(0)=0$ you get $x_0(t) \equiv 0$, i.e.,
$$x(t)=\int_0^t \mathrm{d} t' \frac{1}{\omega_0} \sin[\omega_0(t-t')] A(t').$$
If $A(t)=0$ for $t>T$ of course you get what you expect, because then for $t>T$ you have
$$x(t)=\int_0^T \mathrm{d} t' \sin[\omega_0(t-t')] A(t') = \int_0^T \mathrm{d} t' \frac{1}{\omega_0} \left \{\sin(\omega_0 t) \cos(\omega_0 t') -\cos(\omega_0 t) \sin(\omega' t)]A(t') \right \}=C_1 \cos(\omega_0 t) +C_2 \sin(\omega_0 t)$$
with
$$C_1=-\frac{1}{\omega_0} \int_0^T \mathrm{d}t' \sin(\omega' t) A(t'), \quad C_2=\frac{1}{\omega_0} \int_0^T \mathrm{d}t' \cos(\omega' t) A(t').$$

Thanks for the details, I will go through them carefully, but I might miss-explain my problem (and overcomplicate it). For concreteness, in my setup I have an ion at the center of a Penning trap (I can provide detailed references about this kind of setup if needed). It is common practice (for various purposes), to excite the axial motion of the ion (if we ignore thermal noise, we can assume that the ion is at the center of the trap with zero velocity initially) to a given amplitude then let it oscillate. One can damp the ion motion by coupling the Penning trap to an RLC circuit, such that the induced image current due to the ion motion thermalizes the ion very fast. However, if one changes the RLC resonant frequency enough, the ion will not be damped basically over the time scale of the measurement (ms to seconds) and this second case is where I am at. In all the papers I read about this, they simply assume that the ion motion (after this initial kick which can be done in various ways, including applying a resonant force) is simply oscillating at the resonant frequency i.e. $A_d\cos(\omega_0 t)$. I am not sure where the extra math comes in, as it seems like I am getting in the end the same behaviour as they are getting. I am not sure what am I missing (am I oversimplifing the problem?)

Edit: I am attaching below a figure I have from a PhD thesis (I will try to find the actual thesis). By applying a dipolar excitation they are able to bring the ion to a given axial amplitude of choice (in this case they have thermal noise, so the initial phase affects the final amplitude, but this is exactly what I want to estimate for my setup, too). The point is that once they reach the amplitude they want, they turn off the excitation and the ion simply oscillate as I mentioned above.

 

[sophiecentaur]

So the oscillator would remain stationary.

@vanhees71 Oh yes. That's wrong, the oscillator would move according to its mass and the applied force but there would be a phase delay, according to the sign of the reactance of the oscillator. That would even apply to an oscillator with a finite damping. Well off resonance you just get a phase delay.

There again, the phase shift would still depend on the source impedance of the driver except for a 'voltage source or current source' or equivalent. Much easier when the damping R of the oscillator is much higher than the source impedance.

The situation when the applied waveform is turned on or off, instantly makes it harder because there are components of the square wave modulated driver that will coincide with the natural resonance. But zero energy when there's no bandwidth.

 

[sophiecentaur]

By applying a dipolar excitation they are able to bring the ion to a given axial amplitude of choice

and @vanhees71 I notice that the above graph shows / assumes that the amplitude of the oscillations grows in equal steps. The model assumes the amplitude of the sinusoidal applied force is constant. How to apply a familiar electrical analogue to it? The force from the applied field would have constant amplitude, independent of the motion of the ion (I assume). Thing is, I'm looking for a real damping effect to produce an inverse exponential approach to a real maximum. What would actually be supplying this damping?

 

[vanhees71]

It's hard to say, what exactly is written in the PhD thesis. If it's final and the candidate has defended it, usually it's made public. So you should be able to find it in the internet.

What I could imagine is that the plot shown in #11 is indeed describing the resonance case for the (nearly?) undamped oscillator, because then the amplitude rises linear with time.

We need a particular solution or the inhomogeneous equation for this case. Of course you could use my method with the Green's function, but it's maybe more illuminating to do it from scratch. What does not work, is the "ansatz of the kind of the right-hand side", because you'll not find a solution for the amplitude. That's because of the resonance catastrophe. We use the method with the complexifed right-hand side. At the end we take the imaginary part to solve for the original real equation in #1.

The correct ansatz to solve
$$\ddot{x}+\omega_0^2 x=A \exp(-\mathrm{i} \omega_0 t)$$
is a bit more general than the naive one,
$$x(t)=y(t) \exp(-\mathrm{i} \omega_0 t).$$
Plugging this into the equation gives
$$-2 \mathrm{i} \omega_0 \dot{y} + \ddot{y}=A.$$
We only need one solution, not the general one. That's why we can set $y(t)=B t$, which gives
$$-2 \mathrm{i} \omega_0 B=A \; \Rightarrow \; B=\frac{\mathrm{i} A}{2 \omega_0}.$$
So our particular solution is
$$x(t)=\frac{\mathrm{i} A}{2 \omega_0} t \exp(-\mathrm{i} \omega t).$$
Now taking the imaginary part gives
$$x_{\text{phys}}(t)=-\mathrm{Im} x(t)=-\frac{A}{2 \omega_0} t \cos(\omega t).$$
The general reals solution now is found by adding the general solution of the homogeneous equation,
$$x_{\text{phys}}(t) = \left (-\frac{A}{2 \omega_0} t+C_1 \right) \cos(\omega_0 t) + C_2 \sin(\omega_0 t).$$
For the initial conditions $x_{\text{phys}}(0)=0$ and $\dot{x}_{\text{phys}}(0)=0$ you find $C_1=0$ and $C_2=A/(2 om0^2)$, i.e.,
$$x_{\text{phys}}(t)=\frac{A}{2 \omega_0^2} [\sin(\omega_0 t)-\omega_0 t \cos(\omega_0 t)].$$
For $t \gg 1/\omega$ the 2nd term is dominating, and the asymptotic behavior is thus a harmonic oscillation with linearly growing amplitude.

 

[kelly0303]

It's hard to say, what exactly is written in the PhD thesis. If it's final and the candidate has defended it, usually it's made public. So you should be able to find it in the internet.

What I could imagine is that the plot shown in #11 is indeed describing the resonance case for the (nearly?) undamped oscillator, because then the amplitude rises linear with time.

We need a particular solution or the inhomogeneous equation for this case. Of course you could use my method with the Green's function, but it's maybe more illuminating to do it from scratch. What does not work, is the "ansatz of the kind of the right-hand side", because you'll not find a solution for the amplitude. That's because of the resonance catastrophe. We use the method with the complexifed right-hand side. At the end we take the imaginary part to solve for the original real equation in #1.

The correct ansatz to solve
$$\ddot{x}+\omega_0^2 x=A \exp(-\mathrm{i} \omega_0 t)$$
is a bit more general than the naive one,
$$x(t)=y(t) \exp(-\mathrm{i} \omega_0 t).$$
Plugging this into the equation gives
$$-2 \mathrm{i} \omega_0 \dot{y} + \ddot{y}=A.$$
We only need one solution, not the general one. That's why we can set $y(t)=B t$, which gives
$$-2 \mathrm{i} \omega_0 B=A \; \Rightarrow \; B=\frac{\mathrm{i} A}{2 \omega_0}.$$
So our particular solution is
$$x(t)=\frac{\mathrm{i} A}{2 \omega_0} t \exp(-\mathrm{i} \omega t).$$
Now taking the imaginary part gives
$$x_{\text{phys}}(t)=-\mathrm{Im} x(t)=-\frac{A}{2 \omega_0} t \cos(\omega t).$$
The general reals solution now is found by adding the general solution of the homogeneous equation,
$$x_{\text{phys}}(t) = \left (-\frac{A}{2 \omega_0} t+C_1 \right) \cos(\omega_0 t) + C_2 \sin(\omega_0 t).$$
For the initial conditions $x_{\text{phys}}(0)=0$ and $\dot{x}_{\text{phys}}(0)=0$ you find $C_1=0$ and $C_2=A/(2 om0^2)$, i.e.,
$$x_{\text{phys}}(t)=\frac{A}{2 \omega_0^2} [\sin(\omega_0 t)-\omega_0 t \cos(\omega_0 t)].$$
For $t \gg 1/\omega$ the 2nd term is dominating, and the asymptotic behavior is thus a harmonic oscillation with linearly growing amplitude.

Thanks a lot, that formula makes more sense (in terms of not getting infinite amplitude without damping). So for my case, I need first to know for how long to apply that force such that I can reach a given amplitude $A_d$ (and how would the system behave after that). At that point I don't want to increase the amplitude anymore and just let the ion evolve without any external intervention. So I need at a given time $t_0$ the amplitude to be $A_d$. If we assume that $t_0 = 2\pi / \omega_0$ (which can be done by adjusting the value of $A$ experimentally as needed) we get that $A_d = x(t_0) = -\frac{A}{2\omega_0^2}\frac{2\pi\omega_0}{\omega_0}$, so we get that we need $A_d = \frac{A\pi}{\omega_0^2}$. We also get, as desired, that $\dot{x}(t_0)=0$. At this point I turn off the external force, so if I redefine this $t_0$ as my starting point, I should simply get that the ion moves as $A_d\cos{(\omega_0 t)}$. Is this not right? But this basically gives what I got in my initial post, no? Am I still missing something in my derivation? Thank you!

 

[vanhees71]

For that you just use the Green's function. It's not so easy to solve the equation directly for a force that's not harmonic.

 

[kelly0303]

For that you just use the Green's function. It's not so easy to solve the equation directly for a force that's not harmonic.

But what is wrong with my approach above? The formula you provided is true for a continuous harmonic force. So at a given time $t_0$ the amplitude is given by that formula (the formula doesn't know whether I will turn off the force or not after that). But if I turn it off, the system should evolve as a free harmonic oscillator and in the end I get the result from above. I really don't see what step above is wrong? Also this behavior based on my result is basically what they observe in measurements in Penning traps (simple, undamped harmonic oscillations).

I also looked more carefully at the derivation from here. Equation (25) is exactly what I got in my first post. Why can't I use that in my case?

I guess I don't understand why do I need Green functions in order to get what I need, given that my result appears consistent with experimental observations (and some other derivations I found online). I assume that using Green functions is a more general approach but I don't see what is the difference, in my case, compared to my approach.

 

[vanhees71]

Obviously Eq. (25) doesn't work for $\gamma=0$ and $\omega_d=\omega_0$. That's why you have to treat this case separately as in #14.

If the harmonic motion acts only for a finite time, you need to use the Green's-function formula in #8.