The very first episode in which we introduced Elmer was the Elastic Modes of a Metal Bar episode. In that episode we introduced Elmer by solving a linear elasticity eigenproblem, and mentioned how vibration is an integral part of acoustics, being vibrating bodies one of the principal causes of airborne sound radiation. In this new series of episode we will explore vibration further, integrating in it what we learned so far, and we will explore vibro-acoustic coupling with Elmer.

Project Files

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

• Tuning Fork 1.

A Note on Terminology

In this episode, quite inappropriately, the handle of the tuning fork will be referred as prong while its prongs will be referred as tines. This due to the study having been setup with this terminology originally. Sorry for the confusion…

A Note About Series

As you have probably noticed, in this project there are multiple series going on, for example the Home Studio or the Rigid Walled Room series. Most of these series are not concluded. In fact, each topic could be followed up by many further episodes, as I plan to do. Rather to conclude a series before taking up the next, I decided to develop all series in parallel, so to prevent the project to focus too much on the details of each single problem and instead allowing us to explore Elmer (and possibly other solvers) capabilities more freely. This is why a new series is being introduced now.

Modes of a Tuning Fork

In the Elastic Modes of a Metal Bar episode we used Elmer to solve for the normal modes of vibration of a metal bar. Today we solve for the normal modes of a tuning fork. The analysis we will be performing is largely inspired by those done by Ben Qui and Justin Black (archived here), which I invite you to read beforehand. In fact, we will use the same CAD model for the tuning fork as provided by Justin. I should also mention that tuning fork CAD models are available in other places online, for example on GRABCAD. In this episode we will use Justin’s tuning fork, which is a $512$ $\text{Hz}$ tuning fork. Whilst the problem can be fully setup within ElmerGUI we will set it up by writing the sif file ourselves. As mentioned in other episodes in which we did just that, this is actually beneficial, as not all solvers have a GUI module currently, so becoming confident with the sif file helps us maximising the control over the Elmer solver features.

Model Setup

We will setup this model to search for the first $10$ eigenvalues and eigenfunctions of the tuning fork, with few different boundary conditions:

• Free.
• Prong bottom simply supported.
• Prong simply supported.

The Tuning Fork

The tuning fork itself will be, for this episode, our domain of interest. In order to solve for its mode properly through FEM we need to first understand its properties better.

A tuning fork is a Y shaped piece of metal, were two bars are hold parallel and joined at the base. The fork terminates with a prong (or handle) by which it can be hold (see Figure 1 below).

First of all, we are interested in its nominal frequency which, as we already stated, is $512$ $\text{Hz}$.

Then we are interested in its material, which is aluminium. Justin Black provided the Amazon vendor link for the fork which, unfortunately, only mentions “non-magnetic aluminum alloy” as a material description. We will then follow Justin’s heuristic analysis and use the parameters of Aluminum 6061 for our fork, which are reported below:

Last, but not least (in fact, one can argue we are mostly interested in this, given the deep implications it has for meshing) we are interested in the actual fork geometry. The fork is depicted below with few size annotations. We will refer to these right in the next section, as they directly impact our meshing choices.

Geometry Preprocessing and Meshing

As we typically do, once we have a CAD file we can simply export it as a BREP file to import into Salome. If you use Justin’s CAD file be aware that he scaled the fork down by a factor of $1000$, most likely to avoid having to use the Coordinate scaling keyword within Elmer. The CAD in our repository is instead in physical size and units, as shown in Figure 1 above.

Within Salome we will proceed to explode the geometry in a Solid and Face entities, as we typically do (you can refer to the previous episodes if you need guidance). However, when carrying on this operation in Salome 9.6.0 the behaviour might be slightly different with respect previous Salome version. I recommend that you explode the imported BREP entity first in a solid and then in faces, so that the solid and the faces appear directly below the BREP entity in Salome’s Object Browser, as shown below. This will ensure that the correct groups are created when we mesh the BREP entity. Note that I renamed the faces for the prong (handle) in some meaningful way, so to have them easily located for the application of boundary conditions (as we will see later).

In the Elastic Modes of a Metal Bar episode we used a regular grid mesh. It would be possible to do so also here, by using Salome’s Body Fitting algorithm. However, the mesh produced like that will work within Elmer only with a low Threshold, which in turn produces a “blocky” results as shown below.

