## WaveSolver and Eigen Analysis

Numerical methods and mathematical models of Elmer
CrocoDuck
Posts: 79
Joined: 12 May 2016, 13:15
Antispam: Yes

### WaveSolver and Eigen Analysis

Hello!

I am trying to solve for the acoustic eigenmodes of a 3D rectangular room by using the WaveSolver. I am using this as a guide.

See below my sif file:

Code: Select all

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

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Steady State Max Iterations = 10
Output Intervals = 1
Timestepping Method = BDF
BDF Order = 1
Coordinate Scaling = 0.001
Solver Input File = case.sif
Post File = case.vtu
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

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

Initial Condition 1
Pressure Velocity = Real 1
End

Equation 1
Name = "Wave Equation"
Active Solvers(1) = 1
End

Solver 1

! Equation
Equation = "Wave Equation"
Procedure = "WaveSolver" "WaveSolver"
Variable = "Pressure"
Variable DOFs = 1

! Solver

Steady State Convergence Tolerance = 1e-09

! Nonlinear
!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
!Linear System Scaling = Logical True
!Linear System Symmetric = Logical True
Linear System Solver = string "Iterative"
Linear System Direct Method = Umfpack
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = integer 10000
Linear System Convergence Tolerance = real 1e-8
BiCGstabl polynomial degree = 4
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = True
Linear System Residual Output = 50
Linear System Precondition Recompute = 1

! Special
Eigen Analysis = True
Eigen System Values = 10
Eigen System Convergence Tolerance = Real 1.0e-6
Eigen System Max Iterations = 10000
Eigen System Select = smallest magnitude
!Eigen System Shift = Real 1000
Eigen System Normalize To Unity = Logical True

End

Material 1
Name = "Air (room temperature)"
Heat expansion Coefficient = 3.43e-3
Heat Capacity = 1005.0
Heat Conductivity = 0.0257
Sound speed = 343.0
Viscosity = 1.983e-5
Relative Permittivity = 1.00059
Viscosity = 1.983e-5
Porosity Model = Always saturated
Turbulent Prandtl Number = 0.713
Relative Permittivity = 1.00059
Density = 1.205
Relative Permeability = 1.00000037
Sound Speed = Real 343.0
Sound Damping = Real 0.0
Sound Reaction Damping = Real 0.0
End

Boundary Condition 1
Target Boundaries(6) = 2 3 4 5 6 7
!Pressure = Real 0
!Source Acceleration = Real 0
Pressure Velocity = Real 0
End

``````
My issue is that I always get an infinite norm:

Code: Select all

``````ELMER SOLVER (v 8.4) STARTED AT: 2020/08/18 16:15:10
ParCommInit:  Initialize #PEs:            1
MAIN:
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Version: 8.4 (Rev: unknown, Compiled: 2020-08-17)
MAIN:  Running one task without MPI parallelization.
MAIN: =============================================================
MAIN:
MAIN:
MAIN: -------------------------------------
Model Input:  Unlisted keyword: [pressure velocity] in section: [initial condition 1]
Model Input:  Unlisted keyword: [pressure velocity] in section: [boundary condition 1]
LoadInputFile: Number of Body Forces: 0
LoadInputFile: Number of Initial Conditions: 1
LoadMesh: Scaling coordinates: 1.000E-03 1.000E-03 1.000E-03
LoadMesh: Elapsed REAL time:     0.0863 (s)
MAIN: -------------------------------------
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: wave equation...done.
OptimizeBandwidth: Half bandwidth without optimization: 15311
OptimizeBandwidth:
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 1064
OptimizeBandwidth: ---------------------------------------------------------
WARNING:: InitCond: We can only set timederivative for transients
ElmerSolver: Number of timesteps to be saved: 10
MAIN:
MAIN: -------------------------------------
MAIN: -------------------------------------
MAIN:
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: wave equation
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    241600
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.01
1   Infinity
ERROR:: IterSolve: Numerical Error: System diverged over maximum tolerance.
Note: The following floating-point exceptions are signalling: IEEE_DIVIDE_BY_ZERO
STOP 1

``````
Or, if I try a direct solver:

Code: Select all

