DSPRelated.com
Free Books
   Free Books    Physical Audio Signal Processing  

Wave Digital Mass-Spring Oscillator

Let's look again at the mass-spring oscillator of §F.3.4, but this time without the driving force (which effectively decouples the mass and spring into separate first-order systems). The physical diagram and equivalent circuit are shown in Fig.F.32 and Fig.F.33, respectively.

Figure F.32: Elementary mass-spring oscillator.
\includegraphics{eps/tank}

Figure F.33: Equivalent circuit for the mass-spring oscillator.
\includegraphics{eps/tankec}

Note that the mass and spring can be regarded as being in parallel or in series. Under the parallel interpretation, we have the WDF shown in Fig.F.34 and Fig.F.35.F.5 The reflection coefficient $ \rho$ can be computed, as usual, from the first alpha parameter:

$\displaystyle \rho = \alpha_1 - 1 = \frac{2\Gamma _1}{\Gamma _1+\Gamma _2} - 1 ... ..._1-\Gamma _2}{\Gamma _1+\Gamma _2} = \frac{R_2-R_1}{R_2+R_1} = \frac{m-k}{m+k} $

This result, $ \rho=(m-k)/(m+k)$, is just the ``impedance step over impedance sum'', so no calculation was really necessary.

Figure F.34: Wave digital mass-spring oscillator.
\includegraphics{eps/tankjunc}

Figure F.35: Detailed wave-flow diagram for the wave digital mass-spring oscillator.
\includegraphics{eps/tankwdf}

Oscillation Frequency

From Fig.F.33, we can see that the impedance of the parallel combination of the mass and spring is given by

$\displaystyle R_{m\vert\vert k}(s) \isdef \left.\frac{k}{s} \right\Vert ms = \frac{\frac{k}{s}ms}{\frac{k}{s}+ms} = \frac{ks}{s^2+\frac{k}{m}} \protect$ (F.38)

(using the product-over-sum rule for combining impedances in parallel). The poles of this impedance are given by the roots of the denominator polynomial in $ s$:

$\displaystyle s = \pm j\sqrt{\frac{k}{m}} \protect$ (F.39)

The resonance frequency of the mass-spring oscillator is therefore

$\displaystyle \omega_0 = \sqrt{\frac{k}{m}} \protect$ (F.40)

Since the poles $ s=\pm j\omega_0$ are on the $ j\omega $ axis, there is no damping, as we expect.

We can now write reflection coefficient $ \rho$ (see Fig.F.35) as

$\displaystyle \rho = \frac{m-k}{m+k} = \frac{1-\frac{k}{m}}{1+\frac{k}{m}} = \frac{1-\omega_0^2}{1+\omega_0^2} $

We see that dc ( $ \omega_0=0$) corresponds to $ \rho=1$, and $ \omega_0=\infty$ corresponds to $ \rho=-1$.

DC Analysis of the WD Mass-Spring Oscillator

Considering the dc case first ($ \rho=1$), we see from Fig.F.35 that the state variable $ x_1(n)$ will circulate unchanged in the isolated loop on the left. Let's call this value $ x_1(n)\equiv\message{CHANGE eqv TO equiv IN SOURCE} x_0$. Then the physical force on the spring is always equal to

$\displaystyle f_k(n) = f^{{+}}_k(n) + f^{{-}}_k(n) = 2x_1(n) = 2 x_0. \qquad\hbox{(spring force, dc case)} \protect$ (F.41)

The loop on the right in Fig.F.35 receives $ 2 x_0$ and adds $ x_2(n)$ to that. Since $ x_2(n+1) = 2x_1(n)+x_2(n)$, we see it is linearly growing in amplitude. For example, if $ x_2(0)=0$ (with $ x_1(0)=x_0$), we obtain $ x_2=[0, 2x_0, 4x_0, 6x_0,\ldots]$, or

$\displaystyle x_2(n) = 2 n x_0, \quad n=0,1,2,3,\ldots\,. \protect$ (F.42)

At first, this result might appear to contradict conservation of energy, since the state amplitude seems to be growing without bound. However, the physical force is fortunately better behaved:

$\displaystyle f_m(n) = f^{{+}}_m(n) + f^{{-}}_m(n) = x_2(n+1) - x_2(n) = 2x_0. \protect$ (F.43)

Since the spring and mass are connected in parallel, it must be the true that they are subjected to the same physical force at all times. Comparing Equations (F.41-F.43) verifies this to be the case.

