15,792,217 members
Articles / Programming Languages / C#

# Propagation of Sound in a Room

Rate me:
10 Mar 2023CPOL9 min read 14.9K   270   23   10
Calculate sound transmission between a point source and a receiver point within an enclosed space
This article will go through the steps and limitations of calculating the sound transmission between two points in a room by using modal analysis. This will use Fourier transformation in order to do the calculations based on the wave equation and fully reflective walls. In addition, it will also give two different models for adding dampening to the room, one uses reverberation time and the other will use the absorption coefficient.

## Introduction

Sound transmission between points within a room is one of the most important aspects in practical acoustical problems. In many aspects, it is the cornerstone of what we do, since so many things are dependent on the sound pressure transmission from a source to its listeners. The relation calculated by the geometrical damping and the damping at the boundaries (or wall, ceiling and floor) given a sound source with a reference sound pressure that spreads the sound equally in all directions. The diagram shown is the damping of a flat response of 0 at the source and therefore the negative pressure at the receiver.

The general procedure for calculating the sound pressure level of a receiver have two main approaches: geometrical acoustics which is only valid for high frequency and Fourier acoustics generally only used for the low frequency part due to the number of nodes needed in the calculation increase drastically for increasing frequency.

## Wave equation

We start off with the classical time and space depended wave equation:

$\frac{d^2 \, f(t,x)}{dx^2} = \frac{1}{c^2} \frac{d^2 \, f(t,x)}{dt^2} \tag{1} \label{eq:Wave}$

This equation describes all the details of the sound transmission between all possible sources and receivers, but it is very time consuming to do the calculations. There is a simple way assuming that you do not need all details of the transmission. We perform Fourier (or Laplace) transformation on the differential equation itself and get Helmholtz equation instead:

$\nabla^2 f(x) + k^2 \cdot f(x) = 0 \tag{2} \label{eq:Helmholtz}$

This is what acousticians call a steady state solution to the wave equation, since you integrate over all time to eliminate the time dependency in the equation. For a mathematician, this can lead to another problem: Gibbs phenomena. The Fourier transformation of the wave equation actually eliminates some disjoint initial conditions and solutions that it cannot correct. Actually, the history of Fourier transform begins in 1750 with Euler when he considered the initial condition which is the Fourier transformation of the curve. An English translation of Euler's original work is given here: PDF version of an English translation and there are more translations of more of Euler's work that are also given here. In the article in question, he wrote the initial condition as:

$\alpha \sin(\pi s/a) + \beta \sin(2 \pi s/a) + \gamma \sin(3 \pi s/a) \tag{3} \label{eq:EulerIC}$

Which is exactly a summation of the Fourier coefficients for as also indicated by Daniel Bernoulli for suspended strings, which he published in 1753. Euler had some reservations as he didn't think that the Fourier transformation was general enough, and he is actually correct in his analysis. The Weierstrauss theorem will actually not fix this. In fact, the only reason that Fourier transformation is valid, is an article by Carleson that finally set the question on firm ground. The formulated proof is known as Carleson's theorem states that Fourier transform converges almost everywhere. The issue is that there must be a finite number of disjoint steps in the function, then it will converge. And it will converge almost everywhere except for the disjoint step. This can be illustrated by a step function in blue and the Fourier transformation approximation in red shown below:

What you do to reproduce the error is to give a step function as the initial condition, take the Fourier transformation of it, and at last the inverse Fourier transformation.

The reason for the discrepancy in the figure is simply that the Fourier transform is infinitely smooth, meaning that infinitely many derivatives exist in every point. While the step functions it tries to emulate is missing the derivative at the step point. The reality is that the Fourier transform can accurately replicate any function in $L^2(f(x))$, which is a fancy way of saying that on average it will be correct, but that does not mean that the approximation from the Fourier transform are equal in value at every point.

## Solving Helmholtz equation

We now assume that all the boundaries of the room are infinitely hard. This will mean that the particle velocity at the wall must be zero, which in return will mean that the sound pressure will be the maximum at the same point. Now we are searching for a solution that has its maximum value at the boundaries, so this would be achieved by harmonics of the cos function which is also a solution to the homogenous wave equation.

$p(x,y,z) = C \cdot \cos\left(\frac{M \pi x}{V_x}\right) \cdot \cos\left(\frac{L \pi y}{V_y}\right) \cdot \cos\left(\frac{K \pi z}{V_z}\right) \tag{4} \label{eq:HelmholtzBC}$

