Page  378 ï~~Circulant Feedback Delay Networks for Sound Synthesis and Processing Davide Rocchesso Centro di Sonologia Computazionale Dipartimento di Elettronica e Informatica University degli Studi di Padova Email: Julius 0. Smith Center for Computer Research in Music and Acoustics Stanford University Email: Abstract A special lass of feedback delay networks using circulant matrices is proposed. These structures can be efficiently implemented and allow to control the time and frequency behavior. Several applications of circulant feedback delay networks in audio signal processing are presented. 1 Introduction In the seventies, M. A. Gerzon [19761 introduced feedback delay networks (FDN) as structures well suited for artificial reverberation. These structures are characterized by a set of delay lines connected in a feedback loop through a "feedback matrix" (see Fig. 1). J. Stautner and M. Puckette [1982] have utilized FDNs for multichannel reverberation, and more recently, J. M. Jot [1991, 1992] extensively studied FDNs and developed associated techniques for designing good reverberators. In this work, we explore some of the algebraic properties of FDNs and consider the use of circulant feedback matrices. We show how they can be used to reduce computational complexity and give good control over time-frequency behavior. An FDN having a circulant feedback matrix will be called a circulant feedback delay network (CFDN). In addition to reverberation applications, CFDNs have other uses in sound processing and synthesis. Examples include simulating radiating structures such as instrument bodies, feedback resonators, and digital delay effects in live electronic performance. 2 Basic Formulation An FDN is built using N delay lines, each having a length in seconds given by r = miT, where T = 1/F. is the sampling period. The complete FDN is given by the following relations: N y(n) = Zcisi(n) + dx(n) i=l N si(n + m) = Eai,jsj(n) + bix(n) j=1 (1) where si(n), 1 <i <N, are the delay-line outputs at time sample n. If mi = 1 for each i, we obtain the conventional state-variable description of a discrete-time linear system [Kailath, 1980]. In the present case, the variables si(n) form a subset of the system state at time n, with the remaining state variables being the samples contained within the delay lines at time n. Using the Z-transform, we can rewrite Eq. (1) in the frequency domain as Y(z) = cTS(z) +a dX(z) S(z) = D(z)[AS(z) + bX(z)] (2) where sT(z) = [81(z)... SN(Z)], bT = [b1... bN], cT = [cl... cN]. The diagonal matrix D(z) = diag (z-m1, z-m2,..-. zN) is called the "delay matrix" and A = [a,,j] is called the "feedback matrix". From Eq. (2), the transfer function is easily a2,1 a2,2 a2.3 a3,1 a3,2 a3,3.jL Figure 1: Order 3 Feedback delay network. Sound SynthesisTechniques 378 ICMC Proceedings 1994

Page  379 ï~~found to be H(z) = Y(z)/X(z) = cT[D(z-1) - A]-'b + d. (3) The system poles are the solutions of either det[A - D(z-')] = 0 (4) or det[zI -At]=0, z#0 (5) where At is the state-transition matrix of the linear system of Eq. (2). The matrix At is not uniquely determined by A, but we can express At in such a way that AtAtT = diag (1,n,l,... ImN i, AAT). The system zeros are the solutions of detA-bcT - D(z-')] =0 (6) The formulation of Eq. (1) represents a reference structure, in the sense that, with the appropriate choice of feedback matrix, it is a lossless structure. In practice, we must insert attenuation coefficients and filters in the feedback loop. In general, if one inserts a gain gi at the output of each delay line in the FDN, this corresponds to replacing D(z) with D(z/a) in Eq. (2), where g - = 's(7) In the companion paper we show the relationship between FDNs and waveguide networks [Smith and Rocchesso, 1994]. The feedback matrix can be seen as the scattering matrix of a generalized waveguide junction. We are interested in matrices associated with lossless junctions, because in these cases the signal power is conserved. The set of lossless matrices strictly contains that of unitary matrices. For instance, the scattering matrix of a generic lossless waveguide junction is not necessarily unitary. It is easy to check that A is unitary if and only if the state transition matrix At is unitary. It follows that if A is unitary, its eigenvalues are unit modulus, and therefore, all of the FDN poles are on the unit circle. Thus, choosing a unitary feedback matrix corresponds to choosing undamped, non-decaying eigenmodes for the FDN. Unitary matrices have nice properties in signal processing. In particular, they conserve signal power, even under time-varying conditions, which helps control dynamic range requirements in fixedpoint implementations [Smith, 1987]. 3 Circulant FDN Matrices Consider the class of circulant feedback matrix A, i.e., having the form a(O) a(1)... a(N-1)1 A- a(N-1) a(O)... a(N-2) a(1)... a(0) The following two facts can be proved [Davis, 1979]: Fact 1: If a matrix is circulant, it is normal (i.e., AA = AA). Fact 2:. If a matrix is circulant and has eigenvalues on the unit circle (i.e., it is lossless), then it is unitary. It is well known that every circulant matrix is diagonalized by the DFT matrix. This implies that the eigenvalues of A can be computed by means of the Discrete Fourier Transform of any row or column: {A(A)} = {A(k)} = DFT([a(0)... a(N - 1)]T) where {A(A)} denotes the set of all eigenvalues of A, and {A(k)} denotes the set of complex DFT samples obtained from taking the DFT of {a(n)}. A matrix that is both unitary and circulant has all eigenvalues on the unit circle, and the DFT can be used to compute the eigenvalue phases. In the case of equal-length delay lines, the eigenvalues provide exactly the resonance frequencies. Converselyy we can easily design the circulant matrix to have a desired distribution of eigenvalues. The actual presence of resonance peaks corresponding to eigenvalues depends on the positions of the zeros, as given by Eq. (6). If we are interested in having a maximally flat frequency response, it is natural to put the zeros exactly over the poles. This gives a perfectly flat response for equal-length delay lines, and an almost flat response at low frequencies for slightly different delay lengths. Incidentally, this distribution of zeros is very efficiently obtained by using d = 1, bT = [ 1 1... 1]7 andc having n entries equal to 1, nt entries equal to -1, and zeros for the remaining entries. This result is due to the following fact, which can be proved by Fourier diagonalization: Fact 3: Given a Circulant NxN matrix A, add a constant c to each entry of n rows (columns) and subtract the same constant c from each entry of another n rows (columns). The resultant matrix A' has the same eigenvalues as A. ICMC Proceedings 1994 379 Sound SynthesisTechniques

