14.6: General Analytic Theory for Coupled Linear Oscillators
( \newcommand{\kernel}{\mathrm{null}\,}\)
The discussion of a coupled double-oscillator system in Section 14.5 has shown that it is possible to select symmetric and antisymmetric normal modes that are independent and each have characteristic frequencies. The normal coordinates for these two normal modes correspond to linear superpositions of the spatial amplitudes of the two oscillators and can be obtained by a rotation into the appropriate normal coordinate system. Extension of this to systems comprising n coupled linear oscillators, requires development of a general analytic theory, that is capable of finding the normal modes and their eigenvalues and eigenvectors. As illustrated for the double oscillator, the solution of many coupled linear oscillators is a classic eigenvalue problem where one has to rotate to the principal axis system to project out the normal modes. The following discussion presents a general approach to the problem of finding the normal coordinates for a system of n coupled linear oscillators.
Consider a conservative system of n coupled oscillators, described in terms of generalized coordinates q_k and t with subscript k = 1, 2, 3, \ldots, n for a system with n degrees of freedom. The coupled oscillators are assumed to have a stable equilibrium with generalized coordinates q_{k0} at equilibrium. In addition, it is assumed that the oscillation amplitudes are sufficiently small to ensure that the system is linear.
For the equilibrium position q_k = q_{k0}, the Lagrange equations must satisfy
\begin{align} \dot{q}_k = 0 \label{12.27} \\ \ddot{q}_k = 0 \notag \end{align}
Every non-zero term of the form \frac{d}{dt}\frac{\partial L}{\partial \dot{q}_k} in Lagrange’s equations must contain at least either \dot{q}_k or \ddot{q}_k which are zero at equilibrium; thus all such terms vanish at equilibrium. At equilibrium
\left(\frac{\partial L}{\partial q_k}\right)_0 = \left(\frac{\partial T}{\partial q_k}\right)_0 − \left(\frac{\partial U}{\partial q_k}\right)_0 = 0 \label{12.28}
where the subscript 0 designates at equilibrium.
Kinetic energy tensor T
In chapter 7.6 it was shown that, in terms of fixed rectangular coordinates, the kinetic energy for N bodies, with n generalized coordinates, is expressed as
T = \frac{1}{2} \sum^N_{ \alpha =1} \sum^3_{i=1 } m_{\alpha} \dot{x}^2_{\alpha ,i} \label{12.29}
Expressing these in terms of generalized coordinates x_{\alpha ,i} = x_{\alpha ,i}(q_j , t) where j = 1, 2, ...n, then the generalized velocities are given by
\dot{x}_{\alpha ,i} = \sum^n_{j=1} \frac{\partial x_{\alpha ,i}}{\partial q_j} \dot{q}_j + \frac{\partial x_{\alpha ,i}}{\partial t} \label{12.30}
As discussed in chapter 7.6, if the system is scleronomic then the partial derivative
\frac{\partial x_{\alpha ,i}}{\partial t} = 0 \label{12.31}
Thus the kinetic energy, Equation \ref{12.29}, of a scleronomic system can be written as a homogeneous quadratic function of the generalized velocities
T = \frac{1}{2} \sum^{n}_{j,k} T_{jk} \dot{q}_j \dot{q}_k \label{12.32}
where the components of the kinetic energy tensor \mathbf{T} are
T_{jk} \equiv \sum^{N}_{\alpha} m_{\alpha} \sum^3_i \frac{ \partial x_{\alpha ,i}}{\partial q_j} \frac{ \partial x_{\alpha ,i}}{\partial q_k} \label{12.33}
Note that if the velocities \dot{q} correspond to translational velocity, then the kinetic energy tensor \mathbf{T} corresponds to an effective mass tensor, whereas if the velocities correspond to angular rotational velocities, then the kinetic energy tensor \mathbf{T} corresponds to the inertia tensor.
It is possible to make an expansion of the T_{jk} about the equilibrium values of the form
T_{jk} (q_1, q_2, ..q_n) = T_{jk} (q_{i0}) + \sum_l \left(\frac{\partial T_{jk}}{\partial q_l}\right)_0 q_l + ... \label{12.34}
Only the first-order term will be kept since the second and higher terms are of the same order as the higher order terms ignored in the Taylor expansion of the potential. Thus, at the equilibrium point, assume that \left(\frac{\partial T}{\partial q_k} \right)_0 = 0 where k = 1, 2, 3, ...n.
Potential energy tensor V
Equations \ref{12.28} plus \ref{12.34} imply that
\left(\frac{\partial U}{\partial q_k}\right)_0 = 0 \label{12.35}
where k = 1, 2, 3, ...n.
Make a Taylor expansion about equilibrium for the potential energy, assuming for simplicity that the coordinates have been translated to ensure that q_k = 0 at equilibrium. This gives
U (q_1, q_2, ..q_n) = U_0 + \sum_k \left(\frac{\partial U}{\partial q_k}\right)_0 q_k + \frac{1}{2} \sum_{ j,k} \left( \frac{\partial^2 U}{ \partial q_j \partial q_k } \right)_0 q_j q_k + .. \label{12.36}
The linear term is zero since \left( \frac{\partial U}{\partial q_k} \right)_0 = 0 at the equilibrium point, and without loss of generality, the potential can be measured with respect to U_0. Assume that the amplitudes are small, then the expansion can be restricted to the quadratic term, corresponding to the simple linear oscillator potential
U (q_1, q_2, ..q_n) − U_0 = U^{\prime} (q_1, q_2, ..q_n) = \frac{1}{2} \sum_{ j,k} \left( \frac{\partial^2 U} {\partial q_j\partial q_k} \right)_0 q_j q_k = \frac{1}{2} \sum_{j,k} V_{jk} q_j q_k \label{12.37}
That is
U^{\prime} (q_1, q_2, ..q_n) = \frac{1}{2} \sum_{j,k} V_{jk} q_j q_k \label{12.38}
where the components of the potential energy tensor \mathbf{V} are defined as
V_{jk} \equiv \left( \frac{\partial^2 U^{\prime}}{ \partial q_j\partial q_k} \right)_0 \label{12.39}
Note that the order of differentiation is unimportant and thus the quantity V_{jk} is symmetric
V_{jk} = V_{kj} \label{12.40}
The motion of the system has been specified for small oscillations around the equilibrium position and it has been shown that U^{\prime} (q_1, q_2, ...q_n) has a minimum value at equilibrium which is taken to be zero for convenience.
In conclusion, equations \ref{12.32} and \ref{12.38} give
T = \frac{1}{2} \sum^n_{j,k} T_{jk} \dot{q}_j \dot{q}_k \label{12.41}
U^{\prime} = \frac{1}{2} \sum^n_{j,k} V_{jk} q_j q_k \label{12.42}
where the components of the kinetic energy tensor \mathbf{T} and potential energy tensor \mathbf{V} are
T_{j k} \equiv\left(\sum_{\alpha}^{N} m_{\alpha} \sum_{i}^{3} \frac{\partial x_{\alpha, i}}{\partial q_{j}} \frac{\partial x_{\alpha, i}}{\partial q_{k}}\right)_{0} \label{12.43}
V_{j k} \equiv\left(\frac{\partial^{2} U^{\prime}}{\partial q_{j} \partial q_{k}}\right)_{0} \label{12.44}
Note that q_j and q_k may have different units, but all the terms in the summations for both T and U^{\prime}, have units of energy. The V_{jk} and T_{jk} values are evaluated at the equilibrium point, and thus both V_{jk} and T_{jk} are n \times n arrays of values evaluated at the equilibrium location.
Equations of motion
Both the kinetic energy and potential energy terms are products of the coordinates leading to a set of coupled equations that are complicated to solve. The problem is greatly simplified by selecting a set of normal coordinates for which both T and U are diagonal, then the coupling terms disappear. Thus a coordinate transformation must be found that simultaneously diagonalizes T_{jk} and V_{jk} in order to obtain a set of normal coordinates.
The kinetic energy T is only a function of generalized velocities \dot{q}_k while the conservative potential energy is only a function of the generalized coordinates q_k. Thus the Lagrange equations
\frac{\partial L}{\partial q_k} − \frac{d}{dt} \frac{\partial L}{\partial \dot{q}_k} = 0 \label{12.45}
reduce to
\frac{\partial U}{\partial q_k} + \frac{d}{dt} \frac{\partial T}{\partial \dot{q}_k} = 0 \label{12.46}
But
\frac{\partial U}{\partial q_k} = \sum^n_j V_{jk} q_j \label{12.47}
and
\frac{\partial T}{\partial \dot{q}_k} = \sum^n_j T_{jk} \dot{q}_j \label{12.48}
Thus the Lagrange equations reduce to the following set of equations of motion,
\sum^n_j (V_{jk} q_j + T_{jk} \ddot{q}_j )=0 \label{12.49}
For each k, where 1 \leq k \leq n, there exists a set of n second-order linear homogeneous differential equations with constant coefficients. Since the system is oscillatory, it is natural to try a solution of the form
q_j (t) = a_j e^{i(\omega t−\delta )} \label{12.50}
Assuming that the system is conservative, then this implies that \omega is real, since an imaginary term for \omega would lead to an exponential damping term. The arbitrary constants are the real amplitude a_j and the phase \delta. Substitution of this trial solution for each k leads to a set of equations
\sum_j (V_{jk} − \omega^2 T_{jk} ) a_j = 0 \label{12.51}
where the common factor e^{i(\omega t−\delta )} has been removed. Equation \ref{12.51} corresponds to a set of n linear homogeneous algebraic equations that the a_j amplitudes must satisfy for each k. For a non-trivial solution to exist, the determinant of the coefficients must vanish, that is
\begin{vmatrix} V_{11} − \omega^2 T_{11} & V_{12} − \omega^2 T_{12} & V_{13} − \omega^2 T_{13} & ... \\ V_{12} − \omega^2 T_{12} & V_{22} − \omega^2 T_{22} & V_{23} − \omega^2 T_{23} & ... \\ V_{13} − \omega^2 T_{13} & V_{23} − \omega^2 T_{23} & V_{33} − \omega^2 T_{33} & ... \\ ... & ... & ... & ... \end{vmatrix} = 0 \label{12.52}
where the symmetry V_{jk} = V_{kj} has been included. This is the standard eigenvalue problem for which the above determinant gives the secular equation or the characteristic equation. It is an equation of degree n in \omega^2. The n roots of this equation are \omega^2_r where \omega_r are the characteristic frequencies or eigenfrequencies of the normal modes.
Substitution of \omega^2_r into Equation \ref{12.52} determines the ratio a_{1,r} : a_{2,r} : a_{3,r} : ... : a_{n,r} for this solution which defines the components of the n-dimensional eigenvector \mathbf{a}_r. That is, solution of the secular equations have determined the eigenvalues and eigenvectors of the n solutions of the coupled-channel system.
Superposition
The equations of motion \sum_j (V_{jk} q_j + T_{jk} \ddot{q}_j )=0 are linear equations that satisfy superposition. Thus the most general solution q_j (t) can be a superposition of the n eigenvectors \mathbf{a}_{jr}, that is
q_j (t) = \sum^n_r a_{jr} e^{i(\omega_rt−\delta_r)} \label{12.53}
Only the real part of q_j (t) is meaningful, that is,
q_j (t) = \text{ Re} \sum^n_r a_{jr} e^{i(\omega_r t−\delta_r)} = \sum^n_r a_{jr} \cos (\omega_r t − \delta_r) \label{12.54}
Thus the most general solution of these linear equations involves a sum over the eigenvectors of the system which are cosine functions of the corresponding eigenfrequencies.
Eigenfunction Orthonormality
It can be shown that the eigenvectors are orthogonal. In addition, the above procedure only determines ratios of amplitudes, thus there is an indeterminacy that can be used to normalize the a_{jr}. Thus the eigenvectors form an orthonormal set. Orthonormality of the eigenfunctions for the rank 3 inertia tensor was illustrated in chapter 13.10.2. Similar arguments apply that allow extending orthonormality to higher rank cases such that for n-body coupled oscillators.
The eigenfunction orthogonality for n coupled oscillators can be proved by writing Equation \ref{12.51} for both the s^{th} root and the r^{th} root. That is,
\sum_j V_{jk} a_{ks} = \omega^2_s \sum_j T_{jk} a_{ks} \label{12.55}
\sum_j V_{jk} a_{jr} = \omega^2_r \sum_j T_{jk} a_{jr} \label{12.56}
Multiply Equation \ref{12.55} by a_{jr} and sum over k. Similarly multiply Equation \ref{12.56} by a_{ks} and sum over k. These summations lead to
\sum_{jk} V_{jk} a_{jr}a_{ks} = \omega^2_s \sum_{jk} T_{jk} a_{jr} a_{ks} \label{12.57}
\sum_{jk} V_{jk} a_{jr}a_{ks} = \omega^2_r \sum_{jk} T_{jk} a_{jr}a_{ks} \label{12.58}
Note that the left-hand sides of these two equations are identical. Thus taking the difference between these equations gives
(\omega^2_r − \omega^2_s) \sum_{jk} T_{jk} a_{jr} a_{ks} = 0 \label{12.59}
Note that if (\omega^2_r − \omega^2_s ) \neq 0, that is, assuming that the eigenfrequencies are not degenerate, then to ensure that Equation \ref{12.59} is zero requires that
\sum_{jk} T_{jk} a_{jr} a_{ks} = 0 \quad r \neq s \label{12.60}
This shows that the eigenfunctions are orthogonal. If the eigenfrequencies are degenerate, i.e. \omega^2_r = \omega^2_s, then, with no loss of generality, the axes r and s can be chosen to be orthogonal.
The eigenfunction normalization can be chosen freely since only ratios of the eigenfunction components a_{jr} are determined when \omega_r is used in Equation \ref{12.51}. The kinetic energy, given by Equation \ref{12.32} must be positive, or zero for the case of a static system. That is
T = \frac{1}{2} \sum^n_{j,k} T_{jk} \dot{q}_j \dot{q}_k \geq 0 \label{12.61}
Use the time derivative of Equation \ref{12.54} to determine \dot{q}_r and insert into Equation \ref{12.61} gives that the kinetic energy is
T=\frac{1}{2} \sum_{j, k}^{n} T_{j k} \dot{q}_{j} \dot{q}_{k}=\frac{1}{2} \sum_{j, k}^{n} T_{j k} \sum_{r, s} \omega_{r} \omega_{s} a_{j r} \cos \left(\omega_{r} t-\delta_{r}\right) a_{k s} \cos \left(\omega_{s} t-\delta_{s}\right) \label{12.62}
For the diagonal term r = s
T=\frac{1}{2} \sum_{j, k}^{n} T_{j k} \dot{q}_{j} \dot{q}_{k}=\left[\frac{1}{2} \sum_{r}^{n} \omega_{r}^{2} \cos ^{2}\left(\omega_{r} t-\delta_{r}\right)\right] \sum_{j, k} T_{j k} a_{j r} a_{k r} \geq 0 \label{12.63}
Since the term in the square brackets must be positive, then
\sum_{j,k} T_{jk} a_{jr} a_{kr} \geq 0 \label{12.64}
Since this sum must be a positive number, and the magnitude of the amplitudes can be chosen freely, then it is possible to normalize the eigenfunction amplitudes to unity. That is, choose that
\sum_{j,k} T_{jk} a_{jr} a_{ks} = 1 \label{12.65}
The orthogonality equation, \ref{12.60} and the normalization Equation \ref{12.65} can be combined into a single orthonormalization equation
\sum_{j,k} T_{jk} a_{jr} a_{ks} = \delta_{rs} \label{12.66}
This has shown that the eigenvectors form an orthonormal set.
Since the j^{th} component of the r^{th} eigenvector is a_{jr}, then the r^{th} eigenvector can be written in the form
\mathbf{a}_r = \sum_j a_{jr} \widehat{\mathbf{e}_j} \label{12.67}
where \widehat{\mathbf{e}_j} are the unit vectors for the generalized coordinates.
Normal coordinates
The above general solution of the coupled-oscillator problem is best expressed in terms of the normal coordinates which are independent. It is more transparent if the superposition of the normal modes are written in the form
q_{j}(t)=\sum_{r}^{n} \beta_{r} a_{j r} e^{i \omega_{r} t} \label{12.68}
where the complex factor \beta_r includes the arbitrary scale factor to allow for arbitrary amplitudes q_j as well as the fact that the amplitudes a_{jr} have been normalized and the phase factor \delta_r has been chosen.
Define
\eta_{r} (t) \equiv \beta_{r} e^{i \omega_{r} t} \label{12.69}
then Equation \ref{12.68} can be written as
q_j (t) = \sum^n_r a_{jr}\eta_r (t) \label{12.70}
Equation \ref{12.70} can be expressed schematically as the matrix multiplication
{\bf q = \{a\} \cdot} \boldsymbol{\eta} \label{12.71}
The \eta_r (t) are the normal coordinates which can be expressed in the form
\boldsymbol{\eta} {\bf = \{a\}^{−1} q } \label{12.72}
Each normal mode \eta_r corresponds to a single eigenfrequency, \omega_r which satisfies the linear oscillator equation
\ddot{\eta}_r + \omega^2_r \eta_r = 0 \label{12.73}