Can't I use a 2nd order mesh with the acoustic solver?

[tilbury]

I am using an acoustic solver (Linearized Navier-Stokes equations in the frequency domain), but the following error occurs for the 2nd order mesh.

ERROR:: ComplexLUDecomp: Matrix is singular.
ERROR:: ComplexInvertMatrix: inversion needs successfull LU-decomposition

For linear meshes, the results come out without any problems.

Below is a portion of the Python script used for automatic mesh generation in SALOME.

Mesh_1 = smesh.Mesh(Solid_1,'Mesh_1')
NETGEN_1D_2D_3D = Mesh_1.Tetrahedron(algo=smeshBuilder.NETGEN_1D2D3D)
NETGEN_3D_Parameters_1 = NETGEN_1D_2D_3D.Parameters()
NETGEN_3D_Parameters_1.SetMaxSize( 0.001 )
NETGEN_3D_Parameters_1.SetMinSize( 0.0001 )
NETGEN_3D_Parameters_1.SetSecondOrder( 1 )
NETGEN_3D_Parameters_1.SetOptimize( 1 )
NETGEN_3D_Parameters_1.SetFineness( 5 )

SetSecondOrder(1) creates a 2nd order tetra mesh and SetSecondOrder(0) creates a linear tetra mesh.

The sif file is as follows.

Code: Select all

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

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian 3D
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Steady state
  Steady State Max Iterations = 1
  Steady State Min Iterations = 1
  Output Intervals(1) = 1
  Solver Input File = case.sif
  Post File = case.vtu
  Frequency = Real MATC "f"
$ f = 1000.0
$ p = 1.2
End

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 1
  Target Bodies(1) = 1
  Name = "Body 1"
  Equation = 1
  Material = 1
End

Solver 1
  Equation = Linearized NavierStokes
  Variable = "Flow"
  Variable DOFs = 10
  Procedure = "Acoustics" "AcousticsSolver"
  Element = "p:1 b:1" 
  Bubbles in Global System = False
  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStabl
  BiCGstabl polynomial degree = 2
  Linear System Preconditioning = ILUT
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Max Iterations = 2000
  Linear System Convergence Tolerance = 1e-13
  Linear System Scaling = Logical False
  Linear System Row Equilibration = Logical True 
  Linear System Abort Not Converged = Logical True
  Linear System Residual Output = 10
End

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

Material 1
  Specific Heat = Real 718
  Specific Heat Ratio = 1.4
  Equilibrium Density = Real MATC "p"
  Equilibrium Temperature = 293
  Heat Conductivity = 0.0257
  Viscosity = 1.983e-5
  Bulk Viscosity = 5.98e-06
End

Boundary Condition 1
  Name = "Radiator"
  Target Boundaries(1) = 1 
  Re Surface Traction 1 = 1
!  Calculate Acoustic Impedance = Logical True
End

Boundary Condition 2
  Name = "Symmetry-y"
  Target Boundaries(2) = 2 5 
  Re Velocity 2 = 0
  Im Velocity 2 = 0
End

Boundary Condition 3
  Name = "Symmetry-z"
  Target Boundaries(2) = 3 4 
  Re Velocity 3 = 0
  Im Velocity 3 = 0
End

Boundary Condition 4
  Name = "Wall"
  Target Boundaries(15) = 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 
  Re Velocity 1 = 0
  Im Velocity 1 = 0
  Re Velocity 2 = 0
  Im Velocity 2 = 0
  Re Velocity 3 = 0
  Im Velocity 3 = 0
End

Is it not possible to use a 2nd order mesh in Linearized Navier-Stokes equations? Or is the solver setup incorrect?

[kevinarden]

If you run the test case
https://github.com/ElmerCSC/elmerfem/tr ... HarmonicNS

but try second order elements with
ElmerGrid 1 2 tube.grd -increase

than the same error appears.

[mika]

This is an old solver which was designed to use the lowest-order basis functions together with additional bubble basis functions to get a stable method. Therefore the second-order basis is not supported by this implementation. The given element definition of the form

Element = "p:1 b:n"

make sense only over the lowest-order mesh.

-- Mika

[tilbury]

This is an old solver which was designed to use the lowest-order basis functions together with additional bubble basis functions to get a stable method. Therefore the second-order basis is not supported by this implementation. The given element definition of the form

Element = "p:1 b:n"

make sense only over the lowest-order mesh.

-- Mika

Thank you for your answer.

I have found that sometimes 2nd order meshes work without error in the acoustic solver analysis of simple geometries. The same analysis results as using a linear mesh were obtained.

Should I assume that in the acoustic solver, even if there is a secondary mesh input, it is simply treated as a linear mesh?

[mika]

If you have a mesh consisting of second-order nodal elements (and avoid giving a conflicting element command "Element = p:..."), the solver should create a discretization based on second-order approximation of all fields. Such solution may however be unstable so it may suffer from spurious numerical oscillations especially when the wave number is small. So technically this might work, but the solution may not be good.

On the other hand, in this case combining the lowest-order mesh and trying to create the second-order discretization with the element command "Element = p:2" doesn't work at all, since the solver code hasn't been modernized to work with generic p-element definitions.