Overview of geometrical room acoustic modeling techniques

Computerized room acoustics modeling has been practiced for almost 50 years up to date. These modeling techniques play an important role in room acoustic design nowadays, often including aurali-zation, but can also help in the construction of virtual environments for such applications as computer games, cognitive research, and training. This overview describes the main principles, landmarks in the development, and state-of-the-art for techniques that are based on geometrical acoustics principles. A focus is given to their capabilities to model the different aspects of sound propagation: specular vs diffuse reﬂections, and diffraction. V C 2015 Author(s). All article content, except where otherwise noted, is licensed under a Creative Commons Attribution 3.0 Unported License . [http:/


I. INTRODUCTION
The general motivation for room acoustic modeling is to enable the construction of acoustically better environments.This is especially important in concert halls and other acoustically challenging spaces, such as theaters and studios, but it is also important in classrooms, railway stations, and other public venues and even in homes.Modeling provides acoustic attributes that characterize the acoustic conditions or enable the possibility to listen to the acoustics of the space under design, which is called auralization. 1 The acoustics of a room can be modeled under several different frameworks and for various purposes.The two main approaches are based either on numerically solving the wave equation or on the assumptions of geometrical acoustics (GA).In principle, wave-based modeling is able to provide the most accurate results.However, these techniques are computationally very expensive; thus, it is often more appropriate to resort to faster but less accurate techniques such as those based on GA.In GA, all of the wave properties of sound are neglected, and sound is assumed to propagate as rays.This assumption is valid at high frequencies, where the wavelength of sound is short compared to surface dimensions and the overall dimensions of the space, but at lower frequencies the approximation errors increase as wave phenomena play a larger role.
In practice, GA techniques are widely used in modeling the mid-and high-frequency behaviors of rooms.This has been the case for decades; Maa had already noted in 1941 that rigorous physical acoustics is insisted upon by "research theorists" but that GA is favored among practitioners. 2 Moreover, he stated this preference long before the use of computers for room acoustic simulations was discussed.

A. Scope of the overview
The goal of this overview is to present the history and state-of-the-art of all the main geometrical room acoustic modeling techniques.Our focus is on techniques that can be used to simulate sound propagation inside three-dimensional (3D) spaces to produce impulse responses and/or timeenergy responses.Time-wise, this means that we concentrate on developments after the publication of the landmark paper by Krokstad et al. in 1968, 3 which utilized ray tracing for computing time-energy responses and showed its applicability to practical room acoustic design.However, this overview is not limited to only room acoustic design but also covers GA modeling techniques developed for other areas such as virtual reality applications, cognitive research, and training.From algorithmic and computational viewpoints, these techniques are similar, although the requirements for performance and accuracy are different.
The modeling of structure-borne sound is excluded from this overview because the related noise propagation models are typically quite different from the techniques discussed here.Similarly, the computation of room acoustic attributes and the actual sound rendering component of auralization 1 and related signal processing techniques, such as artificial reverberation 4 or binaural processing, 5 are beyond the scope of this overview.We attempt to provide a comprehensive overview, specifically of geometrical acoustic-based methods, within the limits of a journal paper.The focus will be on the theoretical aspects of the various methods, which means that a wealth of practical aspects regarding input data acquisition, data uncertainty, and use of the techniques, among others, is omitted.In particular, the focus is on the capabilities of the techniques.These capabilities establish the baseline for achievable accuracy, but the final accuracy and ability to match measurement results depend on a sufficiently large number of factors such that a detailed a) Electronic mail: Lauri.Savioja@aalto.fidiscussion on accuracy is beyond the scope of this overview.Similarly, a detailed comparison of computational performance of different techniques has been left out of this article. 6his overview is organized as follows.Section II presents the early history of room acoustic modeling, starting from the pre-computer era and subsequently introducing the first steps toward computerized modeling.In addition, various landmarks in the adoption of GA modeling techniques in practice are discussed.The required conceptual models for sound propagation in air and for reflections from surfaces are introduced in Sec.III.Similarly, a mathematical framework for all GA modeling techniques is presented.There are several different methods of categorizing these techniques, but we divide them into two separate groups, path-and surfacebased techniques.Section IV presents techniques in which the actual sound propagation paths play a major role, whereas in Sec.V, we focus on techniques that propagate the sound energy to surfaces as an intermediate step and then gather the response at the listener from those surfaces.In Sec.VI, we summarize different modeling techniques and compare their capabilities.

II. BACKGROUND A. Basics of geometrical room acoustics
The main objective of room acoustic simulations is to provide time-energy and/or impulse responses than can be used to compute various attributes such as reverberation time, T 60 , and clarity, C 50 .A more challenging application of the modeling results is auralization, 1 which is a term that was coined by Kleiner et al. in 1993.By their definition, auralization is the process of rendering audible, by physical or mathematical modelling, the sound field of a source in a space, in such a way as to simulate the binaural listening experience at a given position in the modeled space.The basic principles of such auralization systems were previously introduced by Schroeder et al. in 1962 (Ref.7) and revisited again by Schroeder in 1970. 8 typical room impulse response, as measured in a room, or obtained by room acoustical simulation, is illustrated in Fig. 1(a).A corresponding time-energy response from a simulation might look like as in Fig. 1(b).Often, the time-energy response is computed as a first step and an impulse response is synthesized from the former.Such responses can be considered to consist of three separate parts: direct sound, early reflections, and late reverberation.Each peak in the figure represents a single reflection path starting from the direct sound that has not undergone any reflection.In the early part of the response most of the reflection paths can be seen separately, although they might overlap temporally and sum together.The reflections that occur in the room can be roughly categorized either as specular or diffuse.A diffuse reflection path can propagate in any direction from the reflecting surface while the specular reflection is as from a mirror.In the response presented in Fig. 1(b), the early specular reflections were modeled up to second order.
In practice, the sound field typically becomes quite diffuse soon after the early reflections for two reasons.The first reason is that reflections tend to scatter the sound instead of being ideally specular.The second reason is that the reflection density increases such that individual reflections are no longer important.This means that while the sound field starts as a superposition of distinguishable individual reflections it turns gradually more and more diffuse.Thus, any reflection path will eventually become statistical instead of representing an ideal specular reflection path. 9From the modeling viewpoint, this means that a technique should be able to address diffuse reflections in addition to specular ones to provide realistic predictions.This was previously observed in the first Round Robin comparison of room acoustic modeling tools, in which the techniques with diffuse reflection models performed better than did the other techniques. 10o find more detailed information on GA, one very good source is the Room Acoustics book by Kuttruff. 11

B. Early history
Sabine had previously used ray-based acoustics in the early 1900s to investigate sound propagation paths using Schlieren photography. 12,13Their impressive visualizations show wavefronts that are augmented with rays that are perpendicular to the wavefronts.Such images from physical scale models can be thought of as the first step toward computerized ray tracing.The concept of a mean free path is also closely linked to GA because it measures the average distance that a ray will travel between successive reflections in a given geometry.This measure serves as one of the underlying terms in the early reverberation time formulae developed by Sabine and Eyring. 14The validity of the GA approach when the wavelength of sound is small compared to the room dimensions was shown in 1929. 15The first computer model used to determine mean free path lengths was presented by Allred and Newhouse in 1958. 16Atal and Schroeder demonstrated the computation of global decay functions in two-dimensional (2D) rooms using ray tracing in 1967, 8,17 while the first method of computing room acoustic responses for different receiver positions was introduced by Krokstad et al. in 1968. 3

C. GA in practical use and earlier reviews
The first review article on room acoustic modeling was published by Schroeder in 1973. 18In his paper, he reports on concert hall auralization that was used to assist in a renovation and states that computers will change the field of concert hall design.At that time, computers were so limited that he deemed all wave-based modeling techniques impossible and that GA was the only practical approach to room acoustic modeling.In general, the deployment of computerized room acoustic modeling has been quite slow.In 1981, Gilbert stated that computing the reverberation time of a hemispherical dome would cost "several dollars" and that "A more realistic auditorium shape… would be feasible but it might cost hundreds of dollars per case." 19However, Krokstad et al. published a summary paper in 1983 in which they described several case studies conducted with their ray tracing system during the preceding 15 years. 20In addition, in 1989, Ondet and Barbry reported ray tracing simulations for noise maps in workshops for which the results were in good agreement with measurements. 21In 1991, Kuttruff published a review paper 22 in which he still saw the practical application of room acoustic modeling, both auralization and the determination of acoustic attributes, to be in the future.
A boom in room acoustic modeling research occurred in the early 1990s, as illustrated in a Special Issue of Applied Acoustics in 1993. 23Most current commercial modeling tools, such as ODEON (Ref.24), http://www.odeon.dk,CATT-Acoustic (Ref.25), http://www.catt.se,and EASE (Ref.26), http://ease.afmg.eu,date back to that period of time.Since then, numerous scientific papers have introduced improvements to the modeling techniques as well as case studies that aim to validate such techniques.The three Round Robins in which the results given by different modeling tools were compared to each other and against measurement results that took place in 1995, 10 2000, 27 29 The use of computers in room acoustic design has become increasingly popular over the years; however, even today, scale models are constructed to support computer simulations in certain concert hall projects.In addition to room acoustic design, computer games have started to utilize these modeling techniques to enable more realistic audio.However, a typical goal in game audio is to efficiently achieve plausible auralization, whereas in room acoustic design, accuracy is of primary concern, and computational performance is a secondary concern.
The only textbook that focuses on room acoustics modeling and auralization is Auralization: Fundamentals of Acoustics, Modelling, Simulation, Algorithms, and Acoustic Virtual Reality by Vorl€ ander. 5One of the most recent reviews was published by the same author in 2013, in which he focused on uncertainties in input data and how they affect simulations. 30A literature review focused on computer graphics techniques applied to room acoustic modeling was published in 2014 by Charalampous and Michael. 31