WD Mass-Spring Oscillator at Half the Sampling Rate

Under the bilinear transform, the $ s=\infty$ maps to $ z=-1$ (half the sampling rate). It is therefore no surprise that given $ \omega_0=\infty$ ($ \rho=-1$), inspection of Fig.F.35 reveals that any alternating sequence (sinusoid sampled at half the sampling rate) will circulate unchanged in the loop on the right, which is now isolated. Let $ x_2(n) = (-1)^n x_0$ denote this alternating sequence. The loop on the left receives $ - 2 x_2(n)$ and adds $ - x_1(n-1)$ to it, i.e., $ x_1(n+1)= - x_1(n) - 2 x_2(n) = -x_1(n) - 2x_2(0)(-1)^n$. If we start out with $ x_1(0)=0$ and $ x_2(0)=x_0$, we obtain $ x_1 = [0,-2x_0, 4x_0, -6x_0, \ldots]$, or

$\displaystyle x_1(n) = (-1)^n 2 n x_0, \quad n=0,1,2,3,\ldots\,. $

However, the physical spring force is well behaved, since

$\displaystyle f_k(n) = f^{{+}}_k(n) + f^{{-}}_k(n) = x_1(n+1) + x_1(n) = (-1)^{n+1}2 x_0 $

As a check, the mass force is found to be

\begin{eqnarray*} f_m(n) &=& x_2(n+1) - x_2(n)\\ &=& (-1)^{n+1}x_0 - (-1)^n x_0\\ &=& (-1)^{n+1}x_0 + (-1)^{n+1}x_0\\ &=& (-1)^{n+1}2 x_0, \end{eqnarray*}

which agrees with the spring, as it must.

Linearly Growing State Variables in WD Mass-Spring Oscillator

It may seem disturbing that such a simple, passive, physically rigorous simulation of a mass-spring oscillator should have to make use of state variables which grow without bound for the limiting cases of simple harmonic motion at frequencies zero and half the sampling rate. This is obviously a valid concern in practice as well. However, it is easy to show that this only happens at dc and $ f_s/2$, and that there is a true degeneracy at these frequencies, even in the physics. For all frequencies in the audio range (e.g., for typical sampling rates), such state variable growth cannot occur. Let's take closer look at this phenomenon, first from a signal processing point of view, and second from a physical point of view.

A Signal Processing Perspective on Repeated Mass-Spring Poles

Going back to the poles of the mass-spring system in Eq.$ \,$(F.39), we see that, as the imaginary part of the two poles, $ \pm\omega_0 = \pm\sqrt{k/m}$, approach zero, they come together at $ s=0$ to create a repeated pole. The same thing happens at $ \omega_0=\infty$ since both poles go to ``the point at infinity''.

It is a well known fact from linear systems theory that two poles at the same point $ s=s_0=\sigma_0$ in the $ s$ plane can correspond to an impulse-response component of the form $ te^{\sigma_0 t}$, in addition to the component $ e^{\sigma_0 t}$ produced by a single pole at $ s=\sigma_0$. In the discrete-time case, a double pole at $ z=r_0$ can give rise to an impulse-response component of the form $ n r_0^n$. This is the fundamental source of the linearly growing internal states of the wave digital sine oscillator at dc and $ f_s/2$. It is interesting to note, however, that such modes are always unobservable at any physical output such as the mass force or spring force that is not actually linearly growing.

Physical Perspective on Repeated Poles in the Mass-Spring System

In the physical system, dc and infinite frequency are in fact strange cases. In the case of dc, for example, a nonzero constant force implies that the mass $ m$ is under constant acceleration. It is therefore the case that its velocity is linearly growing. Our simulation predicts this, since, using Eq.$ \,$(F.43) and Eq.$ \,$(F.42),

\begin{eqnarray*} v_m(n) &=& \frac{f^{{+}}_m(n)}{m} - \frac{f^{{-}}_m(n)}{m} =... ...m} \left[2(n+1) + 2n\right]x_0 = \frac{1}{m} (4 n x_0 + 2 x_0). \end{eqnarray*}

The dc term $ 2x_0/m$ is therefore accompanied by a linearly growing term $ 2nx_0/m$ in the physical mass velocity. It is therefore unavoidable that we have some means of producing an unbounded, linearly growing output variable.

Mass-Spring Boundedness in Reality

