State Space Formulation
In this section, we will summarize and extend the above discussion by means of a state space analysis [220].
FDTD State Space Model
Let
denote the FDTD state for one of the two subgrids at time
, as defined by Eq.
(E.10). The other subgrid is handled
identically and will not be considered explicitly. In fact, the other
subgrid can be dropped altogether to obtain a half-rate,
staggered grid scheme [55,147]. However, boundary
conditions and input signals will couple the subgrids, in general. To
land on the same subgrid after a state update, it is necessary to
advance time by two samples instead of one. The state-space model for
one subgrid of the FDTD model of the ideal string may then be written
as
To avoid the issue of boundary conditions for now, we will continue working with the infinitely long string. As a result, the state vector

When there is a general input signal vector
, it is necessary to
augment the input matrix
to accomodate contributions over both
time steps. This is because inputs to positions
at time
affect position
at time
. Henceforth, we assume
and
have been augmented in this way. Thus, if there are
input
signals
,
, driving the full
string state through weights
,
, the vector
is of dimension
:
![\begin{displaymath}
\underline{u}(n+2) =
\left[\!
\begin{array}{c}
\underline{\upsilon}(n+2)\\
\underline{\upsilon}(n+1)
\end{array}\!\right]
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4608.png)






forms the output signal as an arbitrary linear combination of
states. To obtain the usual displacement output for the subgrid,
is the matrix formed from the identity matrix by deleting every
other row, thereby retaining all displacement samples at time
and
discarding all displacement samples at time
in the state vector
:
![\begin{displaymath}
\underbrace{\left[\!
\begin{array}{c}
\vdots \\
y_{n,m-2} \...
..._{n,m+4}\\
\vdots
\end{array}\!\right]}_{\underline{x}_K(n)}
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4612.png)




The intra-grid state update for even

For odd




