A Finite State Space Model for Representing the Broadband Infinite Series for Acoustic Reverberation between Parallel Reflecting Planes

The acoustic reverberation between two parallel reflecting planes can be represented by an infinite series of the images caused by the planes. To provide a more useful model for analysis and control, the infinite series version of the Green’s function is converted into a finite state space model that retains the high frequency character that enables broadband noise inputs to be examined.The infinite series is first summed into a very accurate, approximate closed form expression in the time domain in terms of a radical function. The radical is then transformed into an expression containing exponentials which have exact Laplace transforms that lead to an overall closed form transfer function for the system.The system transfer function is transformed into a third-order state spacemodel that theoretically contains all of the frequency characteristics of the infinite series representation. The accuracy of the state space model is examined by comparing it to the infinite series solution for three typical types of acoustical inputs: exponential for impulse noise, single frequency sine for harmonic noise, and a multifrequency Schroeder phased harmonic sequence for random noise.


Introduction
Finite state space models for systems have several important advantages.These include both single-input single-output and multi-input multi-output systems which can be formally treated.It is best suited both for the theoretical treatment of control systems (optimization, output feedback, and estimation with observers) and for numerical simulations.The determination of the system initial condition response is straightforward.It gives a better insight into the inner system behavior such as controllability, observability, and system characteristics.The use of finite state space models to represent partially or fully bounded acoustic and/or structural vibration systems for investigating lower frequency behavior has previously received attention.Several noted examples for truncated state space models used for lower frequencies are included below.Five modes from second-order modal differential equations along a bounded string are used in pole placement output and its derivative for state space feedback control design via an observer [1].A few secondorder modal differential equations for waves in a bounded one-dimensional acoustic duct are coupled with speaker dynamics in a model to investigate low frequency behavior [2][3][4].In [4], output feedback control using H ∞ and lag compensation are examined.A discrete model of acoustic waves in a bounded one-dimensional duct using five sections of the duct and ten state equations is presented in [5].Finite element modeling is used in [6] to develop over twenty-four thousand coupled second-order modal differential equations for a two-dimensional bounded structure coupled to a threedimensional bounded acoustic space.A state space model using 214 states is then used to investigate frequencies up to 400 Hz.Finite element modeling is used in [7] to develop many coupled second-order modal differential equations for a one-dimensional bounded structure coupled to the several electromagnetic piezoelectric patches.A state space model using a finite number of states is then used to investigate the performance of optimal LQ controller frequencies for up to 300 Hz.Most of the above studies are targeted for lower frequencies.
The use of the Euler-Maclaurin sum formula to provide an approximate, accurate closed form solution for the infinite series Green's function that describes partially or fully bounded acoustic and structural systems is developed and verified in [8][9][10] by comparing the closed form solution to the infinite series solution for several types of inputs.In [8], the main part and associated variations of the Euler-Maclaurin closed form with the infinite series solutions in a partially bounded three-dimensional acoustic space between two parallel plates are presented for an exponential input.In [9], the Euler-Maclaurin closed form solution is compared with the infinite series solutions in a partially bounded two-dimensional space acoustic space for a narrow duct.The comparison includes both the time and frequency domains for an exponential input.It also includes a broad range (up to 1000 Hz) transfer function based comparison of the Euler-Maclaurin closed form with the infinite series for feed forward control of a disturbance and closed loop output feedback control.In [10], Euler-Maclaurin closed form representations for the Green's functions in both onedimensional acoustic and structural systems are determined in Laplace transform domain in terms of time integrals.The Laplace transform integrals are solved numerically and used in a broad range (up to 1000 Hz.) transfer function based feed forward control system to meet a desired output for an input disturbance.
For all of the Euler-Maclaurin based studies in [8][9][10], state space representations were not considered.In view of the Euler-Maclaurin closed form results given in [8][9][10], in both time domain response and the transfer function based control for frequencies up to 1000 Hz, and the successes in [1][2][3][4][5][6][7] for the truncated state space models for various types of low frequency control in acoustic and structural systems, the work in [11] focuses on the development of a low order but full finite state space model that would naturally include a broadband frequency range and essentially all of the terms required for the convergent infinite series in the Green's function representation of bounded structures.Also, in [11], a one-dimensional structural system of a string is used to develop a suitable state space for a two-dimensional plate.The developed state space model in [11] is shown to give an accurate trend for representing of the series transfer function for high frequencies up to 2000 Hz. and the accurate time response for both a direct problem with an exponential input and for use in an observer based inverse problem of estimating a disturbance input.
In this paper, the ideas and approach in [11] for structural systems are expanded, modified, and improved to develop a low order but full finite space model that would naturally include a broadband frequency range and essentially of all the terms required for a convergent infinite series in the Green's function representation of the acoustic field between reflecting planes commonly known as reverberation.The distinguishing difference between this paper and the development in [11] is that, in [11], the Euler-Maclaurin sum formula is applied to an infinite series of string modes.In this paper, the Euler-Maclaurin sum formula is applied directly to the Laplace transform form of the Green's function infinite series for acoustic reverberation in [8] providing all polynomial type functions in the Laplace transform domain which are directly combined into a system function.This system transfer function in [11] leads directly to an input/output ordinary differential equation which serves as the basis for defining the state variables and determining the state equations.For the acoustic space of the current paper, the Euler-Maclaurin sum formula must be directly applied to the time domain form of the Green's function infinite series as is done in [8].As previously indicated, that direct application of the Euler-Maclaurin sum formula to the time domain form of the Green's function infinite series [8] leads to an integral of a radical function that has no exact analytical solution for providing polynomial type functions in the Laplace transform domain.Thus, in this paper, this radical function is first approximated with a combination of two exponential types to provide a suitable simple system overall analytical transfer function for proceeding to the input/output ordinary differential equation that is used to develop the state space model.Typical noise inputs of impulse, harmonic, and random each including high frequency character up to 3000 Hz. are used to compare the performance state space model with the performance of many terms in the series model.Two extreme numerical cases of plates separation  = 9 m.(30 ft.) for a typical industrial area to  = 0.3 m. (1 ft.) for a typical tight-fitting enclosure such as for a diesel engine system are presented numerically.