To approach the limit of $ \omega_0 = \sqrt{k/m} = 0$, we must either take the spring constant $ k$ to zero, or the mass $ m$ to infinity, or both.

In the case of $ k\to0$, the constant force must approach zero, and we are left with at most a constant mass velocity in the limit (not a linearly growing one, since there can be no dc force at the limit). When the spring force reaches zero, $ x_1(n)=0$, so that only zeros will feed into the loop on the right in Fig.F.35, thus avoiding a linearly growing velocity, as demanded by the physics. (A constant velocity is free to circulate in the loop on the right, but the loop on the left must be zeroed out in the limit.)

In the case of $ m\to\infty$, the mass becomes unaffected by the spring force, so its final velocity must be zero. Otherwise, the attached spring would keep compressing or stretching forever, and this would take infinite energy. (Another way to arrive at this conclusion is to note that the final kinetic energy of the mass would be $ mv^2/2=\infty$.) Since the total energy in an undriven mass-spring oscillator is always constant, the infinite-mass limit must be accompanied by a zero-velocity limit.F.6 This means the mass's state variable $ x_2(n)$ in Fig.F.35 must be forced to zero in the limit so that there will be no linearly growing solution at dc.

In summary, when two or more system poles approach each other to form a repeated pole, care must be taken to ensure that the limit is approached in a physically meaningful way. In the case of the mass-spring oscillator, for example, any change in the spring constant $ k$ or mass $ m$ must be accompanied by the physically appropriate change in the state variables $ x_1(n)$ and/or $ x_2(n)$. It is obviously incorrect, for example, to suddenly set $ k=0$ in the simulation without simultaneously clearing the spring's state variable $ x_1(n)$, since the force across an infinitely compliant spring can only be zero.

Similar remarks apply to repeated poles corresponding to $ \omega_0=\infty$. In this case, the mass and spring basically change places.

Energy-Preserving Parameter Changes (Mass-Spring Oscillator)

If the change in $ k$ or $ m$ is deemed to be ``internal'', that is, involving no external interactions, the appropriate accompanying change in the internal state variables is that which conserves energy. For the mass and its velocity, for example, we must have

$\displaystyle \frac{1}{2} m_1 v_1^2 =\frac{1}{2} m_2 v_2^2 $

where $ m_1,m_2$ denote the mass values before and after the change, respectively, and $ v_1,v_2$ denote the corresponding velocities. The velocity must therefore be scaled according to

$\displaystyle v_2 = v_1\sqrt{\frac{m_1}{m_2}}, $

since this holds the kinetic energy of the mass constant. Note that the momentum of the mass is changed, however, since

$\displaystyle m_2v_2 = m_2 v_1\sqrt{\frac{m_1}{m_2}} =v_1\sqrt{m_1m_2} =m_1v_1\sqrt{\frac{m_2}{m_1}} $

If the spring constant $ k$ is to change from $ k_1$ to $ k_2$, the instantaneous spring displacement $ x$ must satisfy

$\displaystyle \frac{1}{2} k_1 x_1^2 =\frac{1}{2} k_2 x_2^2 $

In a velocity-wave simulation, displacement is the integral of velocity. Therefore, the energy-conserving velocity correction is impulsive in this case.

Exercises in Wave Digital Modeling

  1. Comparing digital and analog frequency formulas. This first exercise verifies that the elementary ``tank circuit'' always resonates at exactly the frequency it should, according to the bilinear transform frequency mapping $ \omega_a = \tan(\omega_d T / 2)$, where $ \omega_a$ denotes ``analog frequency'' and $ \omega_d$ denotes ``digital frequency''.
    1. Find the poles of Fig.F.35 in terms of $ \rho$.

    2. Show that the resonance frequency is given by

      $\displaystyle f_s\arccos\left(\rho\right) $

      where $ f_s$ denotes the sampling rate.

    3. Recall that the mass-spring oscillator resonates at $ \omega_0=\sqrt{k/m}$. Relate these two resonance frequency formulas via the analog-digital frequency map $ \omega_a = \tan(\omega_d T / 2)$.

    4. Show that the trig identity you discovered in this way is true. I.e., show that

      $\displaystyle f_s \arccos\left[\frac{k-m}{k+m}\right] = 2f_s \arctan\left[\sqrt{\frac{m}{k}}\right]. $


Previous Section:
Mass and Dashpot in Series