Clearly this style of meshing is unsatisfactory as it distorts the original shape of the tuning fork too much. A better fit can be achieved by using a higher value for the Threshold parameter for the Body Fitting algorithm. However, this results in the creation of polyhedral elements that neither ElmerGrid or salomeToElmer seem to be very good at dealing with (conversion of the mesh to Elmer format will result in errors). Hence, we will follow a different root. Even though we will not use Body Fitting for this problem you can refer to the Salome study in the repository to see the details of the algorithm settings that produced the mesh above.

As we seen, since our tuning fork has quite a number of round edges Body Fitting is quite not the best algorithm to mesh it unless we want to create hard-to-deal-with polyhedral elements. Hence, we will fallback to our old trustworthy Netgen algorithm. To mesh our geometry we simply select, in Salome’s Mesh module, our top level BREP entity (geometry.brep_1 with reference to Figure 1), then we crate a new mesh with Mesh > Create Mesh, as always. We only need to figure out what mesh sizes we need.

The overall sizes of tuning fork along the various axes are:

What typically happens with eigenmodes is that the first modes tend to develop along the longest dimensions first. The lower the order the lower the number of local peaks. For example, maybe the first mode will involve vibration of the entirety of the fork, with one local peak only. Then its various parts will be able to vibrate with more then one local displacement peak for higher modes. What we want is to have at least $10$ elements between each of these local peaks. But we do not know yet how many peaks we will see, as we still have to solve the study first.

We will have a first guess by considering the biggest size first. If we choose a maximum element size smaller than $1.64$ $\text{mm}$ we will be able to “tile” the longest axis of the fork with more than $100$ elements, which would be effective for up to $10$ local peaks along this axis, which is a very high number of local peaks, most likely not encountered within the first $10$ modes. With reference to the bar sizes reported in Figure 1, this will also produce between $5$ and $6$ elements per upper edge of the bar (these edges being $7.15$ $\text{mm}$ and $9.56$ $\text{mm}$ respectively). This should be a good mesh density to capture the first rotational modes of the bars, which instead we should expect to see in the first $10$ modes. As a result, we will then use a maximum size of $1.5$ $\text{mm}$. The full mesh parameters are reported below together with a picture of the mesh. Note that Salome 9.6.0 will automatically create the necessary mesh groups from geometry: we do not need to follow that step anymore.

Note how the mesh has been kept to first order. This because we will use $p$-elements within Elmer to set the order.

After we get our first results, we will inspect the resulting displacement field and figure out whether our mesh needs additional refinement.

To export the mesh to Elmer format we can right click it from Salome’s Object Browser, select Export and then UNV File. To convert it to Elmer format, we can use the ElmerGrid command as follows:

ElmerGrid 8 2 Netgen.unv

Assuming that the mesh was exported as Netgen.unv. The 8 and 2 arguments simply specify the input and output file formats respectively. For more information, see the ElmerGrid Manual. This will create a folder called Netgen that contains all the needed mesh files for our project.

Elmer Study

We are now ready to setup our study with Elmer. All we need to do is:

• Create a folder.
• Copy the mesh files into the folder.
• Write our sif files.

Create a Project Folder

Simply create a folder and copy all the contents of the folder created by ElmerGrid inside. For example, I named my folder elmerfem. Its contents, for the time being, will then be:

• mesh.boundary
• mesh.elements
• mesh.header
• mesh.names
• mesh.nodes

Writing the sif Files

We now need to write three sif files, one for each of the studies we want to make, each with different boundary conditions:

• Free.
• Prong bottom simply supported.
• Prong simply supported.

In our project folder we can then create the empty text files below, which we will proceed to fill:

• case_free.sif
• case_bottom_prong.sif
• case_whole_prong.sif

Turns out that these sif files are all essentially the same, so we will go through the process to write one and simply point out when there needs to be a difference between them.

Header Section

This section is the same for every sif file. It simply says Elmer to search for the mesh in the current directory and write the results in the current directory:

Header
  CHECK KEYWORDS Warn
  Mesh DB "." "."
  Include Path ""
  Results Directory ""
End

Simulation Section

In this section we define the main simulation parameters. This section is the same for all sif files aside for the values of the Solver Input File and Post File keywords. Note that we set the Coordinate Scaling keyword so that Elmer can correctly interpret the coordinates value of the mesh nodes. Also note that we specify the vtu format for the output. The actual file name will have a timestap appended, for example case_free_t0001.vtu, as always. The sections are shown below.