Page  380 ï~~4 Implementation In an N-th order FDN, the core computations consist of N updates of the delay lines and a matrix by vector multiplication. The delay line operations can proceed in parallel. The matrix by vector multiplication requires in general O(N2) operations (multiplications and additions). If the matrices arise from the scattering coefficients of a waveguide junction, the computations reduce to O(N). The same order of complexity is required by the normalized junction. For the special case of a junction of equal-impedance waveguides, the multiplications can be replaced by shifts when N is a power of 2. In all these efficient cases the eigenvalues of the feedback matrix are constrained to be at +1 or -1. The circulant matrix offers a more general. eigenvalue distribution. Moreover, the matrix by vector multiplication can be implemented very efficiently in hardware. This multiplication can be viewed as a circular convolution of the column vector with the first row of the matrix. Such a convolution can be performed, when N is a power of 2, using two FFTs (one of which can be precomputed), a dot product between two N-vectors, and an inverse FFT. The complexity of this algorithm is O(N log(N)). It is easy to implement this matrix-vector product in VLSI by means of the butterfly or other hypercubic architectures [Leighton, 1992]. These architectures allow computations of the FFT in 0(log(N)) time steps, and the algorithm can be pipelined. 5 Applications The circulant networks has been used for various purposes in sound synthesis and processing. The initial application was artificial reverberation, but other useful applications can be found. 5.1 Digital Reverberation A reverberator may be characterized by its frequency density and its time density [Jot, 1992]. We define the frequency density D1 as the number of resonances per Hertz. In our case, if the delay lines all have the same length m, assuming that all the poles are distinct and no cancellation occurs, we have an average frequency density given by Nm D! = -,(8) The time density Dt is defined as the number of pulses per second found in the impulse response. In actual rooms, Dt is an increasing function of time. In order to obtain dense reverberation after the first reflections (e.g., after 80 msec), we use different delay lengths. The actual positions of frequency peaks depend on the feedback matrix and on the delay lengths. If the delay lengths are fixed, we can vary the time-frequency properties of the structure simply by varying the distribution of eigenvalues of the feedback matrix. A uniform distribution of eigenvalues along the unit circle is optimum for the frequency response in the sense that it minimizes the maximum distance between peaks. However, it produces a highly repetitive time response. Conversely, clustering the eigenvalues around a point on the unit circle can be good for maximizing the length of time patterns, but the clustering of frequency peaks produces a poor reverberator amplitude response vs. frequency. We see from these considerations that there is a time-frequency tradeoff. This tradeoff can be addressed using circulant matrices. The shape of the frequency response depends also on the zeros, and the discussion leading up to Fact 3 considers ways to choose the vectors b and c to obtain a flat amplitude response at low frequencies. An FDN can be interpreted as a diffusing object in a reverberant environment. Consider a square room as in Fig. 2. In the center of the room imagine an object that is the sound source, the listener, and a diffusing element at the same time. A pictorial view could be a musician listening to himself and contributing with his body to the scattering of acoustic waves. Figure 2: FDN as reverberator Let us consider propagation of waves only in certain well-defined trajectories. These trajectories correspond to closed paths from the source to itself, each involving only a small number of reflections. Let N be the number of these differently shaped trajectories. In the interpretation, each closed ray path corresponds to a delay line. The vector of inputweighting coefficients b corresponds to the radiation strength of the sound source along each path. The feedback matrix corresponds to the diffusing characteristics of the object. For instance, the matrix element aid, corresponds to the energy transferred from pattern j to pattern i. The diagonal of the feedback matrix deter Sound SynthesisTechniques 380 ICMC Proceedings 1994

