Rigid Walled Room Revisited - Part 2
Stefano Tronci
In the Rigid Walled Room Revisited - Part 1 episode we seen how to setup an eigen system problem with Elmer by making use of the Wave Equation solver. In this episode we will review the simulation results and check the agreement with the analytical solutions in the Acoustic Modes of a Rectangular Room episode.
Project Files
All the files used for this project are available at the repositories below:
The Solver Output Data and Data Processing
In the Rigid Walled Room Revisited - Part 1 episode we setup two solvers: the Wave Equation solver, which outputs a vtu file, and a Save Scalars solver, which outputs a dat file. The vtu file will contain the eigenmodes (eigenfunctions) while the dat file will contain the eigenvalues.
In order to ease comparison with the exact eigenvalues and eigenmodes these files are read with a Julia script. The Julia script will take care of computing the exact quantities from the equations presented in the Acoustic Modes of a Rectangular Room episode and produce comparison metrics with the Elmer output.
Resonance Frequency Error Metric
To compute the resonance frequency error the FEM eigenvalues are first converted in resonance frequencies by using this formula:
$$\lambda_{n} = \omega_{n}^{2} = \left(2 \pi \nu_{n} \right)^{2}$$
$$\Rightarrow \nu_{n}=\frac{\sqrt[+]{\lambda_{n}}}{2 \pi}$$
with $n$ an index going from $1$ to $10$ used to select the modes, sorted by increasing resonance frequency, $\lambda_{n}$ being the eigenvalue, $\omega_{n}$ being the angular frequency of oscillation of the mode and $\nu_{n}$ its associated frequency. If we call the exact resonance frequency $f_{n}$ we can compute the error between $\nu_{n}$ and $f_{n}$ in cents as follows:
$$E_{n}=1200 \log_{2} \left( \frac{\nu_{n}}{f_{n}} \right)$$
To be noted that the very first eigenvalue is supposed to be $0$ for every system and it is associated with a constant eigenmode. However, FEM algorithms might fail to produce $0$ as a value and instead produce a small negative number. Small negative eiganvalues will just be rounded to $0$.
Eigenmode Error Metric
Whilst it is possible to define many single valued metrics of error it is more interesting to look at the error over the whole domain, so to understand whether there are patterns. The study of the patterns can indicate ways to improve the study. For this reason the error was evaluated by simply subtracting the exact eigenmode from the FEM one, thus obtaining the following error fields:
$$\epsilon_{n} \left(x, y, z\right)=s_{n} p_{n}\left(x, y, z\right) - S_{n}\left(x, y, z\right)$$
Where $p_{n}$ is the FEM eigenmode, $S_{n}$ is the exact one and $s_{n}$ is a coefficient that can be either $1$ or $-1$. The use of $s_{n}$ is made necessary by the fact that the FEM solution can be affected by polarity inversion, that is, the positive and negative parts of the wave can be exchanged. By introducing the $s_{n}$ correction factor the error field $\epsilon_{n}$ will be only sensitive to the error in the field rather than its overall sign, which can be easily corrected for.
To be noted that the field $S_{n}$ is computed at exactly the same coordinates of the nodes at which $p_{n}$ is computed. Then, $\epsilon_{n}$ is stored over a mesh structured as the original FEM mesh. This is the best way to produce an accurate point-to-point comparison.
Implementation
The processing of the data is implemented in the file validate.jl available in the repository linked above. The code makes use of the meshio package to read vtu files.
The vtu files output by Elmer contain $10$ scalar fields, one for each eigenmode. This makes processing with ParaView somehow cumbersome. The code hence separates the fields in the original file into many vtu files that can be imported into ParaView as a time series. In addition, the error fields and corrected FEM fields $s_{n} p_{n}\left(x, y, z\right)$ are also exported in a similar fashion. You can read the code to see how this is done in details.
Results
Resonance Frequencies
The table below shows a comparison of the FEM resonance frequencies with the exact ones (rounded to two decimal places):
| Mode Number | FEM Resonance Frequency $\text{Hz}$ | Exact Resonance Frequency $\text{Hz}$ |
|---|---|---|
| $1$ | $0$ | $0$ |
| $2$ | $34.31$ | $34.30$ |
| $3$ | $42.90$ | $42.88$ |
| $4$ | $54.96$ | $54.91$ |
| $5$ | $57.23$ | $57.17$ |
| $6$ | $66.77$ | $66.67$ |
| $7$ | $68.70$ | $68.60$ |
| $8$ | $71.59$ | $71.46$ |
| $9$ | $79.44$ | $79.26$ |
| $10$ | $81.07$ | $80.90$ |
At first glance we can already see how the FEM resonance frequencies are all slightly higher than the exact ones. The error is best visualised by plotting the error metric $E_{n}$ defined above:
Click Here for Larger Version
It is possible to see that, overall, the error increases with the mode number, and hence with frequency. This has to be expected, given that the mesh was fixed. In fact, as we already noticed in the Dealing with Convergence Issues and Mesh Order and Accuracy episodes the smaller the mesh size is in relation to the size of the spatial features of the field the higher the accuracy will be. Higher frequency fields are more complex and have shorter wavelength, so FEM solutions are expected to be less accurate with respect lower frequency fields if they are computed with the same mesh size.
Whether this level of accuracy on the resonance frequency is acceptable depends on the application. If musical acoustics is the scope of the study, for example, we can compare the errors to the just noticeable pitch differences for humans. One good reference for this is Psychoacoustics - Facts and Models by Fastl, Hugo, Zwicker and Eberhard. At page 186 of the third edition we can read:
… our hearing system is more sensitive to frequency changes if the task is to recognize differences rather than to recognize modulations. The pause between the two sounds to be compared does not reduce the sensitivity, on the contrary, it increases it! Displacing the data given in Fig. 7.8 down by a factor of three produces a reasonable approximation to the results at the two asymptotes: i.e. at frequencies below 500 Hz, we are able to differentiate between two tone bursts with a frequency difference of only about 1 Hz; above 500 Hz, this value increases in proportion to frequency and is approximately 0.002f.
The authors report that the human hearing system is very good at comparing pitches. So, we could imagine to be playing burst of tones at the exact resonance frequency and bursts of tones at the FEM one. If we were listening to them we would be in the conditions for our ear to be maximally sensitive to the frequency difference between them. All of our frequencies are lower than $500$ $\text{Hz}$. Hence, none of the differences in frequency would be audible to humans, being all the FEM frequencies within $1$ $\text{Hz}$ from the exact ones.
In other words, this FEM analysis is more accurate than human hearing, making this level of accuracy acceptable in the context of musical acoustics. Of course, this is by using just noticeable frequency differences as a guide. Referring to other psychoacoustic phenomena might yield to a different outcome. In short: whether your application is musical acoustics or not, you should do through research in order to understand what levels of errors you can accept.
Eigenmodes
First thing first, let’s take a look at the eigenmodes themselves. They are presented in the following video.
The video cycles through all the modes while showing the FEM solution on the left and the exact solution on the right. The fields are presented as surfaces of constant value (remember that the eigenmodes are normalised to $1$ both in value and unit). It is possible to see that the shape of the eigenmodes is essentially always correct. The biggest differences between the two being the sign: certain FEM solutions appear to have inverted sign. This is a common occurrence in FEM and normally it doesn’t come with significant drawbacks. As we mentioned above, we can use the $s_{n}$ coefficients to correct for the sign of the FEM solutions. The $s_{n}$ values are reported in the table below.
| Mode Number | $s_{n}$ |
|---|---|
| $1$ | $-1$ |
| $2$ | $-1$ |
| $3$ | $+1$ |
| $4$ | $+1$ |
| $5$ | $-1$ |
| $6$ | $-1$ |
| $7$ | $-1$ |
| $8$ | $-1$ |
| $9$ | $-1$ |
| $10$ | $+1$ |
With the coefficients above we can construct the error fields $\epsilon_{n}$, which are shown in the video below.
The video shows the error field on the left and the corrected FEM eigenmode $s_{n} p_{n}$ on the right for reference. It is possible to see that the worst case error (the maximum of the scale of the error shown in the video) is $1.4 \cdot 10^{-2}$, that is, the error is limited to $1.4\text{%}$. Again, only the application will indicate whether this error is acceptable or not. However, we can note that an error of $1.4\text{%}$ is associated to a $0.1$ $\text{dB}$ change in magnitude, which is about the highest accuracy that state of the art acoustic measurement equipment can provide. Hence, our FEM results appear comparable to those that a measurement of the exact field with state of the art instrumentation would provide. In this context we can then say ourselves to be satisfied by the results.
However, few things jump to the eye:
- For all of the eigenmodes the errors are higher where the eigenmode magnitude is higher.
- The highest errors are observed for mode number $5$.
- Errors seems to spread more significantly in the quieter zones for the highest frequency mode.
Point 1 is perhaps related to the internals of the FEM algorithm used to solve for the eigenmodes.
Point 2 is interesting as the mode number $5$ varies only along the $z$ axis. Considering that the mesh size is essentially the same along all axes, it means that mode $5$ varies along the axis with the least number of elements. This might be the root cause of the drop of accuracy for mode $5$.
Point 3 is actually a bit tricky, as the quite zones are narrower for high frequency modes, so it is not fully clear whether the errors are higher in the quieter zones or the quieter zones are just smaller. In both cases the end result is that the higher the frequency the more evenly the errors distribute over the domain, which is another good reason to refine mesh quality if high accuracy at high frequency is needed.
Conclusion
In this episode we reviewed the results from the Eigen System study we setup in the Rigid Walled Room Revisited - Part 1 episode. As a reminder, the maximum mesh size was $170$ $\text{mm}$ and the mesh was a first order mesh.
We found good resonance frequency agreement, within $4$ $\text{cent}$ (and within $1$ $\text{Hz}$) and with a trend of increasing resonance frequency error with mode number, as expected.
The eigenmodes were found to be satisfyingly accurate, with errors within $0.1$ $\text{dB}$ distributed mainly where the eigenmode magnitude is higher. Mode number $5$, which varies along the axis with least amount of elements, was the mode affected by the greatest errors.
The results above suggest that Wave Equation eigensolver does a very good job at finding modes of acoustic systems. We can then proceed to use this solver for more complex studies.
License Information
This work is licensed under a Creative Commons Attribution 4.0 International License.