Euler-Maclaurin Representation.
A schematic of the diametrical cross section in the z-r plane for a three-dimensional acoustic space is shown in Figure 1.
The acoustic field between two perfectly reflecting parallel planes can be represented as an infinite sum of images in the Laplace transform domain (see [12,13]).

4𝜋𝐺 (𝑟, 𝑧
where H is the distance between the planes.
z o and z are the source and receiver position between the planes.r is the distance √ 2 +  2 in the x, y plane from source to receive parallel to the planes.
c is the speed of sound (345 m/s, 1130 f/s).
G(r, z, s) is the Laplace transform of the Green's function g(x, y, z, t) for the wave equation In [14], the form of the Euler-Maclaurin sum formula that transforms an infinite series to a closed form expression with a remainder R is used here.It is given by In [8] the Euler-Maclaurin sum formula in (3) is applied to the four infinite series in (1), the inverse Laplace is applied, and the following closed form solution of the Green's function in the time domain is obtained.
where  is the Heaviside step function and O is an order at most.
The Green's function is composed of discrete delta functions  and a group of four dispersion-like terms with Heaviside step function  that come from the integral in the Euler-Maclaurin sum formula.For purposes of showing the concept of a state space representation, (4) is simplified by using  =   = 0 without the loss of generality.The Green's function in (4) is expressed as (, ) = 2 and is given by where and the remainders  are considered small enough to be neglected as shown in [8].
The main reverberation term is the third term in (6) and monotonically decreases from a maximum of   = /2 2 at  =   to 0 at  = ∞.It is conceivable to reason that the sum of several to many exponentials would provide a reasonably accurate fit over a reasonable range of t.Each exponential would increase the system order by one and thus several to many significant states would result.However, the goal of this approximate model is to provide a more simple compact form with just two exponential terms while still providing a significantly influential model, albeit with a lower number of states.Reasonable for this acoustic system may be defined as providing good enough accuracy for use in numerical simulations for various type of inputs, for use in input estimates (see [11]), and for use in some control applications.In a frequency description of reverberant sound in a large industrial space, it is shown in [15] how a function composed of two exponentials, one with constant coefficient and one with time varying coefficient, could be used to accurately represent multiple model floor and ceiling reflections.A similar model is used here and is given by where parameters , , and  are determined to provide a good fit over the time range from   to some suitable value .
Goodness of fit and suitability are defined here as providing a reasonably accurate comparison of the behavior of the developed state space model with the infinite series of image in (1) for a range of practical input sources.The input sources used here are an impact type rapidly decaying exponential source, a single high frequency harmonic input, and a broadband input generated with a modified multifrequency Schroeder phased harmonic sequence as shown in [11].
A four-part procedure is given for determining a, b, and B: (1) and at  =   , ( For the value of a f rom part 1, define The maximim value of   occurs for (3) max is determined numerically and the last two equations in (9) give 1), (2), and (3) give a good set of values for , , and .These values are then used as the initial set for an optimization procedure to obtain a better set of values.The optimization function to minimize is given by for a time  suitable for covering the primary part of the Green function response.
To demonstrate this, two numerical cases of the ,  set representing somewhat of an extreme range are presented.The first is for  = 3 m.(10 ft.),  = 9 m.(30 ft.) which may apply to a large work space or open plan off area and the second is for  = 0.3 m. (1 ft.),  = 0.9 m. (3 ft.) which may apply to the inside of a tight-fitting machine enclosure.Both cases are selected so that the ratio of  to  is the same value of 3.
Table 1 gives the results for the two cases.Units for the table are based on  and  in ft.
Figures 2 and 3 give a comparison of the two sides of (7) for cases 1 and 2, respectively.In the figures, the left hand side of (7) (which comes from the closed form representation using the Euler-Maclaurin sum formula on the infinite series) is defined as yem, the first term in the right side of ( 7) is defined as y1, the second term in the right side of ( 7) is defined as y2, and the total approximation is defined as yapp = y1 + y2.The horizontal axes in the figures are the time  =  −   .The figures clearly show the good accuracy for the most important early part of the approximation with still reasonable accuracy for the latter part of the transformation.The trend is consistent through the entire time.The match in the time scale of the right hand plots is less because the exponential approximation selected for the simple state space model contains two basic terms.Additional accuracy for longer time could be achieved for the right hand time frame by using more exponential type terms but would defeat the purpose of having a simple model.Practically all of the audible sound heard occurs before this latter time frame.The exponential input results are discussed with respect to accuracy for early and later parts of the time response relative to the peak response.Peak response is used as the baseline for errors because of the importance near the peak.

State Space Representation.
Starting with the Green's function of ( 6) with the right side approximation of ( 7) and using () to represent (, ), Applying the Laplace Transform () = ∫ ∞ 0 (, ) −  gives Cross multiplying by To simplify notation, the delta functions with time shifts are defined as Applying the inverse Laplace transform gives a third-order differential equation for (). ...
The input, output equation is obtained by substituting output  for () and input  for .Define the coefficients as The input, output form of ( 17) is given by ...
Three states are defined: where   ,   ,  1 , and  1 are coefficients to be determined in terms of system parameters , , , , , , and   .
The first state equation is determined by taking the derivative of the first state definition and substituting ẏ from the second state definition.It is given by The second state equation is determined by first taking the derivative of the second state definition which brings a ÿ term.
Next substituting ẋ 2 from the third state definition into (22) and solving for ÿ gives The second state equation is then given by The third state equation is determined by first taking the derivative of the third state definition and using the derivative of (22) which brings in a ...  term.
and equating it to ...  from the input, output form of (19) gives , and ÿ from above and collecting terms gives where 's and 's are defined in (18).
To determine the coefficients   ,   ,  1 , and  1 used in the state definitions, like terms in (27) are combined and the coefficients of the six , u , and ü terms on each side of the equation are set equal to give and   are determined by the third equation in (28).
1 and  1 are determined from the first or second equation combined with the third in (28) and are given by Substituting from (28) and ( 29) into (26) gives the third state equation Mathematical Problems in Engineering The finite state space differential equation is given by The system in (31) contains delays in the input.The output  of interest is the first state  1 .The state matrix A has a canonical form with the third-order characteristic equation  () = sin () (34)

Broadband Input.
To simulate a broadband type of input that can provide a response over a typical range of noise frequencies from 100 to 3000 Hz, a multifrequency Schroeder phased harmonic sequence used in [11] with random amplitudes and variable phases is given by where   = 2(100 Hz),   = 30 to provide an input from 100 Hz to 3000 Hz and 0 ≤   ≤ 1 via a uniformly distributed random number generator.This input is randomlike in time (as used in [11]) and provides a suitable model for a broadband multifrequency input for evaluation in this paper.

Numerical Comparison of Series and State Space Model.
The numerical solution to the state space model of differential equations ( 31) is determined with a special fourth-order Runge-Kutta method modified to efficiently handle the delay inputs.This solution is compared to the solution from a direct sum of three hundred terms in a convolution of the infinite series in (1) with the input.From (1), the series representation of Green's function is given by Convolution with the input () gives the output y(t) for comparison Using the sampling property of the delta function, the series output is given by For computations, 300 terms in the sum provide a very accurate series solution for comparing the state space solution.
For the Euler-Maclaurin closed form approximation, convolution of input  and the Green function in (6) give Using the sampling property of the delta function, the Euler-Maclaurin closed form output is given by where the integral in (40) is solved numerically for each time increment.
Both systems in Table 1 are used as relevant examples comparing the state space output to the series output for each of the three types of inputs.
Figures 4 and 5 give the comparison for the exponential input.Input exponential parameters  = 1 and  = 5 are used to show a range for a typical slow impulse noise and a typical fast impulse noise, respectively.In the figures, the "state space" output is the first state for the solution to the finite state space model of (31) which is obtained from the Euler-Maclaurin closed form Green function.The "closed form" output is from (40) which is obtained from direct integration of the Euler-Maclaurin closed form solution.The "series" output is from (38), a summation of the infinite series for  = 300 values.For the case of  = 3 m.(10 ft.) in Figure 4, the state space model is very accurate.For the case of  = 0.3 m. (1 ft.) in Figure 5, the state space model is very accurate for the early part of transient and follows an accurate decaying trend with values a little lower than the series solution.For both cases, the larger value of  gives the best accuracy.This indicates that since larger values of  tend to simulate more of true impulse response which defines natural system behavior, the state space model may be well suited for control applications.The closed form model provides similar accuracy as the state space for both cases and is included to show how well the state space Euler-Maclaurin model agrees with the state space model derived from it.Table 2 gives a numerical assessment of the accuracy in Figures 4 and 5 for key parts of the response.For an exponential input used in to model for a practical impact type sound input, the important parts of the transient response are near the ascent from zero to the maximum value, the maximum value at the time, and the descent from the maximum value to a time that is several times the peak time  peak.As the time increases to this area and beyond, the effect  from this type of impact input on a system becomes much less important.Many practical noise sources are essentially a series of impacts which are space apart enough such that a periodic frequency associated with it is low, usually below about 300 Hz.The noise level and character one audibly hears or measures is usually broadband and primarily comes from a single impact each time it hits and is amplified by how many hits per second fit into a specified time frame for an equivalent time average level.A common example is impact noise in industrial environments from various and many sources of rigid body contacts occurring frequently and randomly causing the familiar audible sound of myriad of superimposed distinct bangs.The effects of these many similar impacts per second provide a broadband mid to high frequency spectrum that are the main contributors to the overall dBA level.Another example is the interaction of cooling system propeller fan blades with still air and/or guide vanes.Generally the fundamental blade passage harmonic frequency is below 300 Hz, dominates the level of higher harmonics, and is only a small contributor to the overall dBA noise.The primary contributor is the turbulent impulses that occur at each interaction.As the industrial impact noise sounding like a series of bangs, each pulse is audible as broadband turbulence.The effects of these many similar impulses per second provide a broadband mid to high frequency spectrum that are the main contributors to the overall dBA level.There are other examples such as combustion engine pressure cycling that may fall into this class.For this class of noise sources, the audible sound is characteristically that of a particular occurrence of the impact/impulse heard over and over again.
To reflect the importance of the response at and near maximum, the percent error between state space and series responses that is used in Table 2  The percent errors for the exponential input are mostly in the 1% to 7% range with a few higher ones in a 12% to 26% range occurring for case 2 during the descent from the peak.This level of error may be considered reasonable for the intent of the model.One reason for large errors for case 2 during the 2 nd descent may be the use of the simple twoterm approximation for yapp which leads to it being a little less accurate for higher time areas as seen in Figures 3 and  4.However, the important monotonically decreasing trend is still present in this higher time range.
Figures 6 and 7 give the comparison for a harmonic input of typical acoustic range frequencies 1000 Hz. and 2000 Hz.In the figures, the "state space" output is the first state for the solution to the finite state space model of (31) which is obtained from the Euler-Maclaurin closed form Green function.The "series" output is from (38), a summation of the infinite series for  = 300 values.For both  value cases, the amplitudes of the two models are off by about 50%, while the phases between the two models are much closer.for the key part of the exponential input.One reason for larger amplitude errors may be the fact that harmonic inputs can potentially cause standing waves to be more prevalent between the two reflecting plane, especially for case 2 ( = .3m) where the wavelengths (.34 m for 1000 Hz and .17m for 2000 Hz) are near full and half wavelengths of .Also noted is that the amplitude errors for case 1 are about one-half of those for case 2, while the phase errors for case 2 are about one-half of those for case 1.
Figure 8 shows a typical input trace () for the multifrequency Schroeder phased harmonic sequence for   = 30 frequencies in a frequency range from   = 100 Hz. to     = 3000 Hz.The randomness in time is seen primarily via peak amplitude and pulse width.The randomness in frequency is less with a tendency of the highest peaks in a range from 300 Hz. to 500 Hz. and the lower peaks in in a range from 2000 Hz. to 2500 Hz.These ranges are indicated on the figure.
Figure 9 gives the comparison for multifrequency Schroeder phased harmonic sequence.In the figures, the "state space" output is the first state for the solution to the finite state space model of (31) which is obtained from the Euler-Maclaurin closed form Green function.The "series" output is from (38), a summation of the infinite series for  = 300 values.For the case of  = 3 m.(10 ft.), there is very close agreement in that an equivalent dominant frequency of approximately 1100 Hz from eleven pulses over the 0.01 seconds shown is evident in both the Series and EMstate responses.A similar comparison is seen for the case of  = 0.3 m. (1 ft.) where there is very close agreement in that an equivalent dominant frequency of approximately 2600 Hz from twenty-six pulses over the 0.01 seconds shown is evident in both the Series and EMstate responses.For both cases, the series amplitudes are larger than the EMstate amplitudes,  while the EMstate pulse widths are larger than the series pulse widths.Using a combination of like frequency in the comparison as well as similar area under a typical pulse to provide a measure that can be expected to significantly affect sound level and character, the accuracy of the EMstate model may be considered reasonable for its intended uses as described previously.The work in [11] uses a similar state space derived model for estimating the random input to a plate via one output measurement and an observer based on the state model.This also may be a suitable area of application for the acoustic state model.Table 4 shows a numerical comparison associated with the measures above.
The comparison of the two primary indicators of dominant frequency (which has essentially a 0% error between Series and EMstate) and the pulse area (of which the EMstate is close to 70% of the series) may be considered to provide a reasonable amount of accuracy for the intended use of the EM state model.
The following summary of the results presented in the previous comparison curves and tables provides clarity of the scope and the applicability of the state space model for analysis and control.For the exponential input as a model of impact noise, the model may be considered accurate enough for use in both the analysis and control of single and multiple impact sources of the type discussed where the behavior up to, at, and just after the peak response is primarily responsible for the audible noise.For the harmonic input, the model may be considered reasonably accurate for analysis to estimate the response amplitude for larger plate separation distances.However, for the harmonic input, even though phase errors are less than amplitude errors, the response has only a single phase and the level of this phase error shown is considered too high for providing confidence in   the application of common control schemes to reduce levels from pure tone noise sources.For the Schroeder random like input model, the reasonable accuracy obtained for both response frequency occurrence and pulse area indicates that the state space model may be used with a reasonable amount of confidence in the analysis of overall noise level and character.Now since the Schroeder is not exactly random but composed of a series of pure tones with various phases and amplitudes, its makeup gives preference to selected frequencies as clearly shown in Figure 8.But Figure 9 shows that the response has varying phase errors on the order of magnitude as the harmonic phase.In the realm of noise issues such as in manufacturing, transportation, and ventilation, single impacts, repeated impacts, and broadband random type acoustic sources are generally much more prevalent than pure tone sources.The numerical results from the exponential and the Schroeder input representations of these types of sources indicate that the scope of the closed form Euler-Maclaurin state space model is reasonably accurate for analysis and control in these areas.The next section gives a frequency response method used by the author in [9] for testing how a Euler-Maclaurin approximate closed form system model can be tested against the corresponding series solution in a feedback control scheme for a purely random input disturbance of white noise which may be used to represent a broadband noise source.

Frequency Response
Figure 10 is a comparison of the magnitude of the transfer function for input () to output () obtained by substituting  =  in (1) and the transfer function for the Euler-Maclaurin state space model is determined by applying the Laplace transform to (31) and obtaining an expression for the transfer function.
The overall magnitudes across frequencies are comparable.The differences in shape comparison may be expected from the nature of the two models.The series is expected to contain many sharp narrow bandwidth spikes of varying amplitude and frequency as it is obtained from an infinite sum all of the images between two rigid plates.The EM state model is obtained from a closed form expression via the Euler-Maclaurin sum formula and derived using the two exponential approximate expression in (7).The finite EM state transfer function result indicates that the many frequencies due to images are consolidated into several frequencies of the same magnitude but with bandwidths increasing with frequency, indicating that the EM state model gives a look of more damping than the series model.
To test the use of the EM state model for feedback control of a broadband random noise source, a feedback control scheme is shown in Figure 11, where   () is the noise source input considered as a disturbance.
It is a model based output feedback control scheme requiring a measurement of the system output.This method was tested in [9] with a white noise disturbance and was shown to be useful for reducing white noise amplitudes in a two-dimensional acoustical space with an EM transfer function obtained directly from the closed form EM approximate function derived in that paper.The same method is applied here using an EM transfer function obtained from the Euler-Maclaurin state space model of (31) and (33).It is used to determine the controller transfer function required  For both Figures 12 and 13, the control drives the highly variable magnitude of the series transfer function over frequency to the desired closed loop magnitude with very approximate model's use in analysis and control.The model gives reasonably good accuracy and confidence for the use of a simple finite state model in state space noise analysis and control applications for impact and random type noise which are prevalent in practice.For a pure tone noise, the accuracy may be reasonably good for approximate analysis of output magnitude but the phase differences in the outputs of series and state models are significant enough to provide concern for its use in controlling single frequency sources.

Figure 1 :
Figure 1: Schematic of a diametrical cross section of the acoustic space.

Figure 11 :
Figure 11: Feedback control scheme for broadband noise disturbance.

Figure 12
is for |  ()| = 0.06 and Figure 13 is for |  ()| = 0.02.Additionally, the figures give the magnitude of the controller transfer function (see (45)) for each |  ()| value to demonstrate the gain and character of the controller magnitude |  ()|.

Table 3
The amplitude errors for the harmonic inputs are much larger than those for the exponential input's amplitude errors

Table 2 :
Numerical assessment of exponential inputs results in Figures4 and 5.