D. GA and computer graphics
The principles of GA are very general, and similar modeling techniques, especially ray tracing and radiosity, are employed in numerous other fields such as optics and underwater acoustics.However, such principles are most prominently used in computer graphics, where the analogous problem is known as global illumination, wherein the behavior of light is computationally simulated.Note that the first paper on ray tracing 32 in the field of computer graphics dates back to the same year as the paper by Krokstad et al. 3 was published.Although these two research areas have mostly been investigated separately, various clear overlaps are apparent; in addition, a certain amount of convergence has been observed recently.This has resulted in publication of some results in high-performance room acoustic modeling in computer graphics forums; see, e.g., the work by Raghuvanshi et al. 33 One commonality between all GA techniques is that these techniques require a substantial number of geometric computations, especially when finding intersections between rays and geometric primitives in a 3D space.The general approach to speed up such operations is to use various spatial data structures, 34 such as bounding volume hierarchies, kdtrees, and binary space partitioning, and they are frequently utilized in modern acoustic modeling systems.However, such techniques are not in the scope of this overview although they are actively utilized in room acoustic modeling as well, see, e.g., Refs.35 and 36.An interested reader is encouraged to seek more information in the computer graphics literature. 34,37

III. FUNDAMENTAL MODELS OF GA
This section presents the fundamental concepts underlying all GA modeling techniques.We start by discussing how sound propagation in air can be modeled and subsequently present several models for a single reflection for various types of surfaces.Combining these two aspects and applying them consecutively yields a complete model of GA.This model can be mathematically formulated as an integral equation, thereby providing a framework that can be used to describe any GA modeling technique.Several models of edge diffraction can be used to extend the validity of geometrical acoustic techniques; these models are described at the end of this section.
A. Sound propagation in air, pressure vs energy Geometrical acoustic methods usually address the energy and propagation delay of each sound wave as well as the visibility of the contribution as observed at the receiver.For a point source, the contribution is where the visibility factor V is 1 when the source-receiver path is unobstructed and 0 for an obstructed path, and the separation distance r leads to a propagation delay t ¼ r=c.In energy-based modeling techniques, such contributions are summed in discretized time slots, while the frequency dependence is addressed by studying separate bands.As further described below, reflections will introduce additional factors, and dissipation in air can be introduced as well.An essential property of these techniques is that they ignore interference effects resulting from this energetic summation.This approach leads to an echogram, or energy-time response, p 2 ðtÞ, which suffices for finding reverberation times and other energy-based acoustic parameters.
In pressure-based modeling techniques, the phase of the sound pressure is restored by constructing a complex sound pressure amplitude, p, for each contribution of the wave numbers k, whereby interference effects between sound waves are introduced because in these techniques, summations can be performed over sound pressures instead of energies.By considering a complex pressure, phase shifts at reflections can also be introduced if one has access to the complex wall impedances.Conversion from an impulse response pðtÞ to an energetic response p 2 ðtÞ is performed in a straightforward manner by squaring the impulse response.Performing the conversion in reverse is not trivial because the energetic response contains less information than does the corresponding impulse response.However, the impulse response is often needed for auralization.To address this issue, Kuttruff, for example, presented an algorithm for such a conversion. 38inally, the term GA is also used interchangeably with the approach where the contribution from a monopole source is 39 This latter approach is obviously equivalent to the phaserestoration technique described above.

B. Models of reflection
An interface between two media will modify the incident sound field in various ways.An infinite surface generates reflected and transmitted waves.A finite object is usually said to generate a scattered wave field, and for the special case of a finite polyhedral object, the scattered sound field might be decomposed into reflections and diffracted waves.Finally, a non-smooth flat surface is typically said to generate a reflected wave and scattered components.

