# Elmer vs FEniCS - Part 2

In the Elmer vs FEniCS - Part 1 episode we solved the Pulsating Sphere problem over a fairly coarse mesh. We solved the problem both with `Elmer`

and `acoupy_helmholtz`

, our `FEniCS`

based Helmholtz solver. We found that `Elmer`

and `FEniCS`

produce the same solution. This solution compares well to the exact solution, but has some appreciable magnitude error. In this episode we will study how this error changes with mesh size, and whether `Elmer`

and `FEniCS`

keep on providing the same solution.

# Project Files

All the files used for this project are available at the repositories below:

# The Problem

We solve exactly the same problem outlined in Elmer vs FEniCS - Part 1 | The Problem, but with meshes of different resolution.

We use a finer mesh with $10$ $\text{mm}$ maximum size and a fine mesh with $5$ $\text{mm}$ maximum size. We will refer to these as `10 mm`

and `5 mm`

throughout.

Again, we solve the problem on all meshes with both `Elmer`

and `acoupy_helmholtz`

. Refer to `simulate_and_compare.py`

from `acoupy_helmholtz`

for reference on how to implement the simulation. Remember that `acoupy_helmholtz`

(as well as `acoupy_meshutil`

) has a documentation with tutorials and examples accessible from their `README.md`

file.

# Results

## Solution

As done in Elmer vs FEniCS - Part 1 | Results, we take each mesh node and evaluate the solution at the node. We compute the distance of the node from the source centre and plot magnitude and phase of the solutions as functions of this distance. We do this because the Pulsating Sphere field depends only on distance. The solutions are shown below.

`10 mm`

Click Here for Larger Version

`5 mm`

Click Here for Larger Version

### Discussion

We can see that both solutions agree fairly well with the expected solution. The agreement in phase is good in both cases. However, the solutions underestimate the field slightly, the `10 mm`

having the higher error as expected. We will do a more quantitative error analysis in the next section.

## Error

As done in the Mesh Order and Accuracy episode, for each mesh and solver we compute the quotient between the solution and the exact field. By taking the decibel value of the magnitude of this quotient we will have the magnitude error in $\text{dB}$. The angle of this quotient is the phase error. The results are shown below.

`10 mm`

Click Here for Larger Version

`5 mm`

Click Here for Larger Version

### Discussion

We can see that the errors for the phase are fairly low for both meshes, being within $\pm 0.017$ $\text{rad}$ ($\pm 1$ $\text{deg}$) for all solutions. However, it is also evident that the `5 mm`

solution has way lower errors.

Magnitude errors on the other hand are somewhat significant for the `10 mm`

mesh. Again, the errors breach the $\pm 0.1$ $\text{dB}$ boundary. In fact, in absolute terms, they reach $0.3$ $\text{dB}$. We observed a similar total spread of error in Elmer vs FEniCS - Part 1 | Results for a $28$ $\text{mm}$ mesh. However, the `10 mm`

mesh error is actually better behaved. This because the error pattern is relatively flat around a lower mean value. The magnitude errors for the `5 mm `

solution are instead very low, well within the $\pm 0.1$ $\text{dB}$ threshold.

As observed in in Elmer vs FEniCS - Part 1 | Results the error patterns are almost identical. This shows that the solvers are providing the same solution and the inaccuracies are rather due to the mesh itself.

Since there is no substantial difference between the `Elmer`

error and the `FEniCS`

error we can pick any of the solutions to further characterise the error as a function of mesh size. One good way to analyse error is by computing its statistics. As evident from the plots above the error is not really random. In this case the errors appear larger the closer we are to the source. However, error statistics are still useful to quickly gauge the overall trends of the error. In the plots below we report the Kernel Density Estimates of the error as computed with `KDEpy`

. Kernel Density Estimation (KDE) is a technique to estimate the probability distribution of a random process. It is mathematically equivalent to a histogram with uniform bins passed through a smoothing filter. Wikipedia has a good overall introduction. We will not discuss the details of KDE in this episode. Rather, we use a KDE as an instrument to visualise spread and mean. We chosen Gaussian kernels and narrow bandwidth to show the trends without excessive bias.

Click Here for Larger Version

The plots show how, by reducing mesh size, the magnitude error reduces in spread as well as gaining an average value closer to $0$. On the other hand, the phase error appears to be centred more closely around $0$, and mainly reducing in spread.

## A Note on Execution Time

`Elmer`

solutions were very fast to compute due to having used the `BiCGStabl`

iterative method. With `FEniCS`

we used the `MUMPS`

direct method, which is slower. `FEniCS`

allows to select multiple solution methods as well. We will investigate those in the next episodes. Many of the methods are also available with `Elmer`

. Hence, you can refer to the Elmer Model and Solver Parameters episode for an overall introduction.

# Conclusion

`Elmer`

and `FEniCS`

output very well matched solutions for meshes of different density when the element parameters are kept consistent between the two solvers. This shows that, in terms of accuracy, `acoupy_helmholtz`

essentially matches `Elmer`

’s Helmholtz solver. `acoupy_helmholtz`

has the same features as the `Elmer`

Helmholtz solver, plus two additional material parameters to describe damping. This shows how versatile writing our own software is: we are now able to expand the `acoupy_helmholtz`

solver further with source terms and adiabatic absorbers (for example). All by not affecting accuracy. However, `acoupy_helmholtz`

solutions were computed with a direct `MUMPS`

solver, which is slower with respect iterative solvers. We will explore iterative `FEniCS`

solvers in the next episodes.

# License Information

This work is licensed under a Creative Commons Attribution 4.0 International License.