B. The expression under the square root in Eq. (92), ( M 11- M 22)2 + 4 M 12 M 21, is negative. In this case, the square root is imaginary, making the real parts of both roots equal, Re = ( M 11 + M 22)/2, and their imaginary parts equal but opposite. As a result, here there can be just two types of fixed points: (i) Stable focus, at ( M 11 + M 22) < 0. The phase plane trajectories are spirals going to the origin (i.e. toward the fixed point) – see Fig. 8c with the solid arrow.
(ii) Unstable focus, taking place at ( M 11 + M 22) > 0, differs from the stable one only by the direction of motion along the phase trajectories – see the dashed arrow in the same Fig. 8c.
C. Frequently, the border case, M 11 + M 22 = 0, corresponding to the orbital (“indifferent”) stability already discussed in Sec. 3.2, is also distinguished, and the corresponding fixed point is referred to as the center (Fig. 8d). Considering centers as a separate category makes sense because such fixed points are typical for Hamiltonian systems, whose first integral of motion may be frequently represented as the distance of the representing point from a certain center. For example, introducing new variables
~
~
~
q q
and
~
q q
m , we may rewrite Eq. (3.12) of a harmonic oscillator without dissipation (again, 1
2
1
with indices “ef” dropped for brevity), as a system of two first-order differential equations: 1
~
~
~
~
q
q ,
q
q
,
(5.93)
1
2
2
1
m
i.e. as a particular case of Eq. (87), with M
2
11 = M 22 = 0, and M 12 M 21 = –/ m –0 < 0, and hence ( M 11-M
2
~ ~
22)2 + 4 M 12 M 21 = –40 < 0, and M 11 + M 22 = 0. On the symmetrized phase plane q , q / Z , where 1
2
the parameter Z ( m)1/2 m0 is the oscillator’s impedance, the sinusoidal oscillations of amplitude A are represented by a circle of radius A about the center-type fixed point A = 0. In the case when q~ q~
1
is the linear coordinate q of an actual mechanical oscillator, so that ~
~
q q
m is its linear momentum
2
1
p q
m , such a circular trajectory corresponds to the conservation of the oscillator’s energy 2
2
~ 2
2
p
q
q
A
~2
2
E T U
q
const .
(5.94)
2 m
2
2
1
Z
2
This is a convenient moment for a brief discussion of the so-called Poincaré (or “slow-variable”, or “stroboscopic”) plane.35 From the point of view of the basic Eq. (41), the sinusoidal oscillations q( t) 34 The term “saddle” is due to the fact that in this case, the system’s dynamics is qualitatively similar to that of a heavily damped motion in a 2D potential U( ~ ~
q , q ) having the shape of a horse saddle (or a mountain pass).
1
2
35Named after Jules Henri Poincaré (1854-1912), who is credited, among many other achievements in physics and mathematics, for his contributions to special relativity (see, e.g., EM Chapter 9), and the basic idea of unstable trajectories responsible for the deterministic chaos – to be discussed in Chapter 9 of this course.
Chapter 5
Page 24 of 38
CM: Classical Mechanics
= A cos( t – ), described by a circular trajectory on the actual (symmetrized) phase plane, correspond to a fixed point { A, }, which may be conveniently represented by a stationary geometric point on the plane with these polar coordinates – see Fig. 9a. (As follows from Eq. (4), the Cartesian coordinates of the point on that plane are just the variables u A cos and v A cos that were used, in particular, in the last section.) The quasi-sinusoidal process (41), with slowly changing A and , may be represented by slow motion of that point on this Poincaré plane.
(a)
(b)
p
m
v( t)
A( t)
v
( t)
Fig. 5.9. (a) Representation of a
0
u( t)
t
q
sinusoidal oscillation (point) and a
0
slow transient process (line) on the
A
Poincaré plane, and (b) the relation
t
between the “fast” phase plane and the
“slow” (Poincaré) plane.
u
Figure 9b shows a convenient way to visualize the relation between the actual phase plane of an oscillator, with the “fast” symmetrized coordinates q and p/ m, and the Poincaré plane with the “slow”
coordinates u and v: the latter plane rotates relative to the former one, about the origin, clockwise, with the angular velocity .36 Another, “stroboscopic” way to generate the Poincaré plane pattern is to have a fast glance at the “real” phase plane just once during the oscillation period T = 2/.
In many cases, the representation on the Poincaré plane is more convenient than that on the
“real” phase plane. In particular, we have already seen that the reduced equations for such important phenomena as the phase locking and the parametric oscillations, whose original differential equations include time explicitly, are time-independent – cf., e.g., (75) and (79) describing the latter effect. This simplification brings the equations into the category considered earlier in this section, and enables an easy classification of their fixed points, which may shed additional light on their dynamic properties.
In particular, Fig. 10 shows the classification of the only (trivial) fixed point A 1 = 0 on the Poincaré plane of the parametric oscillator, which follows from Eq. (83). As the parameter modulation depth is increased, the type of this fixed point changes from a stable focus (pertinent to a simple oscillator with damping) to a stable node and then to a saddle describing the parametric excitation. In the last case, the two directions of the perturbation growth, so prominently featured in Fig. 8b, correspond to the two possible values of the oscillation phase , with the phase choice determined by initial conditions.
This double degeneracy of the parametric oscillation’s phase could already be noticed from Eqs.
(77), because they are evidently invariant with respect to the replacement + . Moreover, the degeneracy is not an artifact of the van der Pol approximation, because the initial equation (75) is already invariant with respect to the corresponding replacement q( t) q( t – /). This invariance 36 This notion of phase plane rotation is the origin of the term “Rotating Wave Approximation”, mentioned above.
(The word “wave” is an artifact of this method’s wide application in classical and quantum optics.) Chapter 5
Page 25 of 38
CM: Classical Mechanics
means that all other characteristics (including the amplitude) of the parametric oscillations excited with either of the two phases are exactly similar. At the dawn of the computer age (in the late 1950s and early 1960s), there were substantial attempts, especially in Japan, to use this property for storage and processing digital information coded in the binary-phase form. Though these attempts have not survived the competition with simpler approaches based on binary-voltage coding, some current trends in the development of prospective reversible and quantum computers may be traced back to that idea.
stable
4
nodes
saddles
1 /
4
2
2
δ ξ 2
4
stable
stable
focuses
focuses
Fig. 5.10. Types of the trivial fixed
point of a parametric oscillator.
0
5.7. Numerical approaches
If the amplitude of oscillations, for whatever reason, becomes so large that nonlinear terms in the equation describing an oscillator become comparable with its linear terms, numerical methods are virtually the only avenue available for their theoretical studies. In Hamiltonian 1D systems, such methods may be applied directly to Eq. (3.26), but dissipative and/or parametric systems typically lack such first integrals of motion, so that the initial differential equation has to be solved.
Let us discuss the general idea of such methods on the example of what mathematicians call the Cauchy problem (finding the solution for all moments of time, starting from the known initial conditions) for the first-order differential equation
q f ( t, q).
(5.95)
(The generalization to a system of several such equations is straightforward.) Breaking the time axis into small equal steps h (Fig. 11) we can reduce the equation integration problem to finding the function’s value at the next time point, qn+1 q( tn+1) q( tn + h) from the previously found value qn = q( tn) – and, if necessary, the values of q at other previous time steps.
q
n 1
q
qn
h
Fig. 5.11. The basic notions used at numerical
h / 2
integration of ordinary differential equations.
t
t
n
n 1
t
In the simplest approach (called the Euler method), qn+1 is found using the following formula: q
q k,
n 1
n
(5.96)
k h f ( t , q ).
n
n
Chapter 5
Page 26 of 38
CM: Classical Mechanics
This approximation is equivalent to the replacement of the genuine function q( t), on the segment [ tn, tn+1], with the two first terms of its Taylor expansion in point tn:
q( t h) q( t ) q( t ) h q( t ) hf ( t , q ).
(5.97)
n
n
n
n
n
n
This approximation has an error proportional to h 2. One could argue that by making the step h sufficiently small, the Euler method’s error might be made arbitrarily small, but even with all the number-crunching power of modern computer platforms, the CPU time necessary to reach sufficient accuracy may be too large for big problems.37 Besides that, the increase of the number of time steps, which is necessary at h 0 at a fixed total time interval, increases the total rounding errors and eventually may cause an increase, rather than the reduction of the overall error of the computed result.
A more efficient way is to modify Eq. (96) to include the terms of the second order in h. There are several ways to do this, for example using the 2nd-order Runge-Kutta method: q
q k ,
n 1
n
2
h
k
(5.98)
k h f t ,
1
q
,
k h f ( t , q ).
2
n
2 n
2
1
n
n
One can readily check that this method gives the exact result if the function q( t) is a quadratic polynomial, and hence in the general case its errors are of the order of h 3. We see that the main idea here is to first break the segment [ tn, tn+1] in half (see Fig. 11 again), evaluate the right-hand side of the differential equation (95) at the point intermediate (in both t and q) between the points number n and ( n
+ 1), and then use this information to evaluate qn+1.
The advantage of the Runge-Kutta approach over other second-order methods is that it may be readily extended to the 4th order, without an additional breakup of the interval [ tn, tn+1]: 1
q
q ( k 2 k 2 k k ),
n 1
n
1
2
3
4
6
(5.99)
h
k
h
k
k h f ( t ,
h q k ), k h f t ,
2
q
,
k h f t ,
1
q
,
k h f ( t , q ).
4
n
n
3
3
n
2 n
2
2
n
2 n
2
1
n
n
This method has a much lower error, O( h 5), without being not too cumbersome. These features have made the 4th-order Runge-Kutta the default method in most numerical libraries. Its extension to higher orders is possible, but requires more complex formulas, and is justified only for some special cases, e.g., very abrupt functions q( t).38 The most frequent enhancement of the method is an automatic adjustment of the step h to reach the pre-specified accuracy, but not make more calculations than necessary.
Figure 12 shows a typical example of an application of that method to the very simple problem of a damped linear oscillator, for two values of fixed time step h (expressed in terms of the number N of such steps per oscillation period). The black straight lines connect the adjacent points obtained by the 37 In addition, the Euler method is not time-reversible – the handicap that may be essential for the integration of Hamiltonian systems described by systems of second-order differential equations. However, this drawback may be partly overcome by the so-called leapfrogging – the overlap of time steps h for a generalized coordinate and the corresponding generalized velocity.
38 The most popular approaches in such cases are the Richardson extrapolation, the Bulirsch-Stoer algorithm, and a set of so-called prediction-correction techniques, e.g. the Adams-Bashforth-Moulton method – see the literature recommended in MA Sec. 16(iii).
Chapter 5
Page 27 of 38
CM: Classical Mechanics
4th-order Runge-Kutta method, while the points connected with the green straight lines represent the exact analytical solution (22). The plots show that a-few-percent errors start to appear only at as few as
~10 time steps per period, so that the method is indeed very efficient.
Let me hope that the discussion in the next section will make the conveniences and the handicaps of the numerical approach to problems of nonlinear dynamics very clear.
(a)
(b)
1
1
N 30
N 6
0.5
0.5
q( t)
k
0
0
0.5
0.5
1
1
0
10
20
30
40
50
60
0
10
20
30
40
50
60
t
t
0
0
Fig. 5.12. Results of the Runge-Kutta solution of Eq. (6) (with /0 = 0.03) for: (a) 30 and (b) 6 points per oscillation period. The results are shown by
points; the black and green lines are only the guides for the eye.
5.8. Higher-harmonic and subharmonic oscillations
Figure 13 shows the numerically calculated39 transient process and stationary oscillations in a linear oscillator and a very representative nonlinear system, the pendulum described by Eq. (42), both with the same 0. Both systems are driven by a sinusoidal external force of the same amplitude and frequency – in this illustration, equal to the own small-oscillation frequency 0 of both systems. The plots show that despite a very substantial amplitude of the pendulum oscillations (the angle amplitude of about one radian), their waveform remains almost exactly sinusoidal.40 On the other hand, the nonlinearity affects the oscillation amplitude very substantially. These results imply that the corresponding reduced equations (60), which are based on the assumption (41), may work very well far beyond its formal restriction q << 1.
Still, the waveform of oscillations in a nonlinear system always differs from that of the applied force – in our case, from the sine function of frequency . This fact is frequently formulated as the generation, by the system, of higher harmonics. Indeed, the Fourier theorem tells us that any non-sinusoidal periodic function of time may be represented as a sum of its basic harmonic of frequency , plus higher harmonics with frequencies n, with integer n > 1.
Note that an effective generation of higher harmonics is only possible with adequate nonlinearity of the system. For example, consider the nonlinear term q 3 used in the equations explored in Secs. 2
39 All numerical results shown in this section have been obtained by the 4th-order Runge-Kutta method with the automatic step adjustment that guarantees the relative error of the order of 10-4 – much smaller than the pixel size in the shown plots.
40 In this particular case, the higher harmonic content is about 0.5%, dominated by the 3rd harmonic, whose amplitude and phase are in a very good agreement with Eq. (50).
Chapter 5
Page 28 of 38
CM: Classical Mechanics
and 3. If the waveform q( t) is sinusoidal, such term will have only the basic (1st) and the 3rd harmonics –
see, e.g., Eq. (50). As another example, the “pendulum nonlinearity” sin q cannot produce, without a time-independent component (”bias”) in q( t), any even harmonic, including the 2nd one. The most efficient generation of harmonics may be achieved using systems with the sharpest nonlinearities – e.g., semiconductor diodes whose current may follow an exponential dependence on the applied voltage through several orders of magnitude. 41
2
2
1
1
q( t) 0
0
1
1
2
2
0
5
10
15
20
28
29
30
2
2
1
1
q( t)
0
0
1
1
2
2
0
5
10
15
20
28
29
30
t /
2
t /
2
0
0
Fig. 5.13. The oscillations induced by a similar sinusoidal external force (turned on at t = 0) in two systems with the same small-oscillation frequency ω 0 and low damping: a linear oscillator (two top panels) and a pendulum (two bottom panels). In all cases, δ/ ω 0 = 0.03, f 0 = 0.1, and ω = ω 0.
Another way to increase the contents of an n th higher harmonic in a nonlinear oscillator is to reduce the excitation frequency to ~0/ n, so that the oscillator resonated at the frequency n 0 of the desired harmonic. For example, Fig. 14a shows the oscillations in a pendulum described by the same Eq. (42), but driven at frequency = 0/3. One can see that the 3rd harmonic amplitude may be comparable with that of the basic harmonic, especially if the external frequency is additionally lowered (Fig. 14b) to accommodate for the deviation of the effective frequency ω 0( A) of own oscillations from its small-oscillation value ω 0 – see Eq. (49), Fig. 4, and their discussion in Sec. 2 above.
However, numerical modeling of nonlinear oscillators, as well as experiments with their physical implementations, bring more surprises. For example, the bottom panels of Fig. 15 show oscillations in a pendulum under the effect of a strong sinusoidal force with a frequency close to 30. One can see that at some parameter values and initial conditions, the system’s oscillation spectrum is heavily contributed (almost dominated) by the 3rd sub harmonic, i.e. the Fourier component of frequency /3 0.
This counter-intuitive phenomenon of such subharmonic generation may be explained as follows. Let us assume that subharmonic oscillations of frequency / 3 0 have somehow appeared, and coexist with the forced oscillations of frequency 3:
41 This method is used in practice, for example, for the generation of electromagnetic waves with frequencies in the terahertz range (1012-1013 Hz), which is still in wait for efficient electronic self-oscillators.
Chapter 5
Page 29 of 38
CM: Classical Mechanics
q( t)
t
A cos Ψ A cos Ψ ,
Ψ
where
t , Ψ
. (5.100)
sub
sub
sub
3
sub
Then the leading nonlinear term, q 3, of the Taylor expansion of the pendulum’s nonlinearity sin q, is proportional to
3
q ( A cos Ψ A cos Ψ )3
sub
sub
(5.101)
3
A cos3 Ψ 3 2
A A cos2 Ψ cos Ψ
3
2
AA cos Ψ cos2 Ψ
3
A cos3 Ψ .
sub
sub
ub
s
sub
sub
sub
2
1
q( t) 0
1
2
0
5
10
15
20
25
30
2
1
0
q( t)
1
2 0
5
10
15
20
25
30
t /
2
0
Fig. 5.14. The oscillations induced in a pendulum, with damping δ/ ω
0 = 0.03, by a sinusoidal external
force of amplitude f 0 = 0.75, and frequencies ω 0/3 (top panel) and 0.8 ω 0/3 (bottom panel).
2
2
1
1
q( t) 0
0
1
1
2
2
0
5
10
15
20
28
29
30
2
2
1
1
q( t) 0
0
1
1
2
2
0
5
10
15
20
28
29
30
t /
2
t /
2
0
0
Fig. 5.15. The oscillations of a pendulum with δ/ ω 0 = 0.03, driven by a sinusoidal external force of amplitude f 0 = 3 and frequency 0.83 ω 0, at initial conditions q(0) = 0 (the top panels) and q(0) = 1
(the bottom panels), with dq/ dt (0) = 0 in both cases.
Chapter 5
Page 30 of 38
CM: Classical Mechanics
While the first and the last terms of the last expression depend only on the amplitudes of the individual components of oscillations, the two middle terms are more interesting, because they produce so-called combinational frequencies of the two components. For our case, the third term, 2
2
3
3 A A cos Ψ cos Ψ
2
AA cos(Ψ Ψ
2
) ... ,
(5.102)
ub
s
sub
4
ub
s
sub
is of special importance, because it produces (besides other combinational frequencies) the subharmonic component with the total phase
t
2
2 .
(5.103)
sub
sub
3
Thus, this nonlinear contribution is synchronous with the subharmonic oscillations, and describes the interaction that can, within a certain range of the mutual phase shift between the Fourier components, deliver to them energy from the external force, so that the oscillations may be sustained. Note, however, that the amplitude of the term describing this energy exchange is proportional to the square of A sub, and vanishes at the linearization of the equations of motion near the trivial fixed point. This means that the point is always stable, i.e., the 3rd subharmonic cannot be self-excited and always needs an initial “kick-off” – compare the two panels of Fig. 15. The same is true for higher-order subharmonics.
Only the second subharmonic is a special case. Indeed, let us make a calculation similar to Eq.
(102), by replacing Eq. (101) with
q( t)
t
A cos Ψ A cos Ψ ,
Ψ
where
t , Ψ
,
(5.104)
sub
sub
sub
2
sub
for a nonlinear term proportional to q 2:
2
q ( A cos Ψ A cos Ψ )2
2
A cos2 Ψ 2 AA cosΨ cosΨ
2
A cos2 Ψ .
(5.105)
sub
sub
sub
sub
sub
sub
Here the combinational-frequency term capable of supporting the 2nd subharmonic,
2 AA cos Ψ cos Ψ
AA cos
AA
t
,
(5.106)
sub
sub
sub
Ψ Ψsub
cos
sub
sub
...
is linear in the subharmonic’s amplitude, i.e. survives the linearization near the trivial fixed point. This means that the second subharmonic may arise spontaneously, from infinitesimal fluctuations.
Moreover, such excitation of the second subharmonic is very similar to the parametric excitation that was discussed in detail in Sec. 5, and this similarity is not coincidental. Indeed, let us redo the expansion (106) making a somewhat different assumption – that the oscillations are a sum of the forced oscillations at the external force’s frequency and an arbitrary but weak perturbation: q t
( ) A cos( t q~
)
t
( ),
q~
with
A.
(5.107)
Then, neglecting the small term proportional to ~ 2
q , we get
2
2
~
q A cos2 ( t ) 2 q ( t) A cos( t ).
(5.108)
Besides the inconsequential phase shift , the second term in the last formula is exactly similar to the term describing the parametric effects in Eq. (75). This fact means that for a weak perturbation, a system with a quadratic nonlinearity in the presence of a strong “pumping” signal of frequency is equivalent to a system with parameters changing in time with frequency . This fact is broadly used for the Chapter 5
Page 31 of 38
CM: Classical Mechanics
parametric excitation at high (e.g., optical) frequencies, where the mechanical means of parameter modulation (see, e.g., Fig. 5) are not practicable. The necessary quadratic nonlinearity at optical frequencies may be provided by a non-centrosymmetric nonlinear crystal, e.g., the -phase barium borate (BaB2O4).
Before finishing this section, let me elaborate a bit on a general topic: the relation between the numerical and analytical approaches to problems of dynamics – and to physics as a whole. We have just seen that sometimes numerical solutions, like those shown in Fig. 15b, may give vital clues for previously unanticipated phenomena such as the excitation of subharmonics. (The phenomenon of deterministic chaos, which will be discussed in Chapter 9 below, presents another example of such
“numerical discoveries”.) One might also argue that for problems without exact analytical solutions, the numerical simulation may be an equally productive theoretical tool. These hopes are, however, muted by the general problem that is frequently called the curse of dimensionality,42 in which the last word refers to the number of parameters of the problem to be solved.43
Indeed, let us have one more look at Fig. 15. OK, we have been lucky to find a new phenomenon, the 3rd subharmonic generation, for a particular set of parameters – in that case, five of them: /0 = 0.03, /0 = 2.4, f 0 = 3, q(0) = 1, and dq/ dt (0) = 0. Could we tell anything about how common this effect is? Are subharmonics with different n possible in this system? The only way to address these questions computationally is to carry out similar numerical simulations at many points of the d-dimensional (in this case, d = 5) space of parameters. Say, we have decided that breaking the reasonable range of each parameter to N = 100 points is sufficient. (For many problems, even more points are necessary – see, e.g., Sec. 9.1.) Then the total number of numerical experiments to carry out is Nd = (102)5 = 1010 – not a simple task even for the powerful modern computing facilities. (Besides the pure number of required CPU cycles, consider the storage and analysis of the results.) For many important problems of nonlinear dynamics, e.g., turbulence, the parameter dimensionality d is substantially larger, and the computer resources necessary even for one numerical experiment, are much greater.
In view of the curse of dimensionality, approximate analytical considerations, like those outlined above for the subharmonic excitation, are invaluable. More generally, physics used to stand on two legs: experiment and analytical theory. The enormous progress of computer performance during the few last decades has provided it with one more support point (a tail? :-) – numerical simulation. This does not mean we can afford to discard any of the legs we are standing on.
5.9. Relaxation oscillations
Such synthesis of the analytical and numerical approaches is also beneficial for the discussion of the last subject of this chapter: nonlinear oscillators with high damping. Perhaps the most interesting effect in such systems is the so-called relaxation oscillations, a type of self-oscillations with highly non-sinusoidal waveforms. Let me demonstrate them using our old friend, Eq. (62) with < 0, whose 42 This term had been coined in 1957 by Richard Bellman in the context of the optimal control theory (where the dimensionality means the number of parameters affecting the system under control), but gradually has spread all over quantitative sciences using numerical methods.
43 In EM Sec. 1.2, I discuss the implications of this “curse” for a different case, when both analytical and numerical solutions to the same problem are possible.
Chapter 5
Page 32 of 38
CM: Classical Mechanics
properties at << 0 were discussed
2
2
in Sec. 4, because it will enable us to
D 2
.
0
follow the crossover from the harmonic
1
oscillations to the relaxation ones.
1
x
x
Figure 16 shows the results of the
numerical solution of this equation for
0
0
three characteristic values of its only
substantial parameter44
1
1
2
D
0 . (5.109)
2
2
0
2
1
0
1
2
0
0.5
1
1.5
y
/
2
(Indeed, if we introduce the natural
2
2
dimensionless time
D 0
.
2
0 t, displacement
x q/ q 0, where q 0 is the scale defined in
1
Eq. (64), and velocity y dx/ d, then the
1
second-order differential equation (62) x
x
may be rewritten as the following system
0
0
of two first-order equations: 45
dx y,
1
1
d
(5.110)
dy D1 2
y y x,
2
2
2
1
0
1
2
0
0.5
1
1.5
d
y
/
2
with
1
1
D being its only parameter.) The
D 20
left panels show phase planes [ x, dx/d]
of the oscillator, with their axes
0.5
0.5
x
x
swapped46 for the comparison with the
right panels showing the displacement x D
D
0
0
y 0 x
as a function of time.
If the damping is low (top two 0.5
0.5
panels), the system, launched from any
initial state, gradually approaches the
“limit
cycle”
of
nearly-harmonic 1
1
2
1
0
1
2
0
2
4
6
8
oscillations. Note that even for this, not
y
/
2
extremely
small
value
D = 0.2,
Fig. 5.16. The phase plane and time evolution of the self-
deviations of the waveform x() from a oscillator described by Eqs. (62) and (110), for three values of the purely sinusoidal function of time are normalized damping (109). The red and blue lines show the very small, its period (in the normalized system’s dynamics for two representative initial conditions, while time ) is very close to the small-the black lines, its asymptotic behavior (the “limit cycles”).
44 As Eq. (11) shows, for positive damping, this parameter is just the reciprocal Q-factor.
45 A somewhat different equation used in 1926 by B. van der Pol to trace the harmonic-to-relaxation-oscillation crossover for the first time, may be also reduced to Eq. (110), using the so-called Liénard's transformation.
46 Note that while on the usual phase plane, the free-oscillation process corresponds to a clockwise rotation of the representation point (see, e.g., Fig. 9), the axes’ swap in Fig. 16 makes the rotation counterclockwise.
Chapter 5
Page 33 of 38
Essential Graduate Physics
CM: Classical Mechanics
oscillation value 2, and its amplitude is also very close to the value 2/3 1.15 predicted by the van der Pol method – see Eq. (64).
As the damping is increased to D = 2 (middle panels), the limit cycle’s deviations from the circle, and hence the deviations of the waveform x() from a sinusoidal function become obvious. Note also that while the oscillation period becomes somewhat longer than its small-oscillation value, the transient processes of approaching the limit cycle become faster.
The trend of these changes becomes evident on the bottom panels, showing case D = 20. (The further increase of the damping does not change the results noticeably, only rescaling the displacements as x D – note the vertical scale of the bottom panels of Fig. 16.) It shows that the oscillation period is dominated by two similar parts, of equal duration. During these two intervals of relatively slow evolution, the limit cycle closely follows the declining branches of the function
x D
2
1 y y ,
(5.111)
0 0
corresponding to the zero value of the first (and nominally, the largest) term in the second equation of the system (110) – see the dashed line on the left bottom panel. During these intervals, the displacement x grows in accordance with the first of these equations, with its right-hand part virtually equal to the y 0
corresponding to Eq. (111). Even without solving the resulting differential equation exactly,47 we see that at these brunches, with y 0 1, x() changes with a speed of the order of D, and hence the path from the initial and final points of each branch, of a length x ~ D, takes a time interval of the order of 1, just as the right panel shows.
As soon as the system reaches the branch’s endpoint x = (2/33)D 0.385D, where the derivative dy 0/ dx diverges, the balance of the terms on the right-hand part of the second Eq. (110) is not more possible, and its magnitude abruptly becomes of the order of D >> 1. As a result, the system jumps from this point to the opposite branch of the curve (111) very rapidly, during a time interval ~ y 0/D
~ 1/D << 1, insufficient for x to change much. (The initial transient processes, i.e. the approaches to the limit cycle from almost arbitrary initial conditions, are equally fast, also with x const.) Upon reaching the new branch, the system “relaxes” to a relatively slow evolution in the opposite direction (“relaxation oscillations”), and the process repeats again and again.
Such oscillations take place in a large number of practical mechanical systems and electronic devices, ranging from the bowed string musical instruments (including those of the violin family), to usual mechanical clocks, to car light blinkers. Many of them allow for simple analyses; to save time/space, let me leave a couple of problems of this type for the reader’s exercise.
5.10. Exercise problems
m
5.1. A body of mass m is connected to its support not only with an elastic
spring but also with a damper (say, an air brake) that provides a drag force obeying Eq. (5) – see the figure on the right.
47 Its integration leads to an elementary function for ( y), but transcendental equations for y() and x().
Chapter 5
Page 34 of 38
CM: Classical Mechanics
(i) How to select the constants and to minimize the body’s vibrations caused by vertical oscillations of its support with frequency ?
(ii)* What if the oscillations are random?
5.2. For a system with the response function given by Eq. (17):
(i) prove Eq. (26), and
(ii) use an approach different from the one used in Sec. 1, to derive Eq. (34).
Hint: You may like to use the Cauchy integral theorem and the Cauchy integral formula for analytical functions of a complex variable.48
f ( t)
5.3. A square-wave pulse of force (see the figure on the right) is exerted
f 0
on a linear oscillator of frequency 0 (with no damping), initially at rest.
Calculate the law of motion q( t), sketch it, and interpret the result.
0
2 /
t
0
5.4. At t = 0, a sinusoidal external force F( t) = F 0cos t, starts to be exerted on a linear oscillator with frequency 0 and damping , which was at rest at t 0.
(i) Derive the general expression for the time evolution of the oscillator’s displacement, and interpret the result.
(ii) Spell out the result for the exact resonance ( = 0) in a system with low damping ( << 0), and, in particular, explore the limit 0.
5.5. A pulse of external force F( t), with a finite duration T, is exerted on a linear oscillator, initially at rest in its equilibrium position. Neglecting dissipation, calculate the change of the oscillator’s energy, using two different approaches, and compare the results.
5.6. A bead may slide, without friction, in a vertical plane along a parabolic curve y = x 2/2, in a uniform gravity field g = – gn y. Calculate the frequency of its free oscillations as a function of their amplitude A, in the first nonvanishing approximation in A 0, using two different approaches.
5.7. For a system with the Lagrangian function
m
2
2
4
L
q q q ,
2
2
with small parameter , use the harmonic balance method to find the frequency of free oscillations as a function of their amplitude.
5.8. Use a different approach to derive Eq. (49) for the frequency of free oscillations of the system described by Eq. (43) with = 0, in the first nonvanishing approximation in the small parameter
A 2/ 2
0 << 1.
5.9. On the plane [ a 1, a 2] of two real, time-independent parameters a 1 and a 2, find the regions in which the fixed point of the following system of equations,
48 See, e.g., MA Eq. (15.1).
Chapter 5
Page 35 of 38
CM: Classical Mechanics
q a ( q q ),
1
1
2
1
q a q q ,
2
2 1
2
is unstable, and sketch the regions of each fixed point type – stable and unstable nodes, focuses, etc.
5.10. Solve Problem 4(ii) using the reduced equations (57), and compare the result with the exact solution.
5.11. Use the reduced equations to analyze forced oscillations in an oscillator with weak nonlinear damping, described by the following equation:
q 2
2
3
q q q f cos t,
0
0
with 0; , > 0; and A 2 << 1. In particular, find the stationary amplitude of the forced oscillations and analyze their stability. Discuss the effect(s) of the nonlinear term on the resonance.
5.12. Within the approach discussed in Sec. 4, calculate the average frequency of a self-oscillator outside of the range of its phase-locking by a small external sinusoidal force.
5.13.* Use the reduced equations to analyze the stability of the forced nonlinear oscillations described by the Duffing equation (43). Relate the result to the slope of resonance curves (Fig. 4).
5.14. Use the van der Pol method to find the condition of
2
( t)
parametric excitation of the oscillator described by the following
0
2
2
equation:
201
q 2
2
q ( t) q 0,
0
2
0
where 2
0 ( t) is the square-wave function shown in the figure on the
t
right, with
201
0.
5.15. Use the van der Pol method to analyze parametric excitation of an oscillator with weak nonlinear damping, described by the following equation:
q 2
3
2
q
q
t q
0 1
cos 2
,
0
with 0; , > 0; and , A 2 << 1. In particular, find the amplitude of stationary oscillations and analyze their stability.
5.16.* Adding nonlinear term q 3 to the left-hand side of Eq. (75),
(i) find the corresponding addition to the reduced equations,
(ii) calculate the stationary amplitude A of the parametric oscillations,
(iii) find the type and stability of each fixed point of the reduced equations,
(iv) sketch the Poincaré phase planes of the system in major parameter regions.
5.17. Use the van der Pol method to find the condition of parametric excitation of an oscillator with weak modulation of both the effective mass m( t) = m 0(1 + m cos 2 t) and the effective spring constant ( t) = 0[1 + cos( 2 t –)], with the same frequency 2 20, at arbitrary modulation depths Chapter 5
Page 36 of 38
Essential Graduate Physics
CM: Classical Mechanics
ratio m/ k and phase shift . Interpret the result in terms of modulation of the instantaneous frequency
( t) [( t)/ m( t)]1/2 and the mechanical impedance Z( t) [( t) m( t)]1/2 of the oscillator.
5.18.* Find the condition of parametric excitation of a nonlinear oscillator described by the following equation:
q 2
2
2
q
q q f cos 2 t,
0
0
with sufficiently small , , f 0, and – 0.
5.19. Find the condition of stability of the equilibrium point q = 0 of a parametric oscillator described by Eq. (75), in the limit when << 0 << and << 1. Use the result to analyze the stability of the Kapitza pendulum mentioned in Sec. 5.
5.20.* Use numerical simulation to explore phase-plane trajectories [ q, q ] of an autonomous pendulum described by Eq. (42) with f 0 = 0, for both low and high damping, and discuss their most significant features.
5.21. Analyze relaxation oscillations of the system shown
m
in the figure on the right: an elastic spring prevents a block of
k
s
mass m from being carried away by a horizontal conveyer belt
moving with a constant velocity u. Assume that the coefficient k
u const
of the kinematic friction between the block and the belt is lower
than the static friction coefficient s.
5.22. The figure on the right shows the electric circuit of the
simplest relaxation oscillator. In it, N is a bistable circuit element that
switches very rapidly from its very-high-resistance state to a very-low- E , R
C
N
resistance state as the voltage across it is increased beyond some value
V t, and switches back as the voltage is decreased below another value
V t ’ < V t.49 Calculate the waveform and the time period of voltage oscillations in the circuit.
Hint: The solution of this problem requires a very basic understanding of electric circuits, including the e.m.f. E and internal resistance R of a dc current source such as an electric battery.
49 This is a good model for many two-terminal gas-discharge devices (such as glow lamps), whose effective resistance may drop by up to 5 orders of magnitude when the discharge has been activated by voltage V > V t. In usual neon glow lamps, V t ’ is about 30% lower than V t.
Chapter 5
Page 37 of 38
CM: Classical Mechanics
This page is
intentionally left
blank
Chapter 5
Page 38 of 38
Essential Graduate Physics
CM: Classical Mechanics
Chapter 6. From Oscillations to Waves
In this chapter, the discussion of oscillations is extended to systems with two and more degrees of freedom. This extension naturally leads to another key notion of physics – waves, so far in simple, mostly 1D systems. (In the next chapter, this discussion will be extended to more complex elastic continua.) However, even the limited scope of the models analyzed in this chapter will still enable us to discuss such important general aspects of waves as their dispersion, phase and group velocities, impedance, reflection, and attenuation.
6.1. Two coupled oscillators
Let us discuss oscillations in systems with several degrees of freedom, starting from the simplest case of two linear (harmonic), dissipation-free, 1D oscillators. If the oscillators are independent of each other, the Lagrangian function of their system may be represented as a sum of two independent terms of the type (5.1):
m
,
1 2
2
,
1 2
2
L L L , L T U
q
q .
(6.1)
1
2
,
1 2
,
1 2
,
1 2
,
1 2
,
1 2
2
2
Correspondingly, Eqs. (2.19) for qj = q 1,2 yields two independent equations of motion of the oscillators, each one being similar to Eq. (5.2):
2
2
,
1 2
m q m q ,
0
where
.
(6.2)
,
1 2
,
1 2
,
1 2
,
1 2
,
1 2
,
1 2
m ,12
(In the context of what follows, 1,2 are sometimes called the partial frequencies.) This means that in this simplest case, an arbitrary motion of the system is just a sum of independent sinusoidal oscillations at two frequencies equal to the partial frequencies (2).
However, as soon as the oscillators are coupled (i.e. interact), the full Lagrangian L contains an additional mixed term L int depending on both generalized coordinates q 1 and q 2 and/or generalized velocities. As a simple example, consider the system shown in Fig. 1, where two small masses m 1,2 are constrained to move in only one direction (shown horizontal), and are kept between two stiff walls with three springs.
L
m
M
m
R
1
2
Fig. 6.1. A simple system of two
coupled linear oscillators.
q
q
1
2
In this case, the kinetic energy is still separable, T = T 1 + T 2, but the total potential energy, consisting of the elastic energies of three springs, is not:
L
2
M
2
R
2
U
q
( q q )
q ,
(6.3a)
1
1
2
2
2
2
2
where q 1.2 are the horizontal displacements of the particles from their equilibrium positions. It is convenient to rewrite this expression as
© K. Likharev
CM: Classical Mechanics
1
2
2
2
U
q
q q q ,
where , , .
(6.3b)
1
2
1 2
1
L
M
2
R
M
M
2
2
This formula shows that the Lagrangian function L = T – U of this system contains, besides the partial terms (1), a bilinear interaction term:
L L L L , L q
q .
(6.4)
1
2
int
int
1 2
The resulting Lagrange equations of motion are as follows:
2
m q m q q
,
Linearly
1 1
1
1 1
2
(6.5)
coupled
2
m q m q q
.
oscillators
2 2
2
2 2
1
Thus the interaction leads to an effective generalized force q 2 exerted on subsystem 1 by subsystem 2, and the reciprocal effective force q 1.
Please note two important aspects of this (otherwise rather simple) system of equations. First, in contrast to the actual physical interaction forces (such as F 12 = – F 21 = M( q 2 – q 1) for our system1) the effective forces on the right-hand sides of Eqs. (5) do not obey the 3rd Newton law. Second, the forces are proportional to the same coefficient ; this feature is a result of the general bilinear structure (4) of the interaction energy, rather than of any special symmetry.
From our prior discussions, we already know how to solve Eqs. (5), because it is still a system of linear and homogeneous differential equations, so that its general solution is a sum of particular solutions of the form similar to Eqs. (5.88),
t
t
q c e ,
q c e ,
(6.6)
1
1
2
2
with all possible values of . These values may be found by plugging Eq. (6) into Eqs. (5), and requiring the resulting system of two linear, homogeneous algebraic equations for the distribution coefficients c 1,2, 2
2
m c m c c
,
1
1
1
1 1
2
(6.7)
2
2
m c m c c
,
2
2
2
2 2
1
to be self-consistent. In our particular case, we get a characteristic equation,
m ( 2
2
)
1
1
,
0
(6.8)
m ( 2
2
)
2
2
that is quadratic in 2, and thus has a simple analytical solution:
1/
2
1
1
2
2 2
Ω Ω
Ω Ω
Ω Ω
1
2
21
2
2
2
2
2
2
1
2
2
4
m m
1
2
1/ 2
1
2
2
1
2
Ω Ω
(6.9)
1
2
Ω2 Ω
1
2 2
2
.
2
4
m m
1
2
1 Using these expressions, Eqs. (5) may be readily obtained from the Newton laws, but the Lagrangian approach used above will make their generalization, in the next section, more straightforward.
Chapter 6
Page 2 of 30
CM: Classical Mechanics
According to Eqs. (2) and (3b), for any positive spring constants, the product 12 = ( L +
M)( R + M)/( m 1 m 2)1/2 is always larger than /( m 1 m 2)1/2 = M/( m 1 m 2)1/2, so that the square root in Eq.
(9) is always smaller than ( 2
2
1 +2 )/2. As a result, both values of 2 are negative, i.e. the general
solution to Eq. (5) is a sum of four terms, each proportional to exp{ i t}, where both normal frequencies (or “natural frequencies”, or “eigenfrequencies”) i are real: 1/
2
1
1
Anticrossing:
2
2
.
(6.10)
2 2
1
2
21
2
2
2
2
example
2
4
m m
1
2
A plot of these eigenfrequencies as a function of one of the partial frequencies (say, 1), with the other partial frequency fixed, gives us the famous anticrossing (also called the “avoided crossing” or
“non-crossing”) diagram – see Fig. 2. One can see that at weak coupling, the normal frequencies are close to the partial frequencies 1,2 everywhere besides a narrow range near the anticrossing point 1 =
2. Most remarkably, at passing through this region, + smoothly “switches” from following 2 to following 1 and vice versa.
2
2
2
Fig. 6.2. The anticrossing diagram for two
values of the normalized coupling strength
/( m
2
1 m 2)1/22 : 0.3 (red lines) and 0.1 (blue
lines). In this plot, 1 is assumed to be changed
by varying 1 rather than m 1, but in the opposite
2
2
case, the diagram is qualitatively similar.
0
2
1
The reason for this counterintuitive behavior may be found by examining the distribution coefficients c 1,2 corresponding to each branch of the diagram, which may be obtained by plugging the corresponding value of = – i back into Eqs. (7). For example, at the anticrossing point 1 = 2 , Eq. (10) is reduced to
2
2
2
(6.11)
m m
1
2
1
1/ 2
1 2 .
1/ 2
Plugging this expression back into any of Eqs. (7), we see that for the two branches of the anticrossing diagram, the distribution coefficient ratio is the same by magnitude but opposite by sign: 1/ 2
c
m
1
2
, at
.
(6.12)
1
2
c
m
2
1
In particular, if the system is symmetric ( m 1 = m 2, L = R), then at the upper branch, corresponding to + > -, we get c 1 = – c 2. This means that in this so-called hard mode,2 masses oscillate 2 In physics, the term “mode” (or “normal mode”) is typically used to describe the distribution of a variable in space, at its oscillations with a single frequency. In our current case, when the notion of space is reduced to two oscillator numbers, each mode is fully specified by the corresponding ratio of two distribution coefficients c 1,2.
Chapter 6
Page 3 of 30
CM: Classical Mechanics
in anti-phase: q 1( t) – q 2( t). The resulting substantial extension/compression of the middle spring (see Fig. 1 again) yields additional returning force which increases the oscillation frequency. On the contrary, at the lower branch, corresponding to –, the particle oscillations are in phase: c 1 = c 2, i.e. q 1( t)
q 2( t), so that the middle spring is neither stretched nor compressed at all. As a result, in this soft mode, the oscillation frequency - is lower than +, and does not depend on M:
2
2
L
R
.
(6.13)
m
m
m
Note that for both modes, the oscillations equally engage both particles.
Far from the anticrossing point, the situation is completely different. Indeed, a similar calculation of c 1,2 shows that on each branch of the diagram, the magnitude of one of the distribution coefficients is much larger than that of its counterpart. Hence, in this limit, any particular mode of oscillations involves virtually only one particle. A slow change of system parameters, bringing it through the anticrossing, leads, first, to a maximal delocalization of each mode at 1 = 2, and then to a restoration of the localization, but in a different partial degree of freedom.
We could readily carry out similar calculations for the case when the systems are coupled via their velocities, L q
m q , where m is a coupling coefficient – not necessarily a certain physical int
1 2
mass.3 The results are generally similar to those discussed above, again with the maximum level splitting at 1 = 2 :
m
2
2
2 1
,
(6.14)
1 m / m m
m m
1
2
1/ 2
1 2 1/2
the last relation being valid for weak coupling. The generalization to the case of simultaneous coordinate and velocity coupling is also straightforward – see the next section.4
One more property of weakly coupled oscillators is a periodic slow transfer of energy between them, especially strong at or near the anticrossing point 1 = 2. Let me leave an analysis of such transfer for the reader’s exercise. (Due to the importance of this effect for quantum mechanics, it will be discussed in detail in the QM part of this series.)
6.2. N coupled oscillators
The calculations of the previous section may be readily generalized to the case of an arbitrary number (say, N) of coupled harmonic oscillators, with an arbitrary type of coupling. It is obvious that in this case Eq. (4) should be replaced with
3 In mechanics, with q 1,2 standing for the actual displacements of particles, such coupling is not very natural, but there are many dynamic systems of non-mechanical nature in which such coupling is the most natural one. The simplest example is the system of two LC (“tank”) circuits, with either capacitive or inductive coupling. Indeed, as was discussed in Sec. 2.2, for such a system, the very notions of the potential and kinetic energies are conditional and interchangeable.
4 Note that the anticrossing diagram shown in Fig. 2, is even more ubiquitous in quantum mechanics, because, due to the time-oscillatory character of the Schrödinger equation solutions, a weak coupling of any two quantum states leads to qualitatively similar behavior of the eigenfrequencies of the system, and hence of its eigenenergies (“energy levels”) E = of the system.
Chapter 6
Page 4 of 30
CM: Classical Mechanics
N
N
L L
L .
(6.15)
j
jj'
j1
j, j'1
Moreover, we can generalize the above expression for the mixed terms Ljj’, taking into account their possible dependence not only on the generalized coordinates but also on the generalized velocities, in a bilinear form similar to Eq. (4). The resulting Lagrangian may be represented in a compact form, N
m
jj'
jj'
L
q q
q q ,
(6.16)
j
j'
j
j'
j, j' 1
2
2
where the off-diagonal terms are index-symmetric: mjj’ = mj’j, jj’ = j’j, and the factors ½ compensate for the double-counting of each term with j j’, at the summation over two independently running indices.
One may argue that Eq. (16) is quite general if we still want to keep the equations of motion linear – as they always are if the oscillations are small enough.
Plugging Eq. (16) into the general form (2.19) of the Lagrange equation, we get N equations of motion of the system, one for each value of the index j’ = 1, 2,…, N:
N
m q q
(6.17)
jj'
j
jj'
j 0.
j 1
Just as in the previous section, let us look for a particular solution to this system in the form
t
q c e .
(6.18)
j
j
As a result, we are getting a system of N linear, homogeneous algebraic equations, N
2
m
c
,
(6.19)
jj'
jj'
0
j
j 1
for the set of N distribution coefficients cj. The condition that this system is self-consistent is that the determinant of its matrix equals zero:
Det
2
m
(6.20)
jj'
jj'
0.
This characteristic equation is an algebraic equation of degree N for 2, and so has N roots (2) n. For any Hamiltonian system with stable equilibrium, the matrices mjj’ and jj’ ensure that all these roots are real and negative. As a result, the general solution to Eq. (17) is the sum of 2 N terms proportional to exp
{ i nt}, n = 1, 2,…, N, where all N normal frequencies ωn are real.
Plugging each of these 2 N values of = i n back into a particular set of linear equations (17), one can find the corresponding sets of distribution coefficients cj. Generally, the coefficients are complex, but to keep qj( t) real, the coefficients cj+ corresponding to = + i n, and cj- corresponding to
= – i n have to be complex-conjugate of each other. Since the sets of the distribution coefficients may be different for each n, they should be marked with two indices, j and n. Thus, at general initial conditions, the time evolution of the j th generalized coordinate may be represented as
N
N
1
q
c exp{ i t} c* exp{ i t}
Re
c
i
exp{ t}.
(6.21)
j
jn
n
jn
n
jn
n
2 n1
n1
Chapter 6
Page 5 of 30
CM: Classical Mechanics
This formula shows very clearly again the physical sense of the distribution coefficients cjn: a set of these coefficients, with different values of index j but the same mode number n, gives the complex amplitudes of oscillations of the coordinates in this mode, i.e. for the special initial conditions that ensure purely sinusoidal motion of all the system, with frequency n.
The calculation of the normal frequencies and the corresponding modes (distribution coefficient sets) of a particular coupled system with many degrees of freedom from Eqs. (19)–(20) is a task that frequently may be only done numerically.5 Let us discuss just two particular but very important cases.
First, let all the coupling coefficients be small in the following sense: mjj’ << mj mjj and jj’ << j
jj, for all j j, and all partial frequencies j ( j/ mj)1/2 be not too close to each other: 2
2
Ω Ω
m
j
j'
jj '
jj'
,
, for all j j’.
(6.22)
2
Ω
m
j
j
j
(Such situation frequently happens if parameters of the system are “random” in the sense that they do not follow any special, simple rule – for example, the one resulting from some simple symmetry of the system.) Results of the previous section imply that in this case, the coupling does not produce a noticeable change in the oscillation frequencies: { ωn} { j}. In this situation, oscillations at each eigenfrequency are heavily concentrated in one degree of freedom, i.e. in each set of the distribution coefficients cjn (for a given n), one coefficient’s magnitude is much larger than all others.
Now let the conditions (22) be valid for all but one pair of partial frequencies, say 1 and 2, while these two frequencies are so close that the coupling of the corresponding partial oscillators becomes essential. In this case, the approximation { ωn} { j} is still valid for all other degrees of freedom, and the corresponding terms may be neglected in Eqs. (19) for j = 1 and 2. As a result, we return to Eqs. (7) (perhaps generalized for velocity coupling) and hence to the anticrossing diagram (Fig.
2) discussed in the previous section. As a result, an extended change of only one partial frequency (say,
1) of a weakly coupled system produces a sequence of eigenfrequency anticrossings – see Fig. 3.
n
ω 4
4
ω 3
3
ω 2
2
ω 1
Fig. 6.3. The level anticrossing in a system of N
0
weakly coupled oscillators – schematically.
1
6.3. 1D waves
The second case when the general results of the last section may be simplified are coupled systems with a considerable degree of symmetry. Perhaps the most important of them are uniform 5 Fortunately, very effective algorithms have been developed for this matrix diagonalization task – see, e.g., references in MA Sec. 16(iii)-(iv). For example, the popular MATLAB software package was initially created exactly for this purpose. (“MAT” in its name stood for “matrix” rather than “mathematics”.) Chapter 6
Page 6 of 30
Essential Graduate Physics
CM: Classical Mechanics
systems that may sustain traveling and standing waves. Figure 4a shows a simple example of such a system – a long uniform chain of particles, of mass m, connected with light, elastic springs, pre-stretched with the tension force T to have equal lengths d. (This system may be understood as a natural generalization of the two-particle system considered in Sec. 1 – cf. Fig. 1.)
(a)
m
m
m
T
.
.
.
. T
d
d
m
m
m
(b)
.
.
.
.
q
q
j
q
1
j
j 1
( j )
1 d
jd
( j )
1 d
z
m
.
.
(c)
Fig. 6.4. (a) A uniform 1D chain
m
T
q
of elastically coupled particles,
T
j 1
.
. m
and their small (b) longitudinal
q
q
j
and (c) transverse displacements
j 1
(much exaggerated for clarity).
( j )
1 d
jd
( j )
1 d
z
The spring’s pre-stretch does not affect small longitudinal oscillations qj of the particles about their equilibrium positions zj = jd (where the integer j numbers the particles sequentially) – see Fig. 4b.6
Indeed, in the 2nd Newton law for such a longitudinal motion of the j th particle, the forces T and –T
exerted by the springs on the right and the left of it, cancel. However, the elastic additions, q, to these forces are generally different:
q
m ( q q ) ( q q ) .
(6.23)
j
j 1
j
j
j 1
On the contrary, for transverse oscillations within one plane (Fig. 4c), the net transverse component of the pre-stretch force exerted on the j th particle, Tt = T(sin+ – sin-), where are the force direction angles, does not vanish. As a result, direct contributions to this force from small longitudinal oscillations, with qj << d, T/, are negligible. Also, due to the first of these strong conditions, the angles are small, and hence may be approximated, respectively, as + ( qj+1 – qj)/ d and - ( qj – qj-1)/ d. Plugging these expressions into a similar approximation, Tt T(+ – -) for the transverse force, we see that it may be expressed as T( qj+1 – qj)/ d – T( qj – qj-1)/ d, i.e. is absolutely similar 6 Note the need for a clear distinction between the equilibrium position zj of the j th point and its deviation qj from it. Such distinction has to be sustained in the continuous limit (see below), where it is frequently called the Eulerian description – named after L. Euler, even though it was introduced to mechanics by J. d’Alembert. In this course, the distinction is emphasized by using different letters – respectively, z and q (in the 3D case, r and q).
Chapter 6
Page 7 of 30
CM: Classical Mechanics
to that in the longitudinal case, just with the replacement T/ d. As a result, we may write the equation of motion of the j th particle for these two cases in the same form:
q
m ( q
q ) ( q q ) ,
(6.24)
j
ef
j 1
j
ef
j
j 1
where ef is the “effective spring constant”, equal to for the longitudinal oscillations, and to T/ d for the transverse oscillations.7
Apart from the (formally) infinite size of the system, Eq. (24) is just a particular case of Eq. (17), and thus its particular solution may be looked for in the form (18), where, in light of our previous experience, we may immediately take 2 –2. With this substitution, Eq. (24) gives the following simple form of the general system of equations (19) for the distribution coefficients cj:
2
m 2 c c
c
.
(6.25)
ef
0
j
ef
j 1
ef
j 1
Now comes the most important conceptual step toward the wave theory. The translational symmetry of Eq. (25), i.e. its invariance with respect to the replacement j j + 1, allows it to have particular solutions of the following form:
i j
c ae
,
(6.26)
j
where the coefficient may depend on (and system’s parameters), but not on the particle number j.
Indeed, plugging Eq. (26) into Eq. (25) and canceling the common factor ei j, we see that this differential equation is indeed identically satisfied, provided that obeys the following algebraic characteristic equation:
2
m
2
i
e
i
e
.
(6.27)
ef
0
ef
ef
The physical sense of the solution (26) becomes clear if we use it and Eq. (18) with = i, to write 1D
q ( t) Re
,
(6.28) traveling
j
a exp i( kz t
j
) Re a exp ik( z v t
j
ph
)
wave
where the wave number k is defined as k / d. Eq. (28) describes a sinusoidal 8 traveling wave of particle displacements, which propagates, depending on the sign before v ph, to the right or the left along the particle chain, with the so-called phase velocity
v
.
(6.29) Phase
ph
k
velocity
Perhaps the most important characteristic of a wave system is the so-called dispersion relation, i.e. the relation between the wave’s frequency and its wave number k – one may say, between the temporal and spatial frequencies of the wave. For our current system, this relation is given by Eq. (27) with kd. Taking into account that (2 – e+ i – e-i) 2(1 – cos) 4sin2(/2), the dispersion relation may be rewritten in a simpler form:
7 The re-derivation of Eq. (24) from the Lagrangian formalism, with the simultaneous strict proof that the small oscillations in the longitudinal direction and the two mutually perpendicular transverse directions are all independent of each other, is a very good exercise, left for the reader.
8 In optics and quantum mechanics, such waves are usually called monochromatic; I will try to avoid this term until the corresponding parts (EM and QM) of my series.
Chapter 6
Page 8 of 30
CM: Classical Mechanics
1/ 2
kd
sin
sin
,
where
2
ef
.
(6.30)
max
2
max
2
max
m
This result, sketched in Fig. 5, is rather remarkable in several aspects. I will discuss them in some detail, because most of these features are typical for waves of any type (including even the “de Broglie waves”, i.e. wavefunctions, in quantum mechanics), propagating in periodic structures.
max
Fig. 6.5. The dispersion
relation (30).
0
2
kd
2
First, at low frequencies, << max, the dispersion relation (31) is linear:
d
d
1/ 2
vk,
v
max
ef
where
d .
(6.31)
dk
2
m
k 0
Plugging Eq. (31) into Eq. (29), we see that the constant v plays, in the low-frequency limit, the role of the same phase velocity for waves of any frequency. Due to its importance, this acoustic-wave 9 limit will with the subject of the special next section.
Second, when the wave frequency is comparable with max, the dispersion relation is not linear, and the system is dispersive. This means that as a wave, whose Fourier spectrum has several essential components with frequencies of the order of max, travels along the structure, its waveform (which may be defined as the shape of the line connecting all points qj( z), at the same time) changes.10 This effect may be analyzed by representing the general solution of Eq. (24) as the sum (more generally, an integral) of the components (28) with different complex amplitudes a:
1D wave
q t
( ) Re a exp
.
(6.32)
packet
j
k i kzj k t dk
This notation emphasizes the possible dependence of the component wave amplitudes ak and frequencies on the wave number k. While the latter dependence is given by the dispersion relation, in our current case by Eq. (30), the function ak is determined by the initial conditions. For applications, the case when ak is substantially different from zero only in a narrow interval, of a width k << k 0 around some central value k 0, is of special importance. The Fourier transform reciprocal to Eq. (32) shows that this is true, in particular, for the so-called wave packet – a sinusoidal (“carrier”) wave modulated by a spatial envelope function of a large width z ~ 1/ k >> 1/ k 0 – see, e.g., Fig. 6.
9 This term is purely historical. Though the usual sound waves in air, which are the subject of acoustics, belong to this class, the waves we are discussing may have frequencies both well below and well above the human ear’s sensitivity range.
10 The waveform’s deformation due to dispersion (which we are considering now) should be clearly distinguished from its possible change due to attenuation due to energy dissipation – which is not taken into account is our current energy-conserving model – cf. Sec. 6 below.
Chapter 6
Page 9 of 30
CM: Classical Mechanics
q j
v gr
0
z
Fig. 6.6. The phase
j
and group velocities
of a wave packet.
v ph
Using the strong inequality k << k 0, the wave packet’s propagation may be analyzed by expending the dispersion relation ( k) into the Taylor series at point k 0, and, in the first approximation in k/ k 0, restricting the expansion to its first two terms:
d
~
~
( k)
k ,
where
k
k k k
(6.33)
0
0
0
and
,
.
0
dk k k 0
In this approximation, Eq. (32) yields
d
q ( t) Re
j
ak i
~
exp
k k z
j
k t dk
0
~
0
k k
dk
0
(6.34)
d