Infinite, flat surfaces
For the simplest case of an infinite, flat, and smooth surface that is impinged upon by an incident wave, reflected and transmitted waves are generated.From our current perspective, only the reflected wave is of interest, and we can safely assume that the medium of the incident wave is air.The second medium at the interface is either of the local reaction or extended reaction category.The local reaction model is a simplification in which the second medium is assumed to be such that the interaction between the two media is a one-way interaction.This means that the sound waves from the second medium never enter the first medium, whereas in the extended reaction model, the interaction can be more complicated.The local reaction model is frequently used in room acoustics 11 and is applicable in GA even with porous materials that fail the assumptions of local reaction.
An explicit expression for the reflected wave is given by the image source (IS) solution, whereby the IS is obtained by reflecting the sound source against the surface. 40(Note that here, the IS solution refers only to a single reflection, whereas later in Sec.IV A, the IS technique is an entire room acoustic modeling technique.)The IS solution is exact for spherical wave incidence and a rigid surface.The IS solution is also asymptotically correct for high frequencies, spherical wave incidence, and any local reaction surface impedance.This asymptotic validity is employed in GA, where the IS solution is used in terms of the equivalent plane-wave reflec- where Z s is the impedance of the surface, h is the incident angle, and qc denotes the impedance of air.The plane-wave reflection coefficient is accurate as long as the source, or receiver, is several wavelengths away from the infinite reflecting plane. 41An exact solution can be calculated using the complex equivalent source method, 42,43 where a semiinfinite, linear distribution of ISs represents the reflecting surface.However, in practice, the plane-wave reflection coefficient is typically sufficient 44 and suitable for pressure-based GA techniques.In particular, the normal incidence reflection coefficients, for h ¼ 0 , are the data that are provided by impedance tube measurements of materials.
For the more commonly applied energy-based GA techniques, it is customary to use the absorption coefficient a aðhÞ ¼ 1 À jRðhÞj 2 ; (5 where the phase information of the reflection process obviously is lost.The angle-dependent absorption coefficient, aðhÞ, or the random incidence absorption coefficient, a rand , averaged across incidence angles h is the quantity that is used by various implementations of the GA techniques.
However, data made available by surface material manufacturers is typically measured according to the reverberation room standard, with octave-band or third-octave-band resolution for the range 100 Hz-5 kHz, 45 and systematic differences between such measurements and the true a rand has long been a topic of study. 46Rindel suggested a method for estimating the angle dependence of the reflection factor from the random incidence a s by assuming a real-valued impedance, 47 and an evaluation of the use of impedance vs absorption coefficient boundary conditions was presented by Jeong. 46

Curved surfaces
Concave surfaces often produce problems in room acoustic applications due to focusing.The modeling of such surfaces requires special attention in GA methods because the IS solution presented above holds only for flat surfaces.A typical procedure in GA techniques is to split the curved surface into small polygons and apply a flat surface reflection model to each polygon.This approach is sufficiently accurate in areas that are not close to the possible focal points; however, the accurate modeling of the sound pressure at the focus requires the use of wave-based modeling. 48n certain GA methods, such as ray tracing, there is no need for the explicit tessellation of curved surfaces; instead, an IS solution can be locally applied at each point of the surface.In practice, this requires that the surface is represented in such a form that the intersection of a ray and the surface, and the corresponding reflection, can be efficiently computed.In the remainder of this paper, curved surfaces are not treated separately from flat surfaces; the applied reflection models are assumed to be sufficiently general to consider methods that directly use tessellated piecewise-flat approximations.

Non-smooth surfaces
Real-world surfaces are not ideally flat and smooth; they include various irregularities.Irregularities that are substantially smaller than the wavelength can be neglected because they do not affect the reflection, whereas irregularities that are much larger than the wavelength can be thought of as curved surfaces or even as a locally flat surface.However, irregularities whose dimensions are of the same order of magnitude as the wavelength of sound can have a large effect on the reflected sound.This means that a single absorption coefficient or complex impedance is insufficient to describe the acoustic properties of such a surface or material, even if the coefficients are frequency dependent.Instead, more complicated measures are needed in modeling such materials.
The scattering from corrugated walls was previously studied by Rayleigh. 49The behavior of irregular surfaces was further investigated in the 1950s by Twersky (see, e.g., Ref. 50) and in the 1950s-1960s by Biot (see, e.g., Ref. 51).From our current viewpoint, major results were published by Biot in 1968 (Ref.51) as Biot's model is very general and addresses several types of surface irregularities.The reflection can be described by a set of radiating dipoles.This provides an elegant addition to the IS solution described above.In practice, this means that a reflection from a nonsmooth surface can spread in time, whereas in the IS solution, an impulse remains an impulse in a reflection.Siltanen et al. 52 have described how to apply this reflection model in practice with beam-tracing.Their results show that although the values of the room acoustic attributes remain within justnoticeable-differences, the model results in audible differences in auralization.
a. Scattering and diffusion coefficients.Instead of using the complete model by Biot, 51 it is customary to apply substantially simpler approximations.The scattering coefficient is the parameter that is often used to describe the roughness of a surface.The coefficient is defined as the ratio of nonspecularly reflected energy and total reflected energy, and it can be used to simply divide the reflected energy into two components, specular and diffuse, whereby a value of 1 corresponds to a fully diffuse reflection in which there is no specular reflection component.It can be noted that a value of 1 tells nothing of the distribution of this diffusely reflected sound.In contrast, a scattering coefficient with a value of 0 corresponds to an ideally specular reflection.The scattering coefficient can be measured 53 and is often applied in GA modeling techniques.The diffuse component is typically interpreted by assuming Lambertian reflection in which the reflected diffuse energy per solid angle of each ray is proportional to the cosine of the angle of reflection, thus having the same radiance when viewed from any angle.Lambertian reflection stems from the study of light reflection off rough surfaces and is consequently a high-frequency model, which has nevertheless found wide usage in acoustics.
The scattering coefficient is not to be confused with the diffusion coefficient, 54 which describes the irregularity of the reflection pattern, whereby a value of 1 corresponds to a completely uniform reflection distribution, in contrast to the common Lambert model, and a value of 0 is obtained when all of the energy is emitted in a single direction.Current modeling techniques typically employ only the scattering coefficient because no clear interpretation of how the diffusion coefficient should be implemented exists.Instead, proper implementation requires more information than simply one plain coefficient.This data can be represented by the bidirectional reflectance distribution functions (BRDFs), as described in the following.
b. BRDFs.Another aspect of complicated reflection characteristics is given by BRDFs that can be applied in GA techniques.Their origin is in the field of optics, 55 and they are frequently used in computer graphics.These functions provide the ratio of incoming and outgoing energy as a function of incoming and outgoing angles for a given frequency; however, in acoustics, a more complete form would include a temporal dimension as well, although it is often neglected in practice.In addition, rotational invariance can be typically assumed, thus reducing the problem by one dimension such that only the elevation angle of the incident sound is considered, and the azimuthal angle can be neglected.This leaves us with a 3D function for a given frequency such that the outgoing energy at a certain elevation and azimuthal angle depends only on the elevation angle of the incoming energy.That reflection model can be used in several different room acoustic modeling techniques. 56pecular reflection results in a delta-like peak in the function such that all of the energy is reflected in one direction, while the function takes on a value of 0 elsewhere.In contrast, the function used to model ideal diffuse reflection is a constant over the entire hemisphere of reflection.In practice, acoustic BRDFs are more complex and contain elements of both extremes.However, no such data are available that would enable the construction of full BRDFs for typical materials.Instead, a typical BRDF is a mixture of ideal specular and ideal diffuse reflection.Figure 2(a) illustrates an ideal BRDF including these two components, and Fig. 2(b) shows a typical discrete representation that transforms the ideally specular component into a more diffuse component.Having ideally specular components in the BRDF necessitates special consideration because it would produce an infinitely narrow beam that can induce numerical challenges.

Finite surfaces
In Secs.III B 1-III B 3, various types of infinite surfaces were discussed; however, for a finite surface, additional effects must be introduced.For a free-hanging finite reflector with no absorption or surface scattering, the reflected sound will have a high-pass filter effect, as well as interference effects above the high-pass cutoff-frequency, and asymptotically tend toward the infinite-wall reflection factor.Rindel suggested a simplified model for this high-pass filter effect based on the Kirchhoff approximation (KA). 57Two equivalent but contrasting views can be used to describe the effect of a finite surface.The reflection strength can be modified according to various models such as the above-mentioned model. 57Alternatively, diffraction waves emanating from the edges, and of opposite polarity, can be added to the reflection and produce the same high-pass filter effect as well as interference effects.The latter approach allows the angles to connecting planes to be considered.Diffraction modeling is further described in Sec.III D. Note that a weak specular reflection at low frequencies could appear as a corresponding large amount of scattering in all directions.However, for a finite surface, the weak reflection corresponds to the incident sound wave largely continuing past the finite reflector and exhibiting minimal scatter in other directions.

C. Model of geometrical room acoustics-Room acoustic rendering equation
A general model of geometrical room acoustics can be elegantly written as an integral equation.The first such formulation dates back to 1971 and is known as Kuttruff's integral equation. 9,11,58This version is still limited to diffuse reflections in a convex room, but it was extended by Joyce in 1978 to address arbitrary reflection functions as well. 59The basic principle of presenting such reflection characteristics is the same as that of the above-mentioned BRDF model, but it is called the reflection matrix.The diffuse reflection variant of the equation was further refined by Nosal et al. 60 to contain a source term that enables the computation of timeenergy responses, whereas the previous versions were only used to study the energy balance in the room.
Seemingly independent, although much later, activities concerning the development of a very similar equation having its roots in computer graphics have been observed.This equation is a rendering equation that covers most of the algorithms that attempt to compute global illumination in a scene. 61The first version of such an equation for room acoustics was presented by Tsingos and Gascuel 62 The applied terminology is quite different than that in Kuttruff's integral equation, although they are essentially the same, with the exception that this new version includes a visibility term that allows the equation to also address non-convex rooms.The main difference between the computer graphics and room acoustic versions is that the acoustic equation can incorporate the temporal dimension and finite propagation speed of sound, whereas in computer graphics, light is typically assumed to propagate instantaneously; thus, the light from a source instantly reaches all surfaces in a space.
A thorough derivation and presentation of this GA model and related equations was presented by Alarcão. 63In this work, the arbitrary reflection characteristics are called wall reflection functions instead of BRDFs, which is a more common term in the literature.The room acoustic rendering equation 56 presented below is yet another formulation of the model.This equation provides a framework with which most of the GA-based room acoustic modeling techniques can be expressed within the limits of local reaction.Below, we present an energy-based GA formulation of this equation. 64This formulation can be generated based on pressure as well, but we use energies because most of the models used in practice are based on energy and, especially, because the general BRDF reflection formulation presented in Sec.III B 3 is practical only for energies.The equation assumes that the underlying geometry G is represented as surface patches, FIG. 2. The BRDF describes how sound energy is reflected in different directions.In (a), an idealized BRDF is composed of a hemisphere to represent diffuse reflection and of a Dirac-like pulse for specular reflection.(b) A practical discretization in which the specular reflection is spread over a given solid angle.
and the formulation recursively calculates the outgoing timedependent sound energy l at a given surface point x 0 in a direction X: where l 0 is the energy emitted by the surface itself; the integral describes the reflected energy whereby C is the outgoing angle from the surface point x toward x 0 , as illustrated in Fig. 3; and R is the reflection kernel: where Vðx; x 0 Þ is the visibility term (see Sec. III A) and qðx 0 ; H; XÞ is the BRDF (see Sec. III B 3) that gives the fraction of energy [¼ 1 À a, see Eq. ( 5)] that is reflected from an incident angle H in the direction X.The geometry term gðx; x 0 Þ considers the propagation delay, air absorption, and form factor, which is defined by the relative orientation, and mutual distance of the surfaces.The form factor between small elements represented by the points x and x 0 in Fig. 3 can be written as A more complete formula considering the finite area of the surfaces has been presented, e.g., by Lewers. 65That version is practical if Eq. ( 6) is to be evaluated only on coarse points.
Most of the terms in the above equations are frequency independent, and only the BRDF and air absorption are functions of frequency.This means that most of the modeling techniques described in this overview can simultaneously address all of the frequency bands by registering energies in each frequency band in parallel.However, in certain techniques, the energy propagation paths might diverge for different frequencies.In such cases, the room acoustic simulation needs to be performed separately for each band.
The main differences between various GA modeling techniques lie in how they sample the integral in Eq. ( 6) and what type of BRDFs they support.
It is worth noting that other general versions of this equation known as the transport equation for sound energy, 63 or the acoustic radiative transfer equation, 66 have been derived.Whereas Eq. ( 6) is based on a boundary formulation, these equations apply a volumetric approach and are even able to address scattering media, which is caused, e.g., by air turbulence; such a treatment is not possible with Eq. ( 6).However, Eq. ( 6) is slightly simpler but still sufficiently expressive to cover the GA modeling techniques discussed in this overview.In certain cases, the models based on a simplified version of the acoustic radiative transfer equation, which is known as the diffusion equation, are considered as GA methods, but these methods are beyond the scope of our discussion because they remain fundamentally different than the methods studied in this survey.The essential difference between the methods is that the diffusion equation methods use a volumetric grid in which the sound energy is propagated 67,68 instead of the ray-like approach common to other GA methods.

D. Diffraction models
The diffraction phenomenon is not handled by classical GA methods.It is, however, possible to modify GA methods, or to add diffraction components to a GA-based solution and thereby take diffraction into account to some degree.Diffraction modeling is a substantial topic in itself so below, a brief overview is given of approaches that are relevant for combination with GA.Rooms are usually modeled as polygonal structures, and one approach to addressing diffraction is to add diffraction waves that appear to emanate from the edges to the geometrical acoustic solution with specular reflections.This decomposition of the sound field originates from the canonical case of an infinite rigid wedge, for which Biot and Tolstoy presented a time-domain solution. 69Later developments have used a secondary edge source approach, which permits the study of finite edges and higher-order diffraction. 70These secondary source models rely on integrations along each edge of the polygons, and multiple diffraction requires second-and higher-order integration.
Alternative diffraction models employ asymptotic formulations wherein an infinite edge's contribution (the diffraction wave emanating from the edge) is described by an explicit expression and therefore, these models require substantially fewer computations than does the secondary source approach.These models include the Geometrical Theory of Diffraction (GTD) 71 and the Uniform Theory of Diffraction (UTD). 72Another high-frequency asymptotic solution is the time-domain secondary source approach by Vanderkooy, 73 which can address polyhedra but produces substantial inaccuracies at low frequencies.
A separate class of diffraction modeling is based on the KA, which leads to the Kirchhoff-Helmholtz Integral Equation being solved very efficiently, as in the Fresnel-Kirchhoff approximation. 74,75A closely related approximation is the Huygens-Fresnel approximation. 76,77As early as the late 1800s, it was shown that the KA could be used to convert the exact Kirchhoff-Helmholtz integral equation into a line integral along the edges of a finite polygonal reflector, FIG. 3. The form factor between two surface patches is defined by their relative angles C and H and by their mutual separation distance.Overall, the room acoustic rendering equation gives the energy emitted in the X direction from point x 0 and integrates the reflected energy over all of the reflecting surfaces represented by points x.
which is the so-called Maggi-Rubinowicz transformation, and that such a line integral also leads to a secondary edge source formulation. 78This approach has been used in room acoustical modeling; 79 however, for certain geometries, this KA-based diffraction approach is not asymptotically correct for higher frequencies 80 and has inherent limitations at low frequencies.
Other diffraction methods include the asymptotic expression by Pierce, which can also be applied to thick noise barriers. 81 very different and computationally efficient approach to diffraction modeling combines ray tracing with an "uncertainty principle" approach and was first presented by Stephenson in 1996. 82,83Rays that pass near an edge are diverted, and the diversion is stronger the closer to the edge the ray passes.This approach implies that the response behind a thin screen quite closely follows the band-averaged diffraction pattern rather than the discontinuous pattern due to pure GA.The original method was derived for a single edge, or slit, but a 3D version has been presented by Schr€ oder and Pohl. 84,85In addition, Pohl has solved the problems related to reciprocity in the original formulation. 85 similar method was suggested by Benedetto and Spagnolo in 1984 86 and in optics in 1971. 87umerous combinations of the above-mentioned diffraction models and the GA modeling techniques have been proposed, as will be discussed in more detail later, especially in Sec.IV D.

E. Source and receiver models
There are both sound sources and a receiver in addition to the acoustic space along a complete acoustic transmission path.However, the treatment of sources and receivers in this overview is limited and is discussed only to a degree that is relevant to actual room acoustic modeling.From this viewpoint, it is important to note that real sources and receivers have certain directivity characteristics.The propagation of sound from a directive source can be express by extending Eq. (1) as follows: where the directivity factor (DF) is generally a function of the two radiation angles h and u.
On the one hand, standardized measurements of room acoustic properties should be performed with omnidirectional loudspeakers; 88 therefore, omnidirectionality is typically used in the computation of acoustic attributes.0][91][92] In the following, we assume that such directivity patterns are available, and their treatment in different modeling techniques is presented.It should, however, be noted that the complexity of the radiation from music instruments might be impossible to describe by fixed directivity functions. 89eceivers in a room acoustic scenario are microphones, each with a given directivity/sensitivity, or human listeners.
In the current scope, we can assume that for both types, it is sufficient to obtain directional information of incoming acoustic energy.This information can then be used for postprocessing in any type of spatial audio reproduction technology.