case_free.sif

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Steady state
  Steady State Max Iterations = 1
  Output Intervals = 1
  Coordinate Scaling = 0.001
  Solver Input File = case_free.sif
  Post File = case_free.vtu
End

case_bottom_prong.sif

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Steady state
  Steady State Max Iterations = 1
  Output Intervals = 1
  Coordinate Scaling = 0.001
  Solver Input File = case_bottom_prong.sif
  Post File = case_bottom_prong.vtu
End

case_whole_prong.sif

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Steady state
  Steady State Max Iterations = 1
  Output Intervals = 1
  Coordinate Scaling = 0.001
  Solver Input File = case_whole_prong.sif
  Post File = case_whole_prong.vtu
End

Constants Section

This section is the same for all sif files. We simply define a number of useful constants in SI units. Most likely these will not be used by the solver, with the possible exception of Gravity. However, it is a good idea to have these always in, just in case one wants to extend the model with additional solvers.

Constants
  Gravity(4) = 0 -1 0 9.82
  Stefan Boltzmann = 5.670374419e-08
  Permittivity of Vacuum = 8.85418781e-12
  Permeability of Vacuum = 1.25663706e-6
  Boltzmann Constant = 1.380649e-23
  Unit Charge = 1.6021766e-19
End

Body Sections

Since we have only one body, we need only one body section. Since there can be more than one body in one simulation, this section needs an ID, which we can specify simply after the Body keyword. We use this section to declare a body in the simulation and declare what its material and governing equation are. This section is the same for all sif files.

Body 1
  Target Bodies(1) = 1
  Name = "Body 1"
  Equation = 1
  Material = 1
End

The value of the Target Bodies array is chosen by referring to the mesh.names file, reported below:

! ----- names for bodies -----
$ Solid_1 = 1
! ----- names for boundaries -----
$ Face_1 = 2
$ Face_2 = 3
$ Face_3 = 4
$ Face_4 = 5
$ Face_5 = 6
$ Face_6 = 7
$ Face_7 = 8
$ Face_8 = 9
$ Face_9 = 10
$ Face_10 = 11
$ Face_11 = 12
$ Face_12 = 13
$ Face_13 = 14
$ Face_14 = 15
$ Handle_Side = 16
$ Handle_Bottom = 17
$ bnry18 = 18

As you can see, there is only one body (our tuning fork), whose ID is 1. Hence, the array Target Bodies has size 1 (hence the (1) just beside it) and it contains only the body with ID 1.

In the next sections we will define the Equation and Material, both with ID 1, which are assigned to this Body.

Solver Sections

We first need a solver to handle the body governing equation. Again, we can have more than one solver, so we need to specify and ID.

We define first the solver that handles the linear elasticity of the body, which is the same for all sif files. The section is reported below:

Solver 1
  Equation = Linear elasticity
  Eigen System Select = Smallest magnitude
  Eigen Analysis = True
  Eigen System Values = 16
  Procedure = "StressSolve" "StressSolver"
  Element = "p:2"
  Variable = -dofs 3 Displacement
  Exec Solver = Always
  Stabilize = True
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-5
  Nonlinear System Convergence Tolerance = 1.0e-7
  Nonlinear System Max Iterations = 1
  Nonlinear System Newton After Iterations = 3
  Nonlinear System Newton After Tolerance = 1.0e-3
  Nonlinear System Relaxation Factor = 1
  Linear System Solver = Direct
  Linear System Direct Method = Umfpack
  Linear System Convergence Tolerance = 10e-10
End

Note that we set the required keywords to enable the eigen analysis (Eigen Analysis = True) and search for the first 16 eigenmodes (Eigen System Select = Smallest magnitude and Eigen System Values = 16). We search for the fist $16$ because, for the free boundary condition, the first six modes are just $0$ $\text{Hz}$ rigid body translations and rotations. The actual interesting modes will then start from the seventh. In order to include all the first $10$ nonzero frequency modes we then compute a total of $16$.

Note also as we will be using $p$-elements of the second order (Element = "p:2"). This will turn our mesh in a second order mesh, increasing accuracy. Finally, we will set the Nonlinear System Max Iterations to 1, being the problem linear, and use a Direct solver method as these work really well for elasticity (in fact I could not get iterative methods to even converge for this problem).