``````ELMER SOLVER (v 8.4) STARTED AT: 2020/08/18 16:20:29
ParCommInit:  Initialize #PEs:            1
MAIN:
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Version: 8.4 (Rev: unknown, Compiled: 2020-08-17)
MAIN:  Running one task without MPI parallelization.
MAIN: =============================================================
MAIN:
MAIN:
MAIN: -------------------------------------
Model Input:  Unlisted keyword: [pressure velocity] in section: [initial condition 1]
Model Input:  Unlisted keyword: [pressure velocity] in section: [boundary condition 1]
LoadInputFile: Number of Body Forces: 0
LoadInputFile: Number of Initial Conditions: 1
LoadMesh: Scaling coordinates: 1.000E-03 1.000E-03 1.000E-03
LoadMesh: Elapsed REAL time:     0.0840 (s)
MAIN: -------------------------------------
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: wave equation...done.
OptimizeBandwidth: Half bandwidth without optimization: 15311
OptimizeBandwidth:
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 1064
OptimizeBandwidth: ---------------------------------------------------------
WARNING:: InitCond: We can only set timederivative for transients
ElmerSolver: Number of timesteps to be saved: 10
MAIN:
MAIN: -------------------------------------
MAIN: -------------------------------------
MAIN:
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: wave equation
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
ERROR:: EigenSolve:  Error with DNAUPD, info =           -9
STOP 1

``````
My feeling is that something is not OK with the boundary conditions. I tried zero and non-zero values for the boundary conditions without any success. I also played with all the Linear System options still without success.

I am using Elmer on Ubuntu 18.04, installed through the PPA.

Any idea?
raback
Posts: 4700
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

### Re: WaveSolver and Eigen Analysis

Hi,

Try with the same keywords as the example. Now you switched at least material parameters to something different. The default is zero if the parameters are not found.

-Peter
CrocoDuck
Posts: 79
Joined: 12 May 2016, 13:15
Antispam: Yes

### Re: WaveSolver and Eigen Analysis

Hi Peter,

Thank you very much for your help. I tried to change the material parameters to this (as in the example):

Code: Select all

``````Material 1
Name = "Ideal"
Wave Speed = Real 1.0
End
``````
But I will get this error:

Code: Select all

``````ERROR:: ListInitElementKeywords: Keyword required but not present: Sound speed
``````
Wave Speed is not a correct keyword according to the manual. The available material parameters appear to be:
• Sound Speed
• Sound Damping
• Sound Reaction Damping
CrocoDuck
Posts: 79
Joined: 12 May 2016, 13:15
Antispam: Yes

### Re: WaveSolver and Eigen Analysis

Also, this line:

Code: Select all

``````Procedure = "WaveEq" "WaveSolver"
``````

Code: Select all

``````Procedure = "WaveSolver" "WaveSolver"
``````
I tried the same sif as in the example, a part the few changes mentioned above:

Code: Select all

``````Check Keywords "Warn"

Mesh DB "." "."
Include Path ""
Results Directory ""
End

Simulation
Coordinate System = "Cartesian"
Steady State Max Iterations = 1
Coordinate Scaling = 0.001
Max Output Level = 5
End

Body 1
Equation = 1
Material = 1
Initial Condition = 1
End

Initial Condition 1
Potential Velocity = Real 1
End

Material 1
Name = "Ideal"
Sound Speed = Real 1.0
End

Equation 1
Active Solvers(1) = 1
End

Solver 1
Equation = "Poisson"

Variable = "Potential"
Variable DOFs = 1
Procedure = "WaveSolver" "WaveSolver"

Steady State Convergence Tolerance = 1e-09

!  Linear System Solver = "Direct"
Linear System Direct Method = Umfpack

Linear System Scaling = Logical True
!  Linear System Symmetric = Logical True
Linear System Solver = string "Iterative"
Linear System Convergence Tolerance = real 1e-8
Linear System Iterative Method = BiCGStabl
BiCGStabl Polynomial Degree = 4
Linear System Residual Output = integer 1
Linear System Max Iterations = integer 2000
Linear System Preconditioning = ILU2
Linear System Residual Output = 50

Eigen System Select = smallest magnitude
!  Eigen System Shift = Real 1.0e3

! For visualization purposes make amplitude one for all
Eigen System Normalize To Unity = Logical True

Eigen Analysis = True
Eigen System Values = 10
Eigen System Convergence Tolerance = Real 1.0e-6
Eigen System Max Iterations = 1
End

Solver 2
!  Exec Solver = never
Equation = "result output"
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = case
Vtu Format = Logical True
Single Precision = Logical True

Eigen Analysis = True
End

Solver 3
Equation = "SaveScalars"
!  Filename = f.dat

Procedure = "SaveData" "SaveScalars"
Save EigenValues = True

! Compare the norm of the 10th eigenvalue
Show Norm Index = 10
End

Boundary Condition 1
Target Boundaries(6) = 2 3 4 5 6 7
Potential = Real 0
!  Potential Velocity = Real 0
End

Solver 3 :: Reference Norm = 1.730641551602E+002

``````
Same results:

Code: Select all