IV. REFLECTION-PATH-BASED GEOMETRICAL ACOUSTIC MODELING TECHNIQUES
The techniques presented in this section attempt to calculate the room responses by determining both specular and diffuse reflection paths.The final response is determined by combining the contributions of all these paths as well as those of possible diffracted paths.Most of the methods in this group are energy-based methods.Exceptions to this are the methods that address only specular reflections and, possibly, edge diffraction, which can also be computed using pressure-based methods.

A. Image-source method
The main principles of the image method, currently called either the mirror image-source method or, most often, the image-source method, were first developed in the study of electrostatics, and early work in acoustics was presented by Carslaw in 1899, 93 which furthered previous work by Sommerfeld.These early studies addressed canonical problems, such as an infinite wedge, and their extension to the shoebox-shaped room appears to have been first presented by Cremer in 1948 (Ref.94) and later by Mintzer in 1950 (Ref.95) for calculating sound pressure in a rectangular room as a function of time.Gibbs and Jones were the first to use computers to determine the image-source positions and sound pressure levels in a rectangular space. 96Some years later, Allen and Berkley presented a similar system and published a FORTRAN code for the computation of impulse responses with the image-source technique. 97The elegance of this method is in its simplicity, and it still provides an exact solution to the wave equation for a parallelepipedical room with rigid walls.In such a space, all reflections are ideally specular, and there is no diffraction because there are only right-angle corners.Aretz et al. demonstrated that the image-source method provides accurate predictions even at low frequencies for more realistic boundary conditions when the complex angular-dependent reflection coefficients of Eq. ( 4) are utilized. 98A similar, but not as detailed, study was conducted earlier by Lam. 44he main principle of the image-source method is to recursively apply the image-source solution, explained in Sec.III B 1, to the sound source.In practice, this means that a sound source is reflected against all surfaces in a model, thus resulting in a set of ISs; then, these ISs are again reflected against all the surfaces.This process is repeated until a termination condition, such as response length or reflection order, is satisfied.The resulting ISs can be considered as secondary sources that each represent a reflection such that the distance from the IS to a receiver corresponds to the actual reflection path length inside the room.Figure 4 illustrates this process in a 2D setup in which there are four first-order reflections, eight second-order unique reflections (six of which are shown), etc.The corresponding 3D case would include 6 and 18 reflections.In a general case, the number of ISs up to the Kth-order reflections is P K k¼1 NðN À 1Þ kÀ1 for N surfaces.However, in the case of a rectangular room, certain ISs are degenerated, and thus, the actual numbers are smaller than those determined using this expression.
The result of the image-source computation can be observed as a hierarchical image-source tree that includes the sound source as the root, and each branch represents an IS.The impulse response can be constructed as the sum of the contributions of all ISs in this tree.The directivities of sound sources can be easily considered by reflecting the directivity patterns of the sources in addition to the source locations. 99he speed of image-source computations in a rectangular room was increased by McGovern. 100In this technique, all of the redundant and unnecessary computations are eliminated using sorted look-up tables.The obtained result is as accurate as the original method, but the technique was reported to be almost nine times faster.

Approximate extensions of the image-source method in rectangular spaces
In room acoustic design, a rectangular geometry is rarely sufficient.However, this simple model is frequently applied in studies wherein room responses are needed to support other research. 101For example, in speech processing, signals are often processed to contain artificial reverberation to make them sound more natural; in this application, very simple room models are applicable.However, such responses still need to be of a length that corresponds to the reverberation time of the space.To this end, various techniques can be used to extrapolate the simulation results of early reflections to cover the entire impulse response duration and avoid summations per IS.Kristiansen et al. presented an extrapolation technique to obtain longer responses at the cost of the method being approximate, which is unlike the original accurate image-source method. 102In 2007, Duraiswami et al. suggested the use of the multipole technique to enable efficient computations. 184This technique can be applied in scenarios in which the sound source or the listener moves or in which there are several listeners.For a single source-receiver combination, the performance remains the same.The diffuse reverberation approach by Lehmann and Johansson uses the standard image-source technique to determine the energy decay curve of the space and then synthesizes the diffuse component based on that curve. 101,103The listening tests showed that the majority of subjects were not able to differentiate between full imagesource auralization and the proposed diffuse reverberation model.However, these models are purely artificial and do not correspond to any realistic spaces.Modeling such rooms requires more accurate techniques such as those described in the following.

Image-source method in arbitrary geometries
The concept of addressing arbitrary polygonal structures with the image-source model was mentioned as early as 1976 by Santon. 104However, this major improvement to the image-source method was first fully described by Borish in 1984. 105The recursive reflection principle is the same as that of rectangular geometries, but it requires various additional checks.First, ISs should be constructed only for sources and ISs in front of the reflecting side of polygons, whereas all of the backsides should be discarded as reflecting surfaces.Similarly, if the reflecting surface is completely behind the previous reflector, no new IS needs to be created.All of these computations can be performed independently from the listener position, and thus, the resulting imagesource tree is valid for the entire space.However, the second test concerns the visibility of an IS and requires information regarding the location of the listener.In this test, a specular reflection path is formed from the source to the listener using the ISs, as illustrated in Fig. 5.This path must hit all of the reflecting surfaces inside their boundaries, and simultaneously, the path may not intersect any other surface in the model; otherwise, the IS is not visible at the given listener position.Figure 6 illustrates the resulting IS tree up to second order in this 2D geometry.Only a fraction of those possible paths are valid for the given receiver position as there are only three valid first-order ISs and six second-order ones for the selected receiver.Valid first-order ISs are denoted with an * with their corresponding reflection paths.ISs with their reflection points outside the reflection polygon are denoted by an Â.An IS with a valid reflection point, but with an obstructed path, is denoted as ٗ, and the obstructed reflection path is denoted with a dashed line.FIG. 4. A shoebox-shaped room (solid line) with the original source, ᭺, first order ISs, *, second order ISs, , and third order ISs, ᭛, in their respective image rooms (dashed lines).
In the traditional image-source method, the imagesource tree is constructed in a breadth-first manner.This means that the process iterates over all surfaces for a given reflection order before entering the next higher-order reflection.One can also proceed in depth-first order, in which each branch of the image-source tree is constructed up to the highest order reflection before proceeding to the next IS, i.e., the next branch of the tree.Lee and Lee suggested this approach because it avoids the significant memory requirements of the original method. 106In addition, they present the reflection as a coordinate system transformation that can be efficiently implemented as a 4 Â 4 matrix with homogeneous coordinates.They also suggested to prune the image-source tree if the IS is not visible.However, such an approach is not correct because it does not consider possible higher-order reflections.For example, while a second-order IS may be invisible to a listener, some of the third-order reflections of this IS may again be visible; thus, halting the search for ISs at the first invisible IS is likely to lead to an error.
The main computational problem with the image-source technique is the exponential growth of the number of ISs as seen in the breadth of the IS tree in Fig. 6.The image-source method very efficiently computes early reflections.This method can even be used in virtual reality systems to provide information for real-time auralization. 99,107However, at higher reflection orders, the technique quickly becomes intractable.
A thorough survey of the image-source theory was presented by Mechel. 41In addition, the paper describes several improvements over the traditional algorithm.A special emphasis is placed on corners and how they can be accurately modeled.

B. Ray tracing techniques
Whereas the image-source method provides an exact geometrical solution that consists of all of the specular reflection paths, the premises of ray tracing methods are different.Instead of looking for all the paths deterministically, ray tracing methods are stochastic and perform a Monte Carlo sampling of possible reflection paths.
The very early history of computerized room acoustic ray tracing dates back to 1958 when Allred and Newhouse published their study on the use of ray tracing to detect mean free path lengths in a closed rectangular space. 16However, the use of ray tracing in obtaining time-energy responses in 3D spaces was described by Krokstad et al. in their seminal paper in 1968. 3In parallel, Schroeder presented similar ideas, but it appears that his technique was more limited and operated only in 2D spaces. 8Soon after, several papers that introduced similar ideas were published, but numerous years passed before the technique became popular in practice.This occurred only after computers became more commonplace and powerful; more detailed geometries could thus be addressed and longer responses computed.

