DSPRelated.com
Free Books

Introduction

The digital waveguide (DW) method has been used for many years to provide highly efficient algorithms for musical sound synthesis based on physical models [433,447,396]. For a much longer time, finite-difference time-domain (FDTD) schemes have been used to simulate more general situations, usually at higher cost [550,392,74,77,45,397]. In recent years, there has been interest in relating these methods to each other [123] and in combining them for more general simulations. For example, modular hybrid methods have been devised which interconnect DW and FDTD simulations by means of a KW converter [223,226]. The basic idea of the KW-converter adaptor is to convert the ``Kirchoff variables'' of the FDTD, such as string displacement, velocity, etc., to ``wave variables'' of the DW. The W variables are regarded as the traveling-wave components of the K variables.

In this appendix, we present an alternative to the KW converter. Instead of converting K variables to W variables, or vice versa, in the time domain, conversion formulas are derived with respect to the current state as a function of spatial coordinates. As a result, it becomes simple to convert any instantaneous state configuration from FDTD to DW form, or vice versa. Thus, instead of providing the necessary time-domain filter to implement a KW converter converting traveling-wave components to physical displacement of a vibrating string, say, one may alternatively set the displacement variables instantaneously to the values corresponding to a given set of traveling-wave components in the string model. Another benefit of the formulation is an exact physical interpretation of arbitrary initial conditions and excitations in the K-variable FDTD method. Since the DW formulation is exact in principle (though bandlimited), while the FDTD is approximate, even in principle, it can be argued that the true physical interpretation of the FDTD method is that given by the DW method. Since both methods generate the same evolution of state from a common starting point, they may only differ in computational expense, numerical sensitivity, and in the details of supplying initial conditions and boundary conditions.

The wave equation for the ideal vibrating string, reviewed in §C.1, can be written as

$\displaystyle Ky''= \epsilon {\ddot y},
$

where the following notation is used:

\begin{displaymath}\begin{array}{rclrcl} K& \isdef & \mbox{string tension} & \qq...
...isdef & \frac{\partial}{\partial x}y(t,x) \nonumber \end{array}\end{displaymath}    

In the following two subsections, we briefly recall finite difference and digital waveguide models for the ideal vibrating string.

Finite Difference Time Domain (FDTD) Scheme

As discussed in §C.2, we may use centered finite difference approximations (FDA) for the second-order partial derivatives in the wave equation to obtain a finite difference scheme for numerically integrating the ideal wave equation [481,311]:

$\displaystyle {\ddot y}(t,x)$ $\displaystyle \approx$ $\displaystyle \frac{y(t+T,x) - 2 y(t,x) + y(t-T,x) }{T^2}$ (E.1)
$\displaystyle y''(t,x)$ $\displaystyle \approx$ $\displaystyle \frac{y(t,x+X) - 2 y(t,x) + y(t,x-X) }{X^2}
\protect$ (E.2)

where $ T$ is the time sampling interval, and $ X$ is a spatial sampling interval.

Substituting the FDA into the wave equation, choosing $ X=cT$, where $ c \isdeftext \sqrt{K/\epsilon }$ is sound speed (normalized to $ c=1$ below), and sampling at times $ t=nT$ and positions $ x=mX$, we obtain the following explicit finite difference scheme for the string displacement:

$\displaystyle y(n+1,m) = y(n,m+1) + y(n,m-1) - y(n-1,m)$ (E.3)

where the sampling intervals $ T$ and $ X$ have been normalized to 1. To initialize the recursion at time $ n=0$, past values are needed for all $ m$ (all points along the string) at time instants $ n=-1$ and $ n=-2$. Then the string position may be computed for all $ m$ by Eq.$ \,$(E.3) for $ n=0,1,2,\ldots\,$. This has been called the FDTD or leapfrog finite difference scheme [127].


Digital Waveguide (DW) Scheme

We now derive the digital waveguide formulation by sampling the traveling-wave solution to the wave equation. It is easily checked that the lossless 1D wave equation $ Ky''=\epsilon {\ddot y}$ is solved by any string shape $ y$ which travels to the left or right with speed $ c \isdeftext \sqrt{K/\epsilon }$ [100]. Denote right-going traveling waves in general by $ y_r(t-x/c)$ and left-going traveling waves by $ y_l(t+x/c)$, where $ y_r$ and $ y_l$ are assumed twice-differentiable. Then, as is well known, the general class of solutions to the lossless, one-dimensional, second-order wave equation can be expressed as

$\displaystyle y(t,x) = y_r\left(t-\frac{x}{c}\right) + y_l\left(t+\frac{x}{c}\right). \protect$ (E.4)