Assuming that we solve Helmholtz equation in each direction using a power series method. We generally want to add a source to the differential equation:

$\nabla^2 f(x) + k^2 \cdot f(x) = - i \cdot \omega \cdot \rho_0 q(r_0) \tag{5} \label{eq:HelmholtzInhomogenous}$

This is a inhomogeneous Helmholtz equation that can be solved using a power series method. The source distribution function. This can be expressed as individual frequency amplitudes:

$q(r_0) \cdot e^{i \omega t} = \sum_{N=0}^{\infty} Q_N \Psi_N(r) e^{i \omega t} \tag{6} \label{eq:InsertedEquation}$

One also needs the property of the form function to follow certain laws. Mainly that the orthogonality property will give:

$V = \iiint\limits_{V} \Psi^2_N (r) \, dV$

This will require some shape property and can be illustrated by each of the energy integrals from 1,2 and 3 nodes contributing to the sound level.

$w_{axial} = \frac{1}{L_{x}} \int_0^{L_x}{\cos^2{ \left( x \frac{L \cdot \pi}{L_x} \right) \, dx}} = \frac{1}{2}$
$w_{tangential} = \frac{1}{L_{x} L_y} \int_0^{L_x} \int_0^{L_y} \cos^2{\left( x \frac{L \cdot \pi}{L_x} \right)} \cos^2{\left( y \frac{M \cdot \pi}{L_x} \right)} \, dx dy = \frac{1}{4}$
$w_{oblique} = \frac{1}{L_x L_y L_z} \int_0^{L_x} \int_0^{L_y} \int_0^{L_z} \cos^2{\left( x \frac{L \cdot \pi}{L_x} \right)} \cos^2{\left( y \frac{M \cdot \pi}{L_y} \right)} \cos^2{\left( z \frac{K \cdot \pi}{L_z} \right)} \, dx dy dz = \frac{1}{8}$

This means that the axial mode will have double the energy than an tangential node, and that a tangential node will have twice the energy of the oblique node. However the same integrals will appear with the orthogonal requirement, but there we demand that each of these three integrals yields the same result. Therefore, each time either L, M or K is zero, we will have to compensate with 2 for each mode that is zero and 1 for each mode that are not zero.

Applying the power series to every element in the differential equation gives you the following resulting equation:

$\sum_{N=0}^\infty A_N \nabla \Psi_N(r) + \sum_{N=0}^\infty A_N k^2 \Psi_N(r) = -i \omega \rho_0 \sum_{N=0}^\infty Q_N \Psi_N(r) \tag{7} \label{eq:TotalHH}$

Rearranging the terms:

$\sum_{N=0}^\infty A_N \Psi_N(r) \left[ k^2 - k_N^2 \right] = -i \omega \rho_0 \sum_{N=0}^\infty Q_N \Psi_N(r) \tag{8} \label{eq:TotalClean}$

Setting up the final solution with regards to the pressure at the receiving point:

$p(x,y,z) = \rho_0 c_0^2 Q \sum_{L=0}^\infty \sum_{M = 0}^\infty \sum_{K = 0}^\infty \frac{\omega \cdot \Psi_{L,M,K}(x,y,z) \cdot \Psi_{L,M,K}(x_0,y_0,z_0)}{V_{L,M,K} \cdot \left( \omega^2 - \omega^2_{L,M,K} \right)} \tag{9} \label{eq:Final}$

## Include damping

In room acoustic the damping of the room from soft surfaces gives a substantial effect on the sound level in a room.

$m \cdot \frac{d^2 x}{dt^2} + k (1 + i \cdot \eta) \cdot x = 0 \tag{10} \label{eq:PendulumWave}$

Assuming that the solution can be written with the complex $s = \sigma + i \tau$ on the form $e^{i s \cdot t}$. Inserting this expression into the previous differential equation will give the following

$s = \sqrt{\frac{k}{m} (1 + i \eta)} \approx \omega_0 \left( 1 + i \frac{\eta}{2} \right) \tag{11}\label{eq:ComplexSolutionWave}$