Basic implementation
In ray tracing techniques, the main principle is to cast rays from a sound source, reflect them, and register valid paths as illustrated in Fig. 7(a).The rays can be emitted from the source either according to a pre-defined distribution or in random directions via Monte Carlo. 3 However, use of  random directions with too few rays introduces fluctuations in the results and thus the use of proper distributions should be preferred.If a source has a known directivity pattern, the ray distribution can be weighted according to this directivity.Subsequently, the propagation of each ray is traced such that whenever a ray hits a wall, a reflection occurs, and the ray continues in a new direction.Simultaneously, intersections of the rays and detectors are computed to track which reflection paths reach the receiver from the sound source.It is worth noting that the rays are not terminated by hitting a receiver but instead they just pass through the receiver along a straight line.
The termination of the rays and the attenuation of the energy are inter-related and can be addressed using two different methods. 108In the first method, each ray carries information about its energy, and whenever the ray is reflected, the energy is attenuated according to the material properties of the reflector.This can be performed in frequency bands such that each band has its own energy information, and correspondingly, the material absorption is defined in the same bands.However, there cannot be any frequency-dependent scattering that would cause the ray paths to diverge as a function of frequency.The ray is eventually terminated when its energy has decayed below a given constant for all of the frequency bands, or a predefined maximum traveling distance is reached.The other method of implementing ray termination is to apply "Russian Roulette."In this method, in each reflection, a ray has a certain probability of becoming annihilated based on the absorption coefficient of the material.However, this approach is not suitable for the simultaneous tracing of different frequency bands, and is not as efficient as the first one.
One of the problems of the basic ray tracing technique is caused by volumetric detectors.In an ideal case, the detectors would be point-like, but because the probability of a ray to intersect a point is infinitely small, the detectors need to have a volume.However, volumetric detectors cause problems as they might register paths that are not valid for the given listener position, 109,110 especially if the receiver intersects the room geometry.In addition, some of the reflection paths might become registered several times and this needs to be detected if the result is to be used for IS construction.It is also possible to have receivers with variable size such that the size of the receiver increases linearly as a function of the propagation distance. 109However, this technique can be seen as a duality of beam tracing, and such techniques are discussed later in Sec.IV C. In practical room acoustic design, it is often computationally advantageous to have multiple simultaneous receivers in one simulation such that there can be, for example, one receiver at each seat in a model of a concert hall, as proposed by Stephenson. 111tephenson has also conducted an extensive comparison of the image-source method and ray tracing technique and demonstrated that the ray tracing technique exhibits a superior performance in almost all cases. 112he most straightforward way to obtain a time-energy response is to cumulatively register hits in the response each time a ray passes through a receiver, and it is the technique originally applied by Krokstad et al. 3 This approach can be directly used without any weighting for spherical omnidirectional receivers whereas the hits should be weighted by the incidence angle with areal receivers.However, more accurate results can be obtained by taking into account the length of the path a ray travels inside the receiver, and this enables the use of receivers of any shape not just spherical ones.This was shown by Stephenson in a study on sound particle simulation; a term that can be considered to be synonymous to ray-tracing on a general level. 112

Ray tracing and diffuse reflections
One fundamental advantage of the ray tracing approach over the IS technique is the ability to consider arbitrary reflection properties, including diffuse reflections-not simply ideally specular reflections.The basic idea of diffuse reflections in ray tracing was previously discussed in the early papers of Krokstad et al. 3 and Schroeder, 8 but Kuttruff was the first to describe an implementation of ideally diffuse reflections. 58Even more general is the case in which a component of the reflected energy is treated specularly and the other component is treated diffusely; the ratio of the components is determined by the scattering coefficient, as discussed in Sec.III B 3. Realizing diffuse reflections in a ray tracer can be achieved in two conceptually different manners.First, the rays can be split at each reflection into two components: a primary specular ray and numerous diffusely reflected rays equally distributed into the half-space pointing away from the surface as illustrated in Fig. 7(b).Although this process provides a relatively good approximation, it is computationally challenging because the number of rays increases exponentially.A computationally lighter, but less accurate, technique is to maintain a single ray, and in each reflection, the ray will travel in a random direction if the scattering coefficient of the reflecting surface is larger than a random number; otherwise, the reflection is specular, as described by Kuttruff and Straßen 113 and illustrated in Fig. 7(c).Another variant of this single-ray approach is the vector mixing model, in which the new direction of the sound ray is a linear combination of the specular direction vector and of a random Lambertian diffuse direction vector that is weighted by the scattering coefficient. 111,114One way to reduce the computational load is to use a technique known as diffuse rain. 115,116ere the visibility of all the diffuse reflections to the receivers is checked and each visible path is recorded to the receivers taking into account the angle of reflection and the solid angle covered by the receiver as illustrated in Fig. 7(d).In computer graphics, the corresponding technique is known as use of "shadow rays." 117By this means it is possible to reach the same result as without it but with much fewer primary rays.
Another principle for scattering was presented by Mehta and Mulholland in 1976. 118Their approach is based on computing the distance of the ray-surface intersection point to the closest edge.If this distance is shorter than two times the wavelength of the sound wave, then the ray will be scattered to approximate edge diffraction into the visible zone.Their experiments show that this simple modification improves the accuracy of the reverberation time estimates.Christensen and Rindel used scattering caused by edges in another manner.They adjusted the specular reflection strength to obtain an asymptotic high-pass filter effect and assigned the reduced energy to scattering. 119

Accuracy of ray tracing
The stability of a ray tracing solution increases as the number of rays is increased such that the results of the image-source method and ray tracing method converge in a scenario in which only the specular reflections are considered.This means that ray tracing always leads to an approximation of the exact result, where some of the reflection paths may be missed because a finite number of rays must be used.However, in practice, satisfactory results can be easily obtained with ray tracing, but it is hard to know how many rays are needed to obtain reliable results.The challenge here is that the number of required rays depends on the geometry and on the materials such that it is not possible to know the number of rays that are needed before performing a simulation.Subsequently, one can compute, for example, the standard deviation between runs to determine the credibility of the results, as suggested by Kulowski. 120Vorl€ ander has studied the relation of the number of rays and obtained accuracy and showed that the number of required rays is proportional to the volume and absorption area of the space. 108In any case, varying the number of rays provides a good foundation to test the validity of an implementation.If performed properly, the results should converge to a single solution when the number of rays is increased.Similarly, the results should be immune to the underlying tessellation of the geometry.It is typical to split surfaces into simple convex polygons, such as triangles, and how this is performed should not affect the results.

C. Beam tracing techniques
In general, beam tracing means tracing volumetric objects in a geometric model to determine possible specular reflection paths.However, there are two different branches of beam tracing techniques, one being closer to ray tracing and the other more related to the image-source method.Although the name beam tracing describes both techniques well, they are fundamentally different in terms of how the beams are formed.For this reason, they are discussed separately in the following.

Beam tracing as an extension of ray tracing
Basic ray tracing techniques sample the space of possible paths using discrete rays.Clearly, certain paths will be missed if the number of rays is too small with respect to the size of the receiver, or the size of the smallest surfaces.The main principle of beam tracing is to expand these rays to be volumetric objects that can be detected by point-like receivers instead of having straight lines representing rays that need to be detected by volumetric receivers.This approach leads to the advantage that the propagation distances are found exactly, and the pressure-based Eq. ( 3) can be used for the sound propagation.In this type of beam-tracing the space around the source is divided into beams that cover roughly equal solid angles.Figure 8(a) shows a 2D case in which the two highlighted beams are defined by their bounding rays, and their first-order reflections are obtained by specularly reflecting those rays.If a beam hits two or more surfaces as in Fig. 8(b), it needs to be split and reflected separately against all such surfaces.
The first paper to introduce the idea of volumetric rays dates back to 1973 (Ref.121) when Haviland and Thanedar computed room responses in a rectangular room using conical beams instead of simple rays.They suggest the idea of using pyramidal beams, which is currently the favored geometry in beam tracing.However, they did not use the term beam but instead favored the term ray tube.The pyramidal beams have a cross-section of a triangle and can be effectively defined using three rays.One advantage in this structure is that the surface of a sphere can be triangulated relatively accurately and easily, and thus, the entire space around the sound source can be covered by pyramidal beams.
3][124] Their beam tracing system was designed as a general solution to the hidden surface problem in computer graphics.In the Godot system, their technique was used for computing reflection paths and for the detection of diffracting edges.With the computational resources available at the time, they were able to divide the space into 384 beams that were individually traced.In the case of an intersection of a beam with the edge of a surface, they split the beam into two beams as outlined in Fig. 8(b).These new beams do not need to have triangular cross-section anymore but can instead have arbitrary polygonal cross-sections, even such that the beams might end up with holes. 124In addition, they report to generate diffracted beams to be handled with GTD although a detailed description on this is missing. 123lmost parallel to the work by Walsh et al., van Maercke and Martin developed another room acoustic beam tracing tool called Epidaure (originally known as EPIKUL). 125,126The main difference in their technique is that the beams were of conical shape and had a circular cross-section.Such beams can be represented by just one ray standing for the cone axis and one scalar value corresponding to the top angle of the cone.However, such an implementation introduces certain challenges because tessellating a sphere with circles is more complicated and results in overlapping beams.van Maercke and Martin handled this by weighting the contributions such that paths parallel to the cone axis got the full weight whereas in the other extreme the paths outside of the beam had zero weight.However, the main difference between cone-shaped beams and beams having some polygonal cross-section is in splitting of the beams.If a cone-shaped beam hits several surfaces at once, it is not possible to split the beam into new cone-shaped beams with the same coverage whereas with polygonal beams this can be generally performed.
A very similar approach for cone tracing was presented by Vorl€ ander in 1989. 109The technique is based on ray tracing with spherical receivers whose radii increase linearly as a function of propagation distance.To overcome the problem of possible false detections, Vorl€ ander suggested the use of the image-source approach such that the actual paths between the point-like source and receiver are reconstructed and validated.Heinz introduced a hybrid method that further extends this approach using ray tracing to add the diffuse component to the response as well. 127In this technique, the same ray tracing engine can be used to obtain both the specular and diffuse reflections, the difference being that the low-order specular reflections are searched at much higher temporal resolution than are the diffuse reflections.
A good insight into beam tracing variants can be achieved by observing their geometrical accuracy because they tend to be fundamentally different.The following class names are borrowed from computer graphics, where they are often used in relation to visibility queries.The use of triangular beams starts as an exact method because they cover the entire sphere around the source exactly once.Instead, the conical beams can be said to be conservative because they are guaranteed to find all of the reflection paths; however, some paths may be detected multiple times due to certain beams overlapping.The situation becomes more complicated when a beam extending over several surfaces becomes reflected.The applied operations may render it approximate when some of the paths might become registered several times, whereas others might become neglected completely; in the worst-case scenario, false detections may occur.Finally, certain algorithms can be considered to be aggressive because they are guaranteed not to suffer from multiple or false detections, but they might completely miss some paths.
Most of the early beam tracers were not able to properly address situations in which a single beam hits several reflecting surfaces because the beams were not split according to the surface boundaries.This style renders all of the above-mentioned approaches, except for the Godot system, as immediately approximate regardless of the shape of the beam.In general, this problem leads to an overestimation of the detected energy at the receivers, and several modeling software packages, including Epidaurus, 126 ODEON, 24 and Ramsete, 128 included correction factors to compensate for this shortcoming.After the work by Walsh et al., the need to split the beams was noticed and suggested again in 1996 by Stephenson 82 and Dalenb€ ack, 129 but only in 2000 were two very similar geometrical solutions to this problem demonstrated.Both the adaptive beam tracing by Drumm and Lam 130 and adaptive pyramid tracing by Campo et al. 131 suggest splitting the beams into smaller beams upon reflections.Both of these techniques retain the exactness of the initial beam tracing.In addition, the adaptive pyramid tracing can be used in various techniques to reduce computation time and memory requirements by merging nearby narrow beams, but this comes at the cost of exactness.However, the proposed splitting operations make these techniques computationally inefficient; instead, it seems that an even higher performance with exact algorithms can be achieved by a two-phase approach that avoids splitting of the frustumshaped beams. 132In the first phase, a conservative beam tree is constructed using potentially visible sets, 132 and in the second phase, the actual ISs are verified as in the original image-source technique.
In addition to the exact beam tracing described above, various computationally lighter approximate techniques can be used, but they do not guarantee that all of the specular reflection paths will be found.These techniques have numerous applications, especially in computer games and other real-time virtual environments.The frustum tracing technique is such a variation, wherein frustum-shaped beams are loosely cut against reflecting polygons such that they retain their original rectangular cross-section. 36,133Moreover, the adaptive frustum tracing is able to create frustums for diffracting edges, thus further enhancing the realism of the acoustic simulation. 36The diffraction computation is based on the UTD principle.
Similar to other geometrical methods that address only specular reflections, it is possible to compute the actual impulse responses instead of echograms.Jeong et al. compared the accuracy of these two methods with their so-called phased beam-tracing technique.They concluded that the early response at mid-frequencies is more accurate if the phase changes due to propagation are included instead of using pure energies. 134This corresponds to the phase restoration technique described in Sec.III A. In addition, they showed that the use of angle-dependent reflection coefficients only marginally improves the accuracy over the use of angle-independent reflection coefficients.Yousefzadeh and Hodgson compared predictions of energy-and pressurebased coefficients as well, especially from the viewpoint of local and extended reactions. 135Their main conclusion was that a pressure-based model with extended reaction is needed if the space under study includes multi-layered surfaces.
While most beam tracing algorithms can only consider specular reflections, Dalenb€ ack's approach 129 is able to address diffuse reflections as well.The technique uses successive cone tracing passes such that in each higher order pass, the diffuse energy of the previous passes is propagated.This approach provides more freedom to the combination of different reflection types such that sound energy that has been diffusely reflected can still become specularly reflected at a later stage.