In addition to this solver, it is useful to create another utility solver, a SaveScalars solver. This solver will simply save the computed eigenvalues to a file, so that we can read the file and avoid parsing the Elmer log to get them. Moreover, the eigenvalues are written to the files with much higher precision with respect the log. Since each solver is solving for a different boundary condition, this second utility solver must be different for each sif file. Note that defining this solver is not quite possible yet with ElmerGUI. Having access to all these features is one of the pros of writing a sif file from scratch.

case_free.sif

Solver 2
  Equation = "SaveScalars"
  Procedure = "SaveData" "SaveScalars"
  Filename = eigenvalues_free.dat
  Save EigenValues = True
End

case_bottom_prong.sif

Solver 2
  Equation = "SaveScalars"
  Procedure = "SaveData" "SaveScalars"
  Filename = eigenvalues_bottom_prong.dat
  Save EigenValues = True
End

case_whole_prong.sif

Solver 2
  Equation = "SaveScalars"
  Procedure = "SaveData" "SaveScalars"
  Filename = eigenvalues_whole_prong.dat
  Save EigenValues = True
End

Equation Sections

We can now finally define the governing equation we assigned to our Body 1 previously:

Equation 1
  Name = "Linear Elasticity"
  Active Solvers(1) = 1
End

The ID of this Equation is 1, to match that we assigned in the Body. This equation is handled by one single solver, Solver 1. Hence, the Active Solvers array has only 1 element, the element 1, which is the ID of the solver required for this equation.

This section is the same for all sif files. We do not need any other equation for this study.

Material Sections

We only need one material, with ID 1 to match the one we assigned in our Body 1 section. This is the same for all sif files and it simply declares the properties of the Aluminium 6061:

Material 1
  Name = "Aluminum 6061"
  Poisson ratio = 0.33
  Density = 2712.63
  Youngs modulus = 68.9e9
  Porosity Model = Always saturated
End

Boundary Condition Sections

Since our sif files all differ for boundary conditions, these sections are different for the three solvers.

case_free.sif

This is a model of a free tuning fork, so no boundary conditions are applied: the tuning fork is floating in the middle of empty space.

case_bottom_prong.sif

In this case the fork is simply supported at the bottom of its prong. This means that the bottom of the prong is fixed, it cannot displace. Hence:

Boundary Condition 1
  Target Boundaries(1) = 17 
  Name = "Simply Supported"
  Displacement 3 = 0
  Displacement 2 = 0
  Displacement 1 = 0
End

All the coordinates of displacement are set to 0. Since there is only one face of our fork defining the bottom of the prong, the Target Boundaries array has simply one element, face 17. To identify which one is the correct face, you can use the mesh.names file, which will make it especially easy if you named the faces to something meaningful within Salome.

case_whole_prong.sif

In this case the entirety of the prong of the fork is simply supported, the bottom and the side. Hence:

Boundary Condition 1
  Target Boundaries(2) = 16 17 
  Name = "Simply Supported"
  Displacement 3 = 0
  Displacement 2 = 0
  Displacement 1 = 0
End

All the coordinates of displacement are set to 0. The Target Boundaries array has now 2 elements, faces 16 and 17, which together define the boundary of the prong. Again, to identify the faces to list in the Target Boundaries array you can use the mesh.names file, which will make it especially easy if you named the faces to something meaningful within Salome.

Running the Studies

Running the studies is simple. Just open a terminal in the project directory, and issue the command ElmerSolver followed by the sif file you want to solve for. The commands below solve the various studies:

ElmerSolver case_free.sif
ElmerSolver case_bottom_prong.sif
ElmerSolver case_whole_prong.sif

The various vtu and dat files will be created in the project directory.

Conclusion

In this episode we setup a tuning fork eigenproblem. We experimented with the meshing and settled on a Netgen mesh, which we then raised to second order within Elmer. We decided to write three sif files to simultaneously solve for three different boundary conditions, and leveraged the flexibility of the sif file approach to define an additional helper solver to export raw data from the simulation. In this new linear elasticity simulation we were able to put together a lot of things we learned from previous simulations, mainly good meshing paradigms, the use of $p$-elements, sif file writing. We were able to expand a little upon sif file writing and we also learned how to implement multiple variations of a simulation.

In the next episode we will be reviewing the results.

License Information

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

• Physics
• Elmer
• FEM
• Linear System
• Single Physics
• Eigensystem
• Series: Tuning Fork