Sampling these traveling-wave solutions yields
$\displaystyle y(nT,mX)$ $\displaystyle =$ $\displaystyle y_r(nT-mX/c) + y_l(nT+mX/c)$  
  $\displaystyle =$ $\displaystyle y_r[(n-m)T] + y_l[(n+m)T]$  
  $\displaystyle \isdef$ $\displaystyle y^{+}(n-m) + y^{-}(n+m)
\protect$ (E.5)

where a ``$ +$'' superscript denotes a ``right-going'' traveling-wave component, and ``$ -$'' denotes propagation to the ``left''. This notation is similar to that used for acoustic-tube modeling of speech [297].

Figure E.1: Digital simulation of the ideal, lossless waveguide with observation points at $ x=0$ and $ x=3X=3cT$.
\includegraphics[width=\twidth]{eps/fidealCopy}

Figure E.1 (initially given as Fig.C.3) shows a signal flow diagram for the computational model of Eq.$ \,$(E.5), termed a digital waveguide model (developed in detail in Appendix C). Recall that, by the sampling theorem, it is an exact model so long as the initial conditions and any ongoing additive excitations are bandlimited to less than half the temporal sampling rate $ f_s=1/T$ [451, Appendix G]. Recall also that the position along the string, $ x_m = mX= m cT$ meters, is laid out from left to right in the diagram, giving a physical interpretation to the horizontal direction in the diagram, even though spatial samples have been eliminated from explicit consideration. (The arguments of $ y^{+}$ and $ y^{-}$ have physical units of time.)

The left- and right-going traveling wave components are summed to produce a physical output according to

$\displaystyle y(nT,mX) = y^{+}(n-m) + y^{-}(n+m) \protect$ (E.6)

In Fig.E.1, ``transverse displacement outputs'' have been arbitrarily placed at $ x=0$ and $ x=3X$.


FDTD and DW Equivalence

The FDTD and DW recursions both compute time updates by forming fixed linear combinations of past state. As a result, each can be described in ``state-space form'' [449, Appendix G] by a constant matrix operator, the ``state transition matrix'', which multiplies the state vector at the current time to produce the state vector for the next time step. The FDTD operator propagates K variables while the DW operator propagates W variables. We may show equivalence by (1) defining a one-to-one transformation which will convert K variables to W variables or vice versa, and (2) showing that given any common initial state for both schemes, the state transition matrices compute the same next state in both cases.

The next section shows that the linear transformation from W to K variables,

$\displaystyle y(n,m) = y^{+}(n-m) + y^{-}(n+m), \protect$ (E.7)

for all $ n$ and $ m$, sets up a one-to-one linear transformation between the K and W variables. Assuming this holds, it only remains to be shown that the DW and FDTD schemes preserve mapping Eq.$ \,$(E.7) after a state transition from one time to the next. While this has been shown previously [442], we repeat the derivation here for completeness. We also provide a state-space analysis reaching the same conclusion in §E.4.

From Fig.E.1, it is clear that the DW scheme preserves mapping Eq.$ \,$(E.7) by definition. For the FDTD scheme, we expand the right-hand of Eq.$ \,$(E.3) using Eq.$ \,$(E.7) and verify that the left-hand side also satisfies the map, i.e., that $ y(n+1,m) = y^{+}(n+1-m) + y^{-}(n+1+m)$ holds:

\begin{eqnarray*}
y(n+1,m) &=& y(n,m+1) + y(n,m-1) - y(n-1,m) \\
&=& y^{+}(n-m...
...^{+}[(n+1)-m] + y^{-}[(n+1)+m] \\
&\isdef & y(n+1,m) \nonumber
\end{eqnarray*}

Since the DW method propagates sampled (bandlimited) solutions to the ideal wave equation without error, it follows that the FDTD method does the same thing, despite the relatively crude approximations made in Eq.$ \,$(E.2). In particular, it is known that the FDA introduces artificial damping when applied to first order partial derivatives arising in lumped, mass-spring systems [447].

The equivalence of the DW and FDTD state transitions extends readily to the DW mesh [518,447] which is essentially a lattice-work of DWs for simulating membranes and volumes. The equivalence is more important in higher dimensions because the FDTD formulation requires less computations per node than the DW approach in higher dimensions (see [33] for some quantitative comparisons).

Even in one dimension, the DW and finite-difference methods have unique advantages in particular situations [223], and as a result they are often combined together to form a hybrid traveling-wave/physical-variable simulation [351,352,222,124,123,224,263,33].


Next Section:
State Transformations
Previous Section:
Von Neumann Analysis