## 2 Theoretical Modelling  
Acoustic scattering occurs when a travelling sound wave encounters an obstacle or inhomogeneity in its path, e.g., a solid object or a change of medium density, and thus breaks into secondary spreads out from it in a variety of directions. The magnitudes and directions of such scattered waves can be described by the Helmholtz-Kirchhoff theory of diffraction, which uses [PERSON]'s theorem to determine the scattered pressure at any given point by specifying the wave field onto the scattering obstacle. However, the latter formulation is generally approximated depending of contextual assumptions regarding the distance of the source and observer from the scattering object compared to the monochromatic wavelength. This leads to two major situations, i.e., (i) the near-field where the source and observer are relatively close to the obstacle, under the Rayleigh critical distance (\(R_{0}=S/\lambda\), where \(S\) is the area of the surface and \(\lambda\) the wavelength of the wave), and (ii) the far-field where the source and observer are considered to be far away from the obstacle, in the limit of infinity. On a first instance, we will focus on the frequency-dependent scattering made by the metasurface, where the acoustic field at a point \(\mathbf{r}=\mathbf{r}(x,y,z)\) scattered by a surface centred at \(\mathbf{r}_{0}=\mathbf{r}(x,y,z=0)\) can be approximated by the Rayleigh-Sommerfeld integral as  
\[P_{s}(\mathbf{r})=-i\frac{k}{2\pi}\int_{S_{0}}\frac{p_{0}(\mathbf{r}_{0})R( \mathbf{r}_{0})e^{ik|\mathbf{r}-\mathbf{r}_{0}|}}{|\ \mathbf{r}-\mathbf{r}_{0}\ |}\mathrm{d}S, \tag{1}\]  
where \(p_{0}(\mathbf{r}_{0})\) is the incident pressure field, and \(R(\mathbf{r}_{0})\) is the spatially dependent reflection coefficient of the locally reacting surface \(S_{0}\), with \(k=\omega/c_{0}\) being the wavenumber in air at the angular frequency \(\omega\) and speed of sound \(c_{0}\). Here, the Fourier transform of \(f(x)\) is \(\hat{f}(k)\), which is given by \(\hat{f}(k)=\int f(x)e^{-ikx}\mathrm{d}x\). It transpires from Equation (1) that defining the state of the spatially dependent reflection coefficient, \(R(\mathbf{r}_{0})\), is ultimately important to determine the directions and magnitudes of the scattered sound energy. In the case of metadiffusers, the spatially dependent surface reflection coefficient is obtained through the Transfer Matrix Method [30] (TMM), which relates the acoustic pressures and normal particle velocities at the extremities of a one-dimensional acoustic system; here, a slit loaded with Helmholtz resonators (HRs). The surface reflection coefficient and characteristic impedance of the \(n\)-th slit, \(R_{slit}^{n}\) and \(Z_{slit}^{n}\), respectively, can be interchangeably deduced one from the other by the relation  
\[Z_{slit}^{n}=Z_{0}\frac{1+R_{slit}^{n}}{1-R_{slit}^{n}}, \tag{2}\]  
where \(Z_{0}=\rho_{0}c_{0}\) is the characteristic impedance of air. The diffusion coefficient of a surface rates the uniformity of the aforementioned scattered sound field. Moving to a spherical coordinate system where \(P_{s}(\mathbf{r})=P_{s}(\theta,\phi,r)\) with \(\theta\) and \(\phi\) being the elevational and azimuthal planes respectively and \(r\) being the distance to the origin, the directional diffusion coefficient [33], \(\delta_{\phi}\), produced when a sound diffuser is radiated by a plane wave at the incident angle \(\psi=(\theta^{\prime},\phi^{\prime})\) (primed superscripts  
Figure 1: Conceptual diagram illustrating (**a**) a render of a Quadratic Residue Metadiffuser (QRM) with \(N=5\) slits, (**b**) a numerical equivalent surface impedance (ZEQ) of the metadiffuser identically composed of \(N=5\) slits, and (**c**) an implementation of the ZEQ into a large computational volume, e.g., an orchestra pit.