![$ [1,-1,1,-1,1]$](http://www.dsprelated.com/josimages_new/pasp/img4624.png)


![$ [1,-1,1]$](http://www.dsprelated.com/josimages_new/pasp/img4625.png)
![\begin{displaymath}
\underbrace{\left[\!
\begin{array}{l}
\qquad\vdots\\
y_{n+1...
...+3}\\
\qquad\vdots
\end{array}\!\right]}_{\underline{x}_K(n)}
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4626.png)


![\begin{displaymath}
\underbrace{\left[\!
\begin{array}{l}
\qquad\vdots\\
y_{n+1...
...ine{\upsilon}(n+1)
\end{array}\!\right]}_{\underline{u}(n+2)}.
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4627.png)
DW State Space Model
As discussed in §E.2, the traveling-wave decomposition
Eq.(E.4) defines a linear transformation Eq.
(E.10) from the DW
state to the FDTD state:
Since



Multiplying through Eq.


![]() |
![]() |
![]() |
|
![]() |
![]() |
![]() |
(E.30) |
where
To verify that the DW model derived in this manner is the computation diagrammed in Fig.E.2, we may write down the state transition matrix for one subgrid from the figure to obtain the permutation matrix

and displacement output matrix

![\begin{displaymath}
\underbrace{\left[\!
\begin{array}{c}
\vdots \\
y_{n,m-2} \...
...+4} \\
\quad\vdots
\end{array}\!\right]}_{\underline{x}_W(n)}
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4644.png)
DW Displacement Inputs
We define general DW inputs as follows:
![]() |
![]() |
![]() |
(E.33) |
![]() |
![]() |
![]() |
(E.34) |
The



![$ [y^{+}_{n+2,m},y^{-}_{n+2,m}]^T$](http://www.dsprelated.com/josimages_new/pasp/img4651.png)
![$ [\underline{\upsilon}(n+2)^T,\underline{\upsilon}(n+1)^T]^T$](http://www.dsprelated.com/josimages_new/pasp/img4652.png)
Typically, input signals are injected equally to the left and right along the string, in which case

![\begin{displaymath}
\underbrace{\left[\!
\begin{array}{c}
\vdots\\
y^{+}_{n+2,m...
...ine{\upsilon}(n+1)
\end{array}\!\right]}_{\underline{u}(n+2)}.
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4655.png)

![\begin{displaymath}
\underline{x}_W(n+2) = \mathbf{A}_W\underline{x}_W(n) +
\un...
...d{array}\!\right]}_{{\mathbf{B}_W}}
\underline{\upsilon}(n+2).
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4657.png)
To show that the directly obtained FDTD and DW state-space models
correspond to the same dynamic system, it remains to verify that
. It is somewhat easier to show that
![\begin{eqnarray*}
\mathbf{T}\,\mathbf{A}_W&=& \mathbf{A}_K\,\mathbf{T}\\
&=&
\l...
...dots & \vdots & \vdots & \vdots & \vdots
\end{array}\!\right].
\end{eqnarray*}](http://www.dsprelated.com/josimages_new/pasp/img4659.png)
A straightforward calculation verifies that the above identity holds,
as expected. One can similarly verify
, as expected.
The relation
provides a recipe for translating any
choice of input signals for the FDTD model to equivalent inputs for
the DW model, or vice versa.
For example, in the scalar input case (
), the DW input-weights
become FDTD input-weights
according to
![\begin{displaymath}
\left[\!
\begin{array}{l}
\qquad\vdots\\
y_{n+1,m-1}\\
y_{...
...psilon}(n+2)\\
\underline{\upsilon}(n+1)
\end{array}\!\right]
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4662.png)

Finally, when




![\begin{displaymath}
\mathbf{B}_K=
\left[\!
\begin{array}{cc}
\vdots & \vdots\\
...
...0 \\
2 & 0 \\
1 & 0 \\
\vdots & \vdots
\end{array}\!\right]
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4668.png)








DW Non-Displacement Inputs
Since a displacement input at position corresponds to
symmetrically exciting the right- and left-going traveling-wave
components
and
, it is of interest to understand what
it means to excite these components antisymmetrically. As
discussed in §E.3.3, an antisymmetric excitation of
traveling-wave components can be interpreted as a velocity
excitation. It was noted that localized velocity excitations in the
FDTD generally correspond to non-localized velocity excitations in the
DW, and that velocity in the DW is proportional to the spatial
derivative of the difference between the left-going and right-going
traveling displacement-wave components (see Eq.
(E.13)). More
generally, the antisymmetric component of displacement-wave excitation
can be expressed in terms of any wave variable which is linearly
independent relative to displacement, such as acceleration, slope,
force, momentum, and so on. Since the state space of a vibrating
string (and other mechanical systems) is traditionally taken to be
position and velocity, it is perhaps most natural to relate the
antisymmetric excitation component to velocity.
In practice, the simplest way to handle a velocity input in a
DW simulation is to first pass it through a first-order integrator of the
form
to convert it to a displacement input. By the equivalence of the DW and FDTD models, this works equally well for the FDTD model. However, in view of §E.3.3, this approach does not take full advantage of the ability of the FDTD scheme to provide localized velocity inputs for applications such as simulating a piano hammer strike. The FDTD provides such velocity inputs for ``free'' while the DW requires the external integrator Eq.

Note, by the way, that these ``integrals'' (both that done internally
by the FDTD and that done by Eq.(E.37)) are merely sums over
discrete time--not true integrals. As a result, they are exact only
at dc (and also trivially at
, where the output amplitude is
zero). Discrete sums can also be considered exact integrators for
impulse-train inputs--a point of view sometimes useful when
interpreting simulation results. For normal bandlimited signals,
discrete sums most accurately approximate integrals in a neighborhood
of dc. The KW-converter filter
has analogous
properties.
Input Locality
The DW state-space model is given in terms of the FDTD state-space
model by Eq.(E.31). The similarity transformation matrix
is
bidiagonal, so that
and
are both approximately
diagonal when the output is string displacement for all
. However,
since
given in Eq.
(E.11) is upper triangular, the input matrix
can replace sparse input matrices
with only
half-sparse
, unless successive columns of
are equally
weighted, as discussed in §E.3. We can say that local
K-variable excitations may correspond to non-local W-variable
excitations. From Eq.
(E.35) and Eq.
(E.36), we see that
displacement inputs are always local in both systems.
Therefore, local FDTD and non-local DW excitations can only occur when
a variable dual to displacement is being excited, such as velocity.
If the external integrator Eq.
(E.37) is used, all inputs are
ultimately displacement inputs, and the distinction disappears.
Boundary Conditions
The relations of the previous section do not hold exactly when the string length is finite. A finite-length string forces consideration of boundary conditions. In this section, we will introduce boundary conditions as perturbations of the state transition matrix. In addition, we will use the DW-FDTD equivalence to obtain physically well behaved boundary conditions for the FDTD method.
Consider an ideal vibrating string with
spatial samples. This is a sufficiently large number to make clear
most of the repeating patterns in the general case. Introducing
boundary conditions is most straightforward in the DW paradigm. We
therefore begin with the order 8 DW model, for which the state vector
(for the 0th subgrid) will be
![\begin{displaymath}
\underline{x}_W(n) =
\left[\!
\begin{array}{l}
y^{+}_{n,0}\...
...}_{n,4}\\
y^{+}_{n,6}\\
y^{-}_{n,6}\\
\end{array}\!\right].
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4678.png)
![\begin{displaymath}
\mathbf{C}_W=
\left[\!
\begin{array}{ccccccccccc}
1 & 1 & ...
...0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 1 & 1
\end{array}\!\right]
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4679.png)



![\begin{displaymath}
{\mathbf{B}_W}
=
\left[\!
\begin{array}{cc}
0 & 0 \\
0 & ...
...
0 & 0 \\
0 & 0 \\
0 & 0 \\
0 & 0
\end{array}\!\right]
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4681.png)


Resistive Terminations
Let's begin with simple ``resistive'' terminations at the string
endpoints, resulting in the reflection coefficient at each end of
the string, where
corresponds to nonnegative (passive)
termination resistances [447]. Inspection of
Fig.E.2 makes it clear that terminating the left endpoint may be
accomplished by setting


![]() |
(E.38) |
The simplest choice of state transformation matrix


![\begin{displaymath}
\mathbf{T}\isdef
\left[\!
\begin{array}{ccccccccccc}
1 & 1...
... 1 & 1 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 1
\end{array}\!\right]
\end{displaymath}](http://www.dsprelated.com/josimages_new/pasp/img4685.png)


![\begin{eqnarray*}
\tilde{\mathbf{A}}_K&\isdef & \mathbf{T}\tilde{\mathbf{A}}_W\m...
...r \\
0 & 0 & 0 & 0 & 0 & 0 & g_r & -g_r
\end{array}\!\right],
\end{eqnarray*}](http://www.dsprelated.com/josimages_new/pasp/img4687.png)
where
and
. We see that the left
FDTD termination is non-local for
, while the right
termination is local (to two adjacent spatial samples) for all
.
This can be viewed as a consequence of having ordered the FDTD state
variables as
instead of
. Choosing the other ordering
interchanges the endpoint behavior. Call these orderings Type I and
Type II, respectively. Then
; that is, the similarity
transformation matrix
is transposed when converting from Type I
to Type II or vice versa. By anechoically coupling a Type I FDTD
simulation on the right with a Type II simulation on the left,
general resistive terminations may be obtained on both ends which are
localized to two spatial samples.
In nearly all musical sound synthesis applications, at least one of
the string endpoints is modeled as rigidly clamped at the ``nut''.
Therefore, since the FDTD, as defined here, most naturally provides
a clamped endpoint on the left, with more general localized terminations
possible on the right, we will proceed with this case for simplicity in what
follows. Thus, we set and obtain
![\begin{eqnarray*}
\mbox{$\stackrel{{\scriptscriptstyle \vdash}}{\mathbf{A}}$}_K&...
..._r \\
0 & 0 & 0 & 0 & 0 & 0 & g_r & -g_r
\end{array}\!\right]
\end{eqnarray*}](http://www.dsprelated.com/josimages_new/pasp/img4695.png)
Boundary Conditions as Perturbations
To study the effect of boundary conditions on the state transition
matrices
and
, it is convenient to write the terminated
transition matrix as the sum of of the ``left-clamped'' case
(for which
) plus a series of one or more rank-one
perturbations. For example, introducing a right termination with
reflectance
can be written
where



In general, when is odd, adding
to
corresponds to a connection from left-going waves to
right-going waves, or vice versa (see Fig.E.2). When
is
odd and
is even, the connection flows from the right-going to the
left-going signal path, thus providing a termination (or partial
termination) on the right. Left terminations flow from the bottom to
the top rail in Fig.E.2, and in such connections
is even
and
is odd. The spatial sample numbers involved in the connection
are
and
, where
denotes the greatest integer less than or equal to
.
The rank-one perturbation of the DW transition matrix Eq.(E.39)
corresponds to the following rank-one perturbation of the FDTD
transition matrix
:




In general, we have
Thus, the general rule is that








![$ [1,-1,1,-1,\ldots]$](http://www.dsprelated.com/josimages_new/pasp/img4718.png)

Reactive Terminations
In typical string models for virtual musical instruments, the ``nut
end'' of the string is rigidly clamped while the ``bridge end'' is
terminated in a passive reflectance . The condition
for passivity of the reflectance is simply that its gain be bounded
by 1 at all frequencies [447]:
A very simple case, used, for example, in the Karplus-Strong plucked-string algorithm, is the two-point-average filter:





![\begin{eqnarray*}
\mbox{$\stackrel{{\scriptscriptstyle \vdash\!\!\dashv}}{\mathb...
... \\
0 & 0 & 0 & 0 & -1/2 & 1/2 & -1 & -1
\end{array}\!\right].
\end{eqnarray*}](http://www.dsprelated.com/josimages_new/pasp/img4722.png)
This gives the desired filter in a half-rate, staggered grid case. In the full-rate case, the termination filter is really


Another often-used string termination filter in digital waveguide models is specified by [447]
![\begin{eqnarray*}
s(n) &=& -g\left[\frac{h}{4}, \frac{1}{2}, \frac{h}{4}\right]\...
...{j\omega T})&=&
-e^{-j\omega T}g\frac{1 + h \cos(\omega T)}{2},
\end{eqnarray*}](http://www.dsprelated.com/josimages_new/pasp/img4724.png)
where is an overall gain factor that affects the decay
rate of all frequencies equally, while
controls the
relative decay rate of low-frequencies and high frequencies. An
advantage of this termination filter is that the delay is
always one sample, for all frequencies and for all parameter settings;
as a result, the tuning of the string is invariant with respect to
termination filtering. In this case, the perturbation is






![\begin{eqnarray*}
\mbox{$\stackrel{{\scriptscriptstyle \vdash\!\!\dashv}}{\mathb...
...d g_2 & \quad -g_2 & \quad g_3 & \quad -g_3
\end{array}\!\right]
\end{eqnarray*}](http://www.dsprelated.com/josimages_new/pasp/img4728.png)
where

The filtered termination examples of this section generalize
immediately to arbitrary finite-impulse response (FIR) termination
filters . Denote the impulse response of the termination filter
by
![$\displaystyle s(n)=[s_0,s_1,s_2,\ldots,s_N],
$](http://www.dsprelated.com/josimages_new/pasp/img4730.png)



Interior Scattering Junctions
A so-called Kelly-Lochbaum scattering junction [297,447] can be introduced into the string at the fourth sample by the following perturbation








![$ [-1,1]$](http://www.dsprelated.com/josimages_new/pasp/img1532.png)
A single time-varying scattering junction provides a reasonable model for plucking, striking, or bowing a string at a point. Several adjacent scattering junctions can model a distributed interaction, such as a piano hammer, finger, or finite-width bow spanning several string samples.
Note that scattering junctions separated by one spatial sample (as
typical in ``digital waveguide filters'' [447]) will
couple the formerly independent subgrids. If scattering junctions are
confined to one subgrid, they are separated by two samples of delay
instead of one, resulting in round-trip transfer functions of the form
(as occurs in the digital waveguide mesh). In the context of
a half-rate staggered-grid scheme, they can provide general IIR
filtering in the form of a ladder digital filter [297,447].
Lossy Vibration
The DW and FDTD state-space models are equivalent with respect to
lossy traveling-wave simulation. Figure E.4 shows the flow diagram
for the case of simple attenuation by per sample of wave
propagation, where
for a passive string.
The DW state update can be written in this case as





State Space Summary
We have seen that the DW and FDTD schemes correspond to state-space models which are related to each other by a simple change of coordinates (similarity transformation). It is well known that such systems exhibit the same transfer functions, have the same modes, and so on. In short, they are the same linear dynamic system. Differences may exist with respect to spatial locality of input signals, initial conditions, and boundary conditions.
State-space analysis was used to translate initial conditions and boundary conditions from one case to the other. Passive terminations in the DW paradigm were translated to passive terminations for the FDTD scheme, and FDTD excitations were translated to the DW case in order to interpret them physically.
Next Section:
Computational Complexity
Previous Section:
Excitation Examples