We now assume that the energy of any sound waves will be negligible after it has decayed by 60 dB. This is a historical picked value as speech typically have a level at 60 dBA so this was assumed to be the level which you could no longer hear the sound $\eqref{eq:ComplexSolutionWave}$ . This was suggested in the 1900's before the advent of electrical equipment that actually made it possible to measure the actual sound level, so this logic was assumed to be correct. Now it is just a convenient way of expressing the damping value of sound and vibrations in a confined space, and the resulting damping value is dB per time unit. In practice signal to noise ration determines the usable range, and this is seldom more than 20 - 30 dB in room acoustics. Its used in a general way for acousticians to describe damping in both rooms and in plates, wall etc.

$\eta = \frac{\log \left(10^{6} \right) }{T_{60} \cdot \omega} = \frac{\log\left(10^{6} \right)}{T_{60} \cdot 2 \pi f} = \frac{2.2}{T_{60} \cdot f} \tag{12} \label{eq:Damping}$

All that is left is to gather this information to the reverberation time or the absorption factor of the surfaces in the room. $\eqref{eq:Damping}$

$\omega_{L,M,K} = \omega_{L,M,K} \cdot \sqrt{1+ i \cdot \eta } = \omega_{L,M,K} \cdot \sqrt{1 + i \frac{4.4 \cdot \pi}{\omega_{L,M,K} \cdot T_{60}}} \tag{13} \label{eq:TotalComplete}$

We now insert this damping into the previous equation:

$p(x,y,z) = \rho_0 c_0^2 Q \sum_{L=0}^\infty \sum_{M = 0}^\infty \sum_{K = 0}^\infty \frac{\omega \cdot \Psi_{L,M,K}(x,y,z) \cdot \Psi_{L,M,K}(x_0,y_0,z_0)}{V_{L,M,K} \cdot \left( \omega^2 - \omega^2_{L,M,K} - \frac{i \cdot \omega_{L,M,K} \cdot 4.4 \cdot \pi}{T_{60}} \right)} \tag{14} \label{eq:TotalCompleteReverberation}$

There is another way of estimating the damping term by simply inserting the experimentally found variable of damping $\eta = \frac{\alpha}{8}$ by using experience from measurements and together with some approximation techniques.

$p(x,y,z) = \rho_0 c_0^2 Q \sum_{L=0}^\infty \sum_{M = 0}^\infty \sum_{K = 0}^\infty \frac{\omega \cdot \Psi_{L,M,K}(x,y,z) \cdot \Psi_{L,M,K}(x_0,y_0,z_0)}{V_{L,M,K} \cdot \left( \omega^2 - \omega^2_{L,M,K} - i 2 \omega_{L,M,K} \left( \epsilon_L \frac{\alpha_{x}}{8 L_x} + \epsilon_M \frac{\alpha_{y}}{8 L_y} + \epsilon_K \frac{\alpha_{z}}{8 L_z} \right) \right)} \tag{15} \label{eq:TotalCompleteWithAlpha}$

This concludes this article. Plans for extending the solutions for boundaries with impedances is on its way, unless you can beat me to it.

## History

• 23rd February, 2023: Initial version

<quillbot-extension-portal><quillbot-extension-portal>

Written By
Chief Technology Officer
Norway
This member has not yet provided a Biography. Assume it's interesting and varied, and probably something to do with programming.

 First Prev Next
 Tough Room Member 1488449217-Mar-23 13:27 Member 14884492 17-Mar-23 13:27
 Re: Tough Room Kenneth Haugland18-Mar-23 7:00 Kenneth Haugland 18-Mar-23 7:00
 Re: Tough Room Member 1488449220-Mar-23 12:27 Member 14884492 20-Mar-23 12:27
 Re: Tough Room Kenneth Haugland21-Mar-23 7:58 Kenneth Haugland 21-Mar-23 7:58
 Re: Tough Room Member 1488449222-Mar-23 9:00 Member 14884492 22-Mar-23 9:00
 Too dirty Sergey Alexandrovich Kryukov23-Feb-23 18:11 Sergey Alexandrovich Kryukov 23-Feb-23 18:11
 Re: Too dirty Jeffamn24-Mar-23 4:55 Jeffamn 24-Mar-23 4:55
 Re: Too dirty Sergey Alexandrovich Kryukov24-Mar-23 5:32 Sergey Alexandrovich Kryukov 24-Mar-23 5:32
 strange program Member 1531459123-Feb-23 9:50 Member 15314591 23-Feb-23 9:50
 Re: strange program Kenneth Haugland23-Feb-23 9:56 Kenneth Haugland 23-Feb-23 9:56
 Last Visit: 31-Dec-99 19:00     Last Update: 1-Dec-23 3:31 Refresh 1