The Finite Difference Approximation

Free Books Physical Audio Signal Processing

In the musical acoustics literature, the normal method for creating a computational model from a differential equation is to apply the so-called finite difference approximation (FDA) in which differentiation is replaced by a finite difference (see Appendix D) [481,311]. For example

$\displaystyle {\dot y}(t,x)\approx \frac{y(t,x)-y(t-T,x)}{T} \protect$

(C.2)

and

$\displaystyle y'(t,x)\approx \frac{y(t,x)-y(t,x-X)}{X} \protect$

(C.3)

where

is the time sampling interval to be used in the simulation, and

is a spatial sampling interval. These approximations can be seen as arising directly from the definitions of the partial derivatives with respect to

and

. The approximations become exact in the limit as

and

approach zero. To avoid a delay error, the second-order finite-differences are defined with a compensating time shift:

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

(C.4)

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

(C.5)

The odd-order derivative approximations suffer a half-sample delay error while all even order cases can be compensated as above.

FDA of the Ideal String

Substituting the FDA into the wave equation gives

$\displaystyle K\frac{y(t,x+X) - 2 y(t,x) + y(t,x-X)}{X^2} = \epsilon \frac{y(t+T,x) - 2 y(t,x) + y(t-T,x)}{T^2}$

which can be solved to yield the following recursion for the string displacement:

$\displaystyle y(t+T,x)$	$\displaystyle =$	$\displaystyle \frac{KT^2}{\epsilon X^2} \left[ y(t,x+X) - 2 y(t,x) + y(t,x-X)\right]$
		$\displaystyle \qquad\qquad\qquad\qquad + 2 y(t,x) - y(t-T,x) \protect$

In a practical implementation, it is common to set $T=1,\, X=(\sqrt{K/\epsilon })T$ , and evaluate on the integers

and

to obtain the difference equation

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

(C.6)

Thus, to update the sampled string displacement, past values are needed for each point along the string at time instants

and

. Then the above recursion can be carried out for time

by iterating over all

along the string.

Perhaps surprisingly, it is shown in Appendix E that the above recursion is exact at the sample points in spite of the apparent crudeness of the finite difference approximation [442]. The FDA approach to numerical simulation was used by Pierre Ruiz in his work on vibrating strings [392], and it is still in use today [74,75].

When more terms are added to the wave equation, corresponding to complex losses and dispersion characteristics, more terms of the form appear in (C.6). These higher-order terms correspond to frequency-dependent losses and/or dispersion characteristics in the FDA. All linear differential equations with constant coefficients give rise to some linear, time-invariant discrete-time system via the FDA. A general subclass of the linear, time-invariant case giving rise to ``filtered waveguides'' is

$\displaystyle \sum_{k=0}^\infty \alpha_k \frac{\partial^k y(t,x)}{\partial t^k} = \sum_{l=0}^\infty \beta_l \frac{\partial^l y(t,x)}{\partial x^l},$

(C.7)

while the fully general linear, time-invariant 2D case is

$\displaystyle \sum_{k=0}^\infty \sum_{l=0}^\infty \alpha_{k,l} \frac{\partial^k... ...nfty \beta_{m,n} \frac{\partial^m\partial^n y(t,x)}{\partial x^m \partial x^n}.$