``````ERROR:: IterSolve: Numerical Error: System diverged over maximum tolerance.
Note: The following floating-point exceptions are signalling: IEEE_DIVIDE_BY_ZERO
STOP 1
``````
I see the example was last updated in 2017. Perhaps it is not up to date?
raback
Posts: 4700
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

### Re: WaveSolver and Eigen Analysis

Hi

It turned out that the updated WaveSolver under modules didn't assemble MASS for harmonic/eigen systems. So it could not have worked. I fixed it and removed the older WaveEq.f90 from the test case.

https://github.com/ElmerCSC/elmerfem/co ... 2f0ed93dd1

-Peter
CrocoDuck
Posts: 79
Joined: 12 May 2016, 13:15
Antispam: Yes

### Re: WaveSolver and Eigen Analysis

I see, thanks for looking into this! I will maybe look into building the devel branch to test it out.
raback
Posts: 4700
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

### Re: WaveSolver and Eigen Analysis

Hi,

You can just compile the updated solver:

Code: Select all

``````elmerf90 -o WaveSolver.so WaveSolver.F90
``````
(or .dll in Windows).

-Peter
CrocoDuck
Posts: 79
Joined: 12 May 2016, 13:15
Antispam: Yes

### Re: WaveSolver and Eigen Analysis

raback wrote: 20 Aug 2020, 10:06 Hi,

You can just compile the updated solver:

Code: Select all

``````elmerf90 -o WaveSolver.so WaveSolver.F90
``````
(or .dll in Windows).

-Peter
Ha! Nice one. That works! I finally got a solution.

Code: Select all

``````ELMER SOLVER (v 8.4) STARTED AT: 2020/08/21 11:59:36
ParCommInit:  Initialize #PEs:            1
MAIN:
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Version: 8.4 (Rev: unknown, Compiled: 2020-08-17)
MAIN:  Running one task without MPI parallelization.
MAIN: =============================================================
MAIN:
MAIN:
MAIN: -------------------------------------
Model Input:  Unlisted keyword: [pressure velocity] in section: [initial condition 1]
Model Input:  Unlisted keyword: [pressure velocity] in section: [boundary condition 1]
LoadInputFile: Number of Body Forces: 0
LoadInputFile: Number of Initial Conditions: 1
LoadMesh: Scaling coordinates: 1.000E-03 1.000E-03 1.000E-03
LoadMesh: Elapsed REAL time:     0.0841 (s)
MAIN: -------------------------------------
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: wave equation...done.
OptimizeBandwidth: Half bandwidth without optimization: 15311
OptimizeBandwidth:
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 1064
OptimizeBandwidth: ---------------------------------------------------------
WARNING:: InitCond: We can only set timederivative for transients
ElmerSolver: Number of timesteps to be saved: 10
MAIN:
MAIN: -------------------------------------
MAIN: -------------------------------------
MAIN:
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: wave equation
WaveSolver: Solving the divergence pressure wave
WaveSolver: We have a harmonic or eigenmode system
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
CRS_IncompleteLU: ILU(2) (Real), Starting Factorization:
CRS_IncompleLU: Recursive round: 1
CRS_IncompleteLU: ILU(2) (Real), NOF nonzeros:   1102560
CRS_IncompleteLU: ILU(2) (Real), filling (%) :       456
CRS_IncompleteLU: ILU(2) (Real), Factorization ready at (s):     0.17
12 0.1228E-14 0.9235E-09
22 0.7611E-12 0.3690E-08
5 0.5304E-12 0.2564E-08
4 0.1043E-11 0.5115E-08
4 0.1915E-12 0.9567E-09
4 0.1027E-11 0.5092E-08
4 0.6112E-13 0.2997E-09
4 0.1125E-12 0.5514E-09
4 0.3809E-13 0.1910E-09
4 0.4597E-12 0.2274E-08
4 0.1939E-11 0.9609E-08
4 0.1317E-11 0.6608E-08
4 0.1958E-11 0.9892E-08
4 0.4634E-13 0.2300E-09
4 0.7754E-12 0.3925E-08
4 0.5258E-13 0.2682E-09
4 0.4891E-13 0.2444E-09
4 0.3098E-12 0.1575E-08
4 0.5467E-12 0.2729E-08
4 0.1925E-11 0.9782E-08
4 0.1083E-11 0.5455E-08
4 0.1909E-11 0.9570E-08
4 0.1098E-11 0.5615E-08
4 0.1266E-11 0.6381E-08
4 0.1625E-11 0.8228E-08
4 0.1846E-11 0.9409E-08
4 0.3496E-13 0.1789E-09
4 0.1573E-11 0.8068E-08
4 0.5069E-12 0.2593E-08
4 0.1596E-11 0.8240E-08
4 0.1605E-11 0.8212E-08
4 0.9570E-12 0.4932E-08
EigenSolve:
EigenSolve: Eigen system solution complete:
EigenSolve:
EigenSolve: Convergence criterion is:    1.000E-06
EigenSolve: Number of converged Ritz values is: 10
EigenSolve: Number of update iterations taken: 1
EigenSolve:
EigenSolve: Computed 10 Eigen Values
EigenSolve: --------------------------------
EigenSolve: 1:   -4.206542E-10   0.000000E+00
EigenSolve: 2:    4.648218E+04   0.000000E+00
EigenSolve: 3:    7.265766E+04   0.000000E+00
EigenSolve: 4:    1.192635E+05   0.000000E+00
EigenSolve: 5:    1.293377E+05   0.000000E+00
EigenSolve: 6:    1.760333E+05   0.000000E+00
EigenSolve: 7:    1.863612E+05   0.000000E+00
EigenSolve: 8:    2.023335E+05   0.000000E+00
EigenSolve: 9:    2.491492E+05   0.000000E+00
EigenSolve: 10:    2.595073E+05   0.000000E+00
EigenSolve: --------------------------------
ComputeChange: SS (ITER=1) (NRM,RELC): (  0.0000000      0.0000000     ) :: wave equation
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: internalvtuoutputsolver
ResultOutputSolver: -------------------------------------
ResultOutputSolver: Saving with prefix: case
ResultOutputSolver: Creating list for saving - if not present
CreateListForSaving: Field Variables for Saving
ResultOutputSolver: Saving in unstructured VTK XML (.vtu) format
VtuOutputSolver: Saving results in VTK XML format with prefix: case
VtuOutputSolver: Saving number of partitions: 1
ResultOutputSolver: -------------------------------------
ElmerSolver: *** Elmer Solver: ALL DONE ***
ElmerSolver: The end
SOLVER TOTAL TIME(CPU,REAL):        11.13        7.39
ELMER SOLVER FINISHED AT: 2020/08/21 11:59:44
``````
I will play with it some more and maybe be back with some questions.