Beam tracing as an optimization of the image-source method
The goal in this class of beam tracing techniques is to optimize the performance of the image-source technique by limiting the growth in the number of ISs.This is achieved by culling ISs that are unable to provide any valid reflection paths.The basic image-source method suggests numerous such sources, but it is possible to obtain a remarkable performance gain by pruning the image-source tree as early as possible.The main principle of this type of beam tracing is to minimize the number of beams used.The main difference from the previously described beam tracing technique is that the number of beams increases at each reflection, whereas in the ray tracing style, the number of beams remains constant over time as long as there is no beam splitting.
Whereas the space was divided into equally sized beams that were traced in the previous group of beam tracing techniques, in this family of beam tracers, the beams are formed by the reflecting geometry.The first-order ISs are generated as in the original IS technique, but a beam is subsequently formed for each IS by the edges of the surfaces that were used to create the IS.This implicates that the achieved beam tree is of minimal size whereas in the ray tracing style the beam tree typically starts with plenty of abundant beams thus reducing the efficiency of that approach.After computing the first-order reflections, each beam-not just the IS-is reflected only against surfaces that are at least partially inside the beam as opposed to reflecting it against all surfaces, as in the original image-source technique.The surfaces outside of the beam are not visible to this IS, and thus, the surfaces cannot form valid reflection paths.Therefore, the surfaces can be neglected to save computation time and memory, as outlined in Fig. 9.The resulting IS tree, as seen in Fig. 10, is much narrower that the one by the IS technique (see Fig. 6).This is seen already with the second-order reflections although the savings would be much larger at higher reflection orders.Similarly, if the geometry would consist of smaller surfaces the benefits would be more notable.
The basic principles of this type of beam tracing technique for room acoustics were presented in 1996 by Stephenson. 82The proposed Quantized Pyramidal Beam Tracing addresses also diffraction modeling by introducing new beams at diffracting edges of the room geometry.The key principle in this technique is quantization and reunification of beams.In this process, diffracted and reflected beams are spatially grouped and unified in given time slots.By this means the number of beams to be traced can be reduced.This comes at the cost of geometric accuracy, but as it is applied only at higher reflection orders, this should have no perceptual effect in the result.This reunification process is essential from a computational performance viewpoint especially in such cases in which diffraction or diffusion modeling is taken into account as those will cause continuous beam splittings, whereas with only specular reflections beams will get more and more narrow at higher reflection orders and beam split events are rare.
Funkhouser et al. were the first researchers to present a functional implementation of an image-source-based beam tracer in 1998. 136,137The beam tracer was inspired by the earlier Godot beam tracer described above but more importantly by a computer graphics beam tracer aimed at visual rendering. 138The goal of Funkhouser's work was to construct an auralization system that would be able to consider both specular reflections and edge diffraction.Diffraction, according to the UTD model, was demonstrated using this beam tracing approach by Tsingos et al. 139 Note that the improved image-source technique by Mechel 41 includes beam-like constructs, but in this algorithm, the constructs remain as the actual ISs that are mirror reflected instead of the actual beams.
The goal of both of these beam tracers is to be exact, whereby they can provide all possible specular reflection paths similar to the IS method.This means that the beams have to be accurately cut whenever hit by an occluding FIG. 9. Beams, drawn with dotted lines, are defined by an IS and the edges of the reflecting polygon.Each beam is reflected only against those surfaces that are at least partially inside the beam, whereas the plain image-source method would have created many more ISs.FIG. 10.The IS tree for the geometry in Fig. 5 obtained by beam-tracing.The first-order reflections are the same as with the IS technique but the number of second-order image-sources has decreased from 59 to 39.
surface.This easily makes the beam irregular in shape, and properly cutting them and efficiently computing the intersections becomes tedious.Therefore, variations in which the beam cutting is more relaxed have been developed.This simplifies the consideration of beams but requires another computation pass to verify that the obtained reflection paths can be realized.This operation is similar to the last step of the computation in the above image-source method.Using this approach, it is possible to achieve a high computational performance.However, using the approach combined with utilization of the temporal coherence of the scene, it is possible to construct a real-time auralization system in which a user can freely move and the soundscape is interactively updated. 140ne alternative approach to efficient beam-tree computation was presented by Antonacci et al. 141,142 The principle of their technique is to present the beam tracing from the viewpoint of visibility and to reparametrize rays and lines in a dual space, which is called the ray space, such that a beam in 2D space is represented as a line segment in the ray space.This approach reduces the dimensionality of the problem such that efficient beam-tree construction is enabled for both a moving source and a moving receiver, and the essential visibility queries are reduced to one-dimensional lookups.However, this technique remains limited to only 2D and 2.5D (extruded 2D) geometries, and full 3D implementation has not been presented so far.

D. Diffraction modeling in ray tracing and image-source methods
The extension to include diffraction modeling in GA-based methods has various degrees of importance for different geometries.For outdoor sound propagation, noise barriers and city geometries can certainly benefit from including diffraction modeling.In rooms, the primary objects for which diffraction effects can have a major impact include stage ceiling reflectors, balcony edges, orchestra pits, the presence of pillars, and openings between subvolumes.
Diffraction modeling has lately been incorporated into most popular commercial GA implementations, such as CATT-Acoustic and Odeon, using various approaches. 143ractically all software packages for room acoustical modeling use ray tracing and, in addition, many use the IS method for the efficient handling of low-order specular reflections.
The combination of the image-source method and edge diffraction modeling is typically based on a separate detection of diffracting edges.When all such edges are known, they can be considered as edge sources.The secondary source approach by Svensson et al. was applied to room acoustical geometries, 70,[144][145][146] whereas the computationally more efficient but less accurate UTD formulation has been employed for complex geometries. 139,147Subsequent developments presented the efficient implementation of secondorder diffraction to facilitate complex geometries. 148ombining ray tracing and diffraction is a different challenge.As presented in Sec.III D, approaches to considering this combination in 3D have been suggested. 84,85,149In these techniques, it is necessary to use transparent "portals" at important edges, detect the distance of the ray-portal intersection point and diffracting edges, and modify the ray paths according to that distance.

V. SURFACE-BASED GEOMETRICAL ACOUSTIC MODELING TECHNIQUES
The techniques presented in this section utilize surfaces as an intermediate storage of acoustic energy.This means that techniques of this family can include multiple passes.In the first pass, the sound energy is propagated from a source to the surfaces, and in subsequent passes, the energy is further propagated among surfaces.Finally, the energy is propagated to the receiver.As a comparison, the path-based techniques typically gather the response from the source to the receiver in a single pass one path at a time.This surfacebased approach has certain benefits over the path-based approach.One major benefit is that much of the computation can be performed independently from the receiver position because in these multi-pass techniques, only the last pass involves sound propagation to the receiver.This is especially useful in interactive applications in which all of the previous passes can be precomputed before the actual interactive component is performed, in which the sound energy can be gathered in real time from surfaces to the receiver for auralization purposes.This also benefits room acoustic designs in which the responses are computed from a source to several receiver positions; following the initial passes, it is computationally cheap to obtain the response at any listener position.The downside of this approach is that the geometrical accuracy of an exact reflection path is degraded.Collecting the energy at surfaces always results in an approximation, both in time and space, when compared to actual ray-like paths.
The relation between the methods in this category and the room acoustic rendering Eq. ( 6) is clear because in both of them, the surfaces of the space play a key role.