Page  381 ï~~mines the strength of "standing waves" which set up along each pattern. To model the frequency dependent nature of specular versus diffuse reflection, the matrix elements aij may be replaced by digital filters Ai j(z) in which case the diagonal elements would be gentle lowpass filters while the off-diagonal elements would tend to be highpass. It is possible to share lowpass filtering among all "standing waves" and to share highpass filtering among all "diffusing waves." This modeling point of view is limited by the fact that only N "closed ray paths" in the room are being simulated, and all non-specular reflections are being forced to enter some subset of the supported ray paths. A circulant matrix can be interpreted as a diffusing object where the transfer of energy between two patterns depends only on the relative position of the patterns (under a certain ordering). This physical interpretation is useful for computing the lengths of the delay lines, according to the geometrical proportions of the room and to the fact that each path corresponds to a normal mode. For this purpose, a three-dimensional model may be used. 5.2 Resonators FDN's with short delay lines may be used to produce resonances irregularly spread over frequency. A possible application could be the simulation of resonances in the body of a string instrument. M. Mathews and J. Kohut [1973] showed that in this kind of simulation, the exact position and height of resonances is not important; on the contrary, they stated that the Q's of the resonances must be sufficiently large and the- peaks must be sufficiently close together. With CFDNs we can easily achieve these goals, and we can vary the distribution of peaks by acting on the delay lengths and/or the feedback matrix. This means we have the possibility of spanning a large number of resonances using a small number of parameters. Another interesting application of CFDNs is as resonators in a feedback loop for pseudo-physical sound-synthesis techniques. By exciting these structures with bursts of white noise we obtain a multidimensional extension of the Karplus-Strong algorithm [Karplus and Strong, 1983], that is very effective for simulating membranes and bars. Alternatively, we can couple these resonators with nonlinear exciters and explore new families of sustained sounds. 6 Summary Feedback delay networks, and special cases using circulant matrices, have been discussed. We find that circulant networks provide good efficiency and versatility, and they can be used in a variety of applications. We continue to study the algebraic properties of circulant networks and to look for further practical applications. References [Davis, 1979] P. J. Davis. Circulant Matrices, John Wiley & Sons, New York, 1979. [Gerzon, 1976] M. A. Gerzon. Unitary (Energy Preserving) Multichannel Networks with Feedbacks, Electronics Letters V 12(11) 278 -279, 1976. [Jot and Chaigne, 1991] J. M. Jot, A. Chaigne. Digital Delay Networks for Designing Artificial Reverberators, 90th Convention of the Audio Engineering Society, Paris, Feb. 1991. [Jot, 1992] J. M. Jot. Etude et Redlisation d'un spatialisateur de son, PhD Thesis. TELECOM Paris 92 E 019, Paris, Sept. 1992. [Kailath, 1980] T. Kailath, Linear Systems, Prentice-Hall, 1980. [Karplus and Strong, 1983] K. Karplus, A. Strong. Digital Synthesis of Plucked-String and Drum Timbres, Computer Music Journal 7(2): 43-55, 1983. [Leighton, 1992] F. T. Leighton. Introduction to parallel algorithms and architectures: arrays, trees, hypercubes, M. Kaufmann Publishers, San Mateo, California, 1992. [Mathews and Kohut, 1973] M. Mathews, J. Kohut. Electronic Simulation of Violin Resonances, The Journal of the Acoustical Society of America 53(6): 1620-1626, 1973. [Perlis, 1952] S. Perlis. Theory of Matrices, Addison-Wesley, Reading, Mass.,1952. Dover Publications, New York 1991. [Smith, 1987] J. 0. Smith. Music Applications of Digital Waveguides, Report No. STAN-M-39, CCRMA - Stanford University, 1987. [Smith and Rocchesso, 1994] J. 0. Smith, D. Rocchesso. Connections between Feedback Delay Networks and Waveguide Networks for Digital Reverberation, International Computer Music Conference, Aarhus, Denmark, Sept. 1994. [Stautner and Puckette, 1982] J. Stautner, M. Puckette. Designing Multichannel Reverberators, Computer Music Journal 6(1): 52-65, 1982. ICMC Proceedings 1994 381 Sound SynthesisTechniques