Thanks!
CrocoDuck
Posts: 79
Joined: 12 May 2016, 13:15
Antispam: Yes

### Re: WaveSolver and Eigen Analysis

So, the eigenvalues appear to be angular frequency squared. If I calculate the associated frequency in Hz I get these values (omitting the first small negative value clearly associated with the 0 order mode, a constant pressure field with 0 eigenvalue):

Code: Select all

`````` 34.313370192521894
42.900349974920196
54.96344008822581
57.237769116156144
66.77556124149527
68.70650873280937
71.59027467208264
79.44194711791363
81.0764867267743
``````
They really agree well with what I expect:

Code: Select all

`````` 34.3
42.87500000000001
54.90679033598667
57.16666666666666
66.66721666439793
68.6
71.45833333333333
79.26401076641136
80.89638820738537
``````
I guess I will be able to raise the agreement even further with finer and/or higher order mesh.

I had a cursory look at the the associated eigenmodes and they look reasonable, but I still have to do a proper accuracy analysis like this one. I would say that everything is looking good so far.

I was wondering though why the manual does not mention the Eigen Analysis keywords. Is it maybe because the feature is still experimental?

Also, the manual does not mention the possibility of defining boundary conditions of this kind (which I used to get the solution above):

Code: Select all

``````Boundary Condition 1
Target Boundaries(6) = 2 3 4 5 6 7
Pressure Velocity = Real 0
End
``````
Where Pressure is the name I gave to the solver variable:

Code: Select all

``````Solver 1

! Equation
Equation = "Wave Equation"
Procedure = "WaveSolver" "WaveSolver"
Variable = "Pressure"

...
``````
I assume that the Pressure Velocity boundary condition is used to specify the particle velocity associated with the pressure field, and the results seems consistent with this interpretation. Is this correct?
raback
Posts: 4700
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

### Re: WaveSolver and Eigen Analysis

Hi

Great that you found success.

Now the Models Manual does not mention this possibility because it is a library feature mainly. As long as we take care that the time derivative multiplier matrix is registered even though the case is not transient, any time-dependent equation may become eigenvalue or harmonic system in this way. I would guess the "Eigen Analysis" will be found on the SolverManual side.

The same applies to "Pressure Velocity". That being said, I don't really recall what it would do in harmonic case. Does it play a role whether you have it or not? This is also taken care by the library.

From documentation point of view the challenge is that there are N modules and M library features and we don't really feel like documenting each of NxM combinations but try to get by with N+M.

-Peter