A. Radiosity
The radiosity method is in stark contrast to the imagesource method presented in Sec.IV A because it is surfacebased and assumes ideally diffuse reflections, whereas the image-source method is path-based and models only ideally specular reflections.The radiosity technique is the most straightforward application of Eq. ( 6) whereby the BRDF qðx 0 ; H; XÞ of each surface has a constant value that is independent of incoming or outgoing angles.
The radiosity technique was first applied in the field of thermodynamics in the 1950s, although the main radiosity equations can be found in an optics paper by Yamauti in 1926. 150Since the 1980s, the radiosity technique has been a popular technique in computer graphics for computing global illumination. 151,152In such applications, the main principle is to search for an energy equilibrium in which the outgoing energy from an element equals the incoming energy multiplied by the absorption of the element.In principle, this can be achieved by representing the rendering equation in matrix form and by inverting the matrix that is formed by the reflection kernels.
The radiosity techniques in optics and computer graphics provide a solution that represents the steady-state, time-invariant solution and that has only limited, although important, use in acoustics.The basic equation in the acoustic radiosity method was presented in 1971 by Kuttruff, 58 as mentioned in Sec.III C.This equation and its descendants were used to determine the reverberation time of a space. 19,59,153,154In these techniques, there is no need to compute the actual response, but the integral equation can be used to determine the decay constant in a more straightforward manner.However, the general goal of a simulation is to obtain impulse or energy responses in more realistic rooms.In the following, we focus on techniques that are able to provide such time-dependent solutions instead of plain scalar values.These are based on the iterative propagation of sound energy that starts from the source and is propagated until most of the sound energy has been absorbed by the surfaces, as schematically illustrated in a 2D setup in Fig. 11.Obtaining such responses requires use of some iterative solver instead of matrix inversion.
The first study to discuss the principles of radiosity for the determination of room acoustic responses was the doctoral thesis of Moore in 1984. 155This work provides a thorough description of the then state-of-the-art in geometrical room acoustic techniques and introduces the radiosity technique.The work is adapted from the heat radiation literature (see, e.g., the work by Sparrow et al. 156 ), wherein the technique is called radiative exchange.The technique utilizes form factors as presented in Sec.III C and also introduces the concept of an exchange factor, which is a cumulative factor for higher-order reflection paths.This hints at the technique not being limited to memory-less diffuse reflections; the technique can address specular reflections of a diffuse reflection.This ability is achieved using the principle of ISs.Instead of propagating the sound to a given surface, it is also propagated to its ISs corresponding to specular reflections.The thesis by Moore provides a very detailed description of the implementation by first defining a new Pascal-like programming language called Sabine, in which the modeling tool is coded.As such, the technique is robust but memory intensive, and it has not gained much popularity in practice.In addition, one needs to pay special attention to not recording the same reflection path multiple times.This mistake was made in the original work by Sparrow, as noted by Moore. 155The same principle was again used in acoustics by Korany et al. 157 However, their formulation suffers from the same defect of multiple detection.
The first implementation of the acoustic radiosity method to gain more publicity is the implementation by Lewers in 1993, which was introduced as part of a hybrid method to compute late reverberation. 65Another radiosity method was briefly explained by Tsingos and Gascuel in 1997. 62The novelty in their approach is the hierarchical nature of the computation, wherein surfaces are divided into finer detail patches only when necessary.
Le Bot and Bocquillet demonstrated, both numerically and theoretically, that the Monte Carlo-based ray-tracing (see Sec. IV B) and collocation-based radiosity methods are equivalent and give the same results when computing reverberation times in an enclosure with ideally diffusely reflecting boundaries. 158 major step in acoustic radiosity was taken in 2004 when Nosal et al. presented, in detail, how the technique can be used to predict the sound field in non-rectangular rooms with ideally diffusely reflecting boundaries. 60The paper also includes validation against analytical solutions in a spherical enclosure.The algorithm computes the radiation density, which is also known as the radiance in later papers, and follows the same iterative propagation principle, as illustrated in Fig. 11.The presented technique was limited to convex FIG.11.In the radiosity technique, the sound energy is (a) first shot from the sound source to all surfaces and then (b) iteratively propagated between surfaces such that eventually, (c) the final response to the receiver is gathered from the surfaces.The actual sound propagation can be performed as pre-processing and is valid for all listener positions.
geometries wherein no visibility computations were needed.However, in practice, computation of the mutual visibility of surface patches is essential to ensure that energy is propagated only to such surfaces that have a line-of-sight to the patch emitting the sound energy.
One elegant feature of the radiosity technique is that sound energy radiating from a surface is angle independent, as implied by the basic assumption, thus resulting in reduced memory requirements.This makes the technique exceptionally suitable for real-time applications where the radiosity simulation can be conducted as a pre-processing step, and the sound energy at the listener position can be gathered at interactive rates (i.e., several times per second), even though the listener is moving.Similarly, the simulation can be performed in a reciprocal manner, whereby the sound source moves and the listener remains stationary.However, if both the source and listener move, a complete recomputation is needed.

B. Acoustic Radiance Transfer
The basic assumption for the radiosity of all surfaces being ideally diffusing is quite limiting and the use of more advanced BRDFs is needed to obtain more accurate simulations.One approach to extending the radiosity technique is to enable the use of arbitrary reflection properties, as proposed by Siltanen et al. in the acoustic radiance transfer (ART) method. 56Such properties make the radiosity technique more practical; however, they simultaneously increase the computation and memory requirements based on the granularity of the spatial discretization.The basic principles of ART remain the same as in the radiosity, and the simulation result is the sound radiances at surfaces; however, the results are now a function of the outgoing direction.This is achieved by spatially discretizing the responses at surface patches such that each discretized direction can have its own response.Thus the reflections are not anymore memory-less, but instead, emitted energy in a given direction depends on the incident angle.If the results are stored as a function of time, they can also be readily used for auralization purposes.For this purpose, performing the simulation in the frequency domain is most efficient because convolution with an anechoic signal is more efficient in this domain compared to the time domain. 64ne application area in which this type of efficient auralization is important is computer games.However, in a typical game scenario, both the sound source and listener can move.The first improvement was presented by Antani et al. 159 Their main principle was to reduce the complexity of the scene such that the scene can be described using a limited number of points, and for each pair of such points, an acoustic transfer function can be precomputed.This means that the sound transport operator is decoupled from the geometry and from the sound source.However, this requires heavy compression of the data and thus reduces the obtained accuracy.The sound from a source is projected onto those points in real time, and the transfer functions are applied to perform the room acoustic simulation.Eventually, the final gathering is performed as in the basic ART and radiosity techniques to provide impulse responses for auralization.

C. Advanced ray-tracing techniques
The basic ray tracing techniques discussed in Sec.IV B attempted to search for specular and diffuse reflection paths between the source and the receiver, whereas the techniques described here utilize ray tracing as a tool to compute the sound energies at the surfaces.A rendering technique called photon mapping is used in computer graphics. 160This technique is based on bidirectional ray tracing, in which light energy is propagated from light sources to surfaces and stored on a photon map; subsequently, an additional ray tracing pass is performed from the viewer to gather the global illumination data from these photon maps.A similar two-phase approach has also been proposed in room acoustic modeling.In the phonon mapping by Bertram et al., 161 sound energy on surfaces is registered on a phonon map, from which the energy is collected to the receiver position in another pass.They used the gathered energy information to form filters for auralization.They mentioned the use of arbitrary reflection properties, but their implementation is only for specular reflections.However, the phonon map does not store any directional information; therefore, the last reflection is always diffuse.Another variation on the same principle is the sonel mapping by Kapralos et al. 77,162 Their technique incorporates both specular and diffuse reflections and is also able to consider frequencydependent edge diffraction such that the results are physically more accurate than are those using the proposed phonon mapping.
The sound particle radiosity method is another variant in this category. 163This method focuses on diffraction modeling and carefully considers computation time and memory usage based on reunification process similarly as in the sound particle simulation of the Quantized Pyramidal Beam Tracing. 82espite this, the huge memory requirement is a bottleneck with this technique, which is still under development.

D. Hybrid methods
The techniques currently applied in practice are typically hybrids of path-and surface-based techniques, whereby the early parts of responses are computed path-based and the later parts are computed surface-based.
The basic beam tracing described in Sec.IV C can only address specular reflections.However, in any realistic scenario, some of the reflected energy will be scattered in a diffuse manner.In 1993, Lewers presented a hybrid room acoustic modeling technique in which pyramidal beam tracing was used to compute the specular reflections, and for each reflection, the diffuse energy was separately recorded and processed using the radiosity principle in another pass. 65his limits the reflections such that following a diffuse reflection, all consecutive reflections will also be diffuse.
The main trend in the development of hybrid methods has been to utilize image-source or beam tracing techniques to detect the early specular reflection paths and various radiosity-based techniques for the late reverberation as there is no need for the accurate detection of the specular reflections to obtain the later component of the impulse response (see Sec. II A); instead, a computationally more efficient but less accurate technique can be used.For example, Lewers used radiosity, 65 and Dalenb€ ack used cone tracing, 129 as previously discussed.An early hybrid technique that combined the image-source method and basic radiosity principles in rectangular rooms was presented in 1983 by Baines. 164A similar hybrid that combines ray-tracing of the early reflections with a surface-based energy transmission algorithm has been presented also by Tenenbaum et al. 165,166 In that model, special emphasis is put on how to achieve as uniform as possible distribution of rays over a sphere.
A recent hybrid method in this family was presented by Koutsouris et al. in 2013. 167The main advantage of this method is that it includes a data-driven, well-defined technique for switching from a specular mode to a diffuse mode, which is separately defined for each reflection path.This allows the technique to capture various phenomena, such as flutter echo, that are easily missed if the transition is performed prematurely.
There are also such hybrid models that make frequency division such that the lowest frequencies are simulated with wave-based techniques whereas the higher frequencies are simulated with GA techniques, see e.g., Refs.168 and 169.

Approximate hybrid techniques for real-time use and auralization
Real-time hybrid techniques often sacrifice substantial physical accuracy for performance.In addition, such systems typically model complex spaces, which may include numerous rooms and corridors.This produces quite different requirements for simulations compared to in the acoustic design of single spaces such as concert halls or auditoria.For such environments, it is advantageous to split the geometric model into separate rooms that are connected to each other by portals as suggested by Schr€ oder et al. 170,171 They conduct separate path tracing in each room and the resulting impulse responses are fed into a filtering network whose topology is derived from the geometric model.
Stavrakis et al. introduced the similar concept of reverberation graphs in 2008 (Ref.172) for this purpose.In that technique, path tracing is used to compute transport operators from sources to portals, between portals, and from portals to listeners.Altogether these operators form a reverberation graph in which the sound can be efficiently propagated from a source to a receiver in any room, thus enabling realistic reverberation in auralization.
One of the most complete physically based auralization system is the RAVEN system developed at RWTH Aachen.The system is well documented in the doctoral thesis by Schr€ oder. 146The system is a hybrid approach of both ISs and ray tracing and also considers edge diffraction.It contains numerous performance optimizations to enable its use in real time.One such topic is the use of different level-ofdetail (LOD) models such that geometrically more detailed models can be used close to the listener whereas some other parts can be simulated with more coarse geometries. 173ntani et al. continued their above-mentioned work 159 by building a hybrid method that uses path tracing to detect early specular and diffuse paths as well as edge diffraction in real time. 174This is combined with pre-computed ART to obtain full-length impulse responses in real time to be used in interactive applications.The technique is able to consider moving sources and listeners and even moving occluders.The latest work in this direction provides further increased efficiency and is able to consider higher-order diffractions at interactive rates. 175Their technique is based on ray tracing and utilizes the temporal coherence of rays such that if there is a change in the scene, most of the rays will not be affected and thus will not require recomputation.

VI. SUMMARY AND FUTURE
A. Summary of GA techniques Different geometrical-acoustic-based modeling techniques have significantly varying capabilities.The strings in Table I show the types of reflection paths that each technique can consider, wherein B is the most general case and covers everything from ideally specular S to ideally diffuse D reflections, while E denotes the possibility of considering edge diffraction.The syntax of the strings is the same as that in regular expressions, in which some parts are alternative and some sections may repeat.For example, the string S{DjS}*D would mean that the technique can address reflection paths that start with a specular reflection and after that there can be both diffuse and specular reflections but the last reflection must be a diffuse one.Most of the beam tracing techniques are denoted by the shape of the beam, whereby the cross-section of a beam is a circle in cone tracing; the cross-section is a triangle in pyramid tracing, and the beams have a convex quadrilateral cross-section in frustum tracing.
Table I shows that no technique can consider the most general expression {BjE}*, in which there could be arbitrary reflection and edge diffraction events in any order.This would require specular reflections being accurately computed such that edge diffraction could also be accurately modeled.From this perspective, the ray tracing technique appears to be the most attractive option because it is essentially able to consider arbitrary reflection properties and because it can be combined with different diffraction models.Similarly, capable techniques are the ART and surface-based ray tracing techniques.The only difference here is that these techniques are unable to accurately include edge diffraction and must resort to energybased approximations.However, this is important only for specular reflections in which phase is essential, whereas for diffuse reflections, the phase is irrelevant.In addition, the construction of fully specular reflection paths from the ray tracing results requires a separate validation procedure for the ISs.In practice, the specular paths are relevant only at very low reflection orders, while most of the response should be obtained with energy-based techniques.
The geometrical acoustic modeling techniques exhibit different computational performances.In Table II, we have selected a subset of techniques to represent the most advanced current modeling techniques.The table lists their basic performance characteristics in terms of efficiency.This ranking is slightly arbitrary because the required modeling accuracy and efficiency are closely related and application dependent.Regardless, Table II is more informative for applications whereby a plausible result is sufficient, whereas from the room acoustic design perspective, it might not be significantly relevant.Here, the main conclusion is that if the scene is fully dynamic, then path-based ray tracing techniques are the most appropriate choice.In such conditions, pre-computation is not performed; only temporal coherence can reduce the computation load, as suggested by Schissler et al. 175 Instead, if both the geometry and either the sound sources or the listener are static, then the most efficient method is to utilize a surface-based technique such as radiance transfer.In such techniques, it is possible to pre-compute the sound propagation to all surfaces, and the resulting response to the listener can be efficiently obtained from these surfaces.

B. Accuracy of GA modeling techniques
The accuracy requirements of a room acoustic simulation heavily depend on the application.These requirements are the highest in room acoustic design, and thus, software tools with high reliability are required.Simultaneously, there is a substantial amount of responsibility on the user of the software whereby the main principles and limitations of the GA model must be acknowledged.It is important to note that a GA model is only an approximation of the complex behavior of sound in a room.Within such limitations, the accuracy of the modeling techniques still depends on several factors.These factors can basically be divided into two groups.
First, there are weaknesses in the models and in their implementation.Of these weaknesses, the basic assumptions of GA might be the most limiting factor.At high frequencies, the GA model is relatively accurate, but at lower frequencies, those assumptions fail, and the lack of wave phenomena starts to play a major role and renders the results inaccurate.Edge diffraction modeling can only remedy this problem to a certain degree.Another cause of inaccuracy is the GA modeling techniques being limited to the local reaction reflection model, whereas in reality, non-local, extended reactions are common.The implementation of different techniques is not always trivial, which is an issue that can easily result in inaccuracies in the modeling results.One such factor is the treatment of the details of the geometry.The basic principle in GA is that all geometry details should be an order of magnitude larger than the longest wavelength of sound in the simulation.Whether the finer details should be smoothed out and how they should be addressed remain open questions.
The second set of factors that affect the practical accuracy is related to the quality of the input data.The available material data are typically substantially coarser than what is usable by modern modeling techniques; for example, real-world materials do not have angle-dependent absorption data available.Similarly, the accuracy and the LOD in the geometry data often impose challenges in modeling.If the geometry is too coarse, then the errors can be significant.Similarly, if the geometry is too detailed, certain modeling techniques may experience problems, as previously [S ¼ Specular, D ¼ Diffuse, B ¼ Arbitrary BRDF, E ¼ Edge diffraction, j ¼ or, x Ã repeats, x þ occurs 0 or 1 times, and x q occurs a low number (q) of times.]The abbreviations for the diffraction models are as follows: HF ¼ Huygens-Fresnel, BTM ¼ Biot-Tolstoy-Medwin, UTD ¼ Universal Theory of Diffraction, and HU ¼ Heisenberg's uncertainty principle.
mentioned.In practice, the most useful simulation results in acoustic design are typically obtained with simplified geometric models 181,182 of the space, whereas in computer games, the use of the same geometry for visual and aural rendering is preferred.The availability of widely accepted benchmark cases can facilitate the development of more accurate modeling techniques.Such well-documented cases with geometric models, material information, and publicly available measurement results would facilitate validation for those who develop these techniques.

C. Future trends
Currently, there appears to be a gap between the communities that investigate accurate room acoustic modeling techniques for acoustic design and those who focus on the entertainment industry.We believe that these areas would benefit from collaboration.In particular, the highperformance techniques developed for computer games could benefit the more serious modeling area as well.Although high performance does not correlate with increased accuracy, various improvements that do not sacrifice accuracy can be implemented.For example, the use of graphics processing units (GPUs) for ray tracing 179,183 enables much higher ray counts and higher reflection orders compared to what is customary on current central processing units (CPUs).Similarly, advanced beam tracing techniques for exactly finding all purely specular reflection paths are worth studying for use in room acoustic modeling tools.Finally, we would like to highlight the need for collaboration with researchers in the field of visual rendering in computer graphics.Although there are fundamental differences between these fields, there are numerous similar techniques that are studied in both areas, and the synergy between the two fields is clear.

FIG. 1 .
FIG. 1.(a) An impulse response of a room represents the propagation of sound pressure from a source to a receiver, while (b) the time-energy response represents the propagation of sound energy.

FIG. 5 .
FIG.5.Example of a room with an original source, ᭺, and a receiver, þ.Valid first-order ISs are denoted with an * with their corresponding reflection paths.ISs with their reflection points outside the reflection polygon are denoted by an Â.An IS with a valid reflection point, but with an obstructed path, is denoted as ٗ, and the obstructed reflection path is denoted with a dashed line.

FIG. 7 .
FIG. 7. (a) Specular ray tracing will eventually find all the same reflection paths as the IS technique although the volumetric receivers may cause some false detections.(b) Diffuse reflections can be modeled either by splitting a ray into several new rays in a reflection, or (c) by redirecting rays to nonspecular directions.(d) In the diffuse rain technique a shadow-ray is cast from each diffuse reflection point to the receiver to speed-up convergence of the simulation.

FIG. 6 .
FIG.6.The IS tree for the geometry in Fig.5.The boxes represent possible ISs generated by the surfaces labeled in the box.For example, the first-order IS produced by surface 2 will create 8 second-order ISs.The growth-rate in the number of possible ISs is very high.There are 9 possible first-order ISs, but the number of second-order ISs is already 59.The circles denote valid ISs for the selected receiver position.

FIG. 8 .
FIG. 8. (a)In ray-tracing style of beam-tracing, the space around a source is split into equispaced beams that are reflected at surfaces.Two such beams are illustrated here.(b) If a beam hits several surfaces when being reflected it should be split into smaller beams such that there will be one beam per each separate reflecting surface.
and 2005 (Ref.28) represent milestones in this field.A thorough review of different room acoustic modeling techniques was presented in 2002 by Svensson and Kristiansen.

TABLE I .
Possible reflection paths for each modeling technique.