Acoustic Analysis Using the Boundary Element Method (BEM)
Theory and Physics
What is Acoustic BEM?
Professor, what is acoustic BEM?
A method that transforms the Helmholtz equation into a Boundary Integral Equation and calculates the sound field using only surface meshes. Its greatest advantage is that it naturally handles external sound fields (infinite domains) because volume meshes are unnecessary.
Governing Equation
Helmholtz equation for the acoustic field:
$k = \omega/c$: wave number. This is transformed into a boundary integral form using Green's theorem:
$G$: free-space Green's function:
$c(\mathbf{x})$: $1/2$ on the boundary, $1$ inside the domain, $0$ outside the domain.
What are its advantages compared to FEM?
Summary
BEM's Origins are in 1960s Elasticity Theory
The mathematical foundation of the Boundary Element Method was laid by Jaswon (1963) and Sympson from Stanford University. They published a method for discretizing integral equations for elasticity problems, but at the time, applications to acoustics were not envisioned. The adaptation to acoustic BEM came about a decade later, with Schenck's CHEFS method published in 1968 as an IBM technical report being the pioneer, influencing later NASTRAN and commercial solvers.
Physical Meaning of Each Term
- Inertia Term (Mass Term): $\rho \ddot{u}$, i.e., "mass × acceleration". Have you ever experienced being thrown forward during sudden braking? That "feeling of being carried away" is precisely the inertial force. Heavier objects are harder to set in motion and harder to stop once moving. Buildings shake during earthquakes because the ground moves suddenly while the building's mass "gets left behind". In static analysis, this term is set to zero, which assumes "forces are applied slowly enough that acceleration can be ignored". It absolutely cannot be omitted for impact loads or vibration problems.
- Stiffness Term (Elastic Restoring Force): $Ku$ or $\nabla \cdot \sigma$. When you stretch a spring, you feel a "force trying to return it", right? That's Hooke's law $F=kx$, and it's the essence of the stiffness term. So here's a question—if you pull an iron rod and a rubber band with the same force, which stretches more? Obviously the rubber. This "resistance to stretching" is the Young's modulus $E$, which determines stiffness. A common misconception: "high stiffness = strong". Stiffness is "resistance to deformation", strength is "resistance to failure"—they are different concepts.
- External Force Term (Load Term): Body force $f_b$ (e.g., gravity) and surface force $f_s$ (e.g., pressure, contact force). Think of it this way—the weight of a truck on a bridge is a "force acting on the entire volume" (body force), while the force of the tires pushing on the road surface is a "force acting only on the surface" (surface force). Wind pressure, water pressure, bolt tightening force... all are external forces. A common mistake here: getting the load direction wrong. Intending "tension" but ending up with "compression"—it sounds like a joke, but it actually happens when coordinate systems are rotated in 3D space.
- Damping Term: Rayleigh damping $C\dot{u} = (\alpha M + \beta K)\dot{u}$. Try plucking a guitar string. Does the sound continue forever? No, it gradually fades away. That's because the vibration energy is converted into heat by air resistance and internal friction in the string. Car shock absorbers work on the same principle—they deliberately absorb vibration energy to improve ride comfort. What if damping were zero? Buildings would continue shaking forever after an earthquake. Since that doesn't happen in reality, setting appropriate damping is crucial.
Assumptions and Applicability Limits
- Continuum assumption: Treats material as a continuous medium, ignoring microscopic inhomogeneity.
- Small deformation assumption (for linear analysis): Deformation is sufficiently small compared to initial dimensions, and the stress-strain relationship is linear.
- Isotropic material (unless otherwise specified): Material properties are independent of direction (anisotropic materials require separate tensor definitions).
- Quasi-static assumption (for static analysis): Ignores inertial and damping forces, considering only the balance between external and internal forces.
- Non-applicable cases: Large deformation/large rotation problems require geometric nonlinearity. Nonlinear material behavior such as plasticity or creep requires constitutive law extensions.
Dimensional Analysis and Unit Systems
| Variable | SI Unit | Notes / Conversion Memo |
|---|---|---|
| Displacement $u$ | m (meter) | When inputting in mm, unify load and elastic modulus to MPa/N system |
| Stress $\sigma$ | Pa (Pascal) = N/m² | MPa = 10⁶ Pa. Be careful of unit system inconsistency when comparing with yield stress |
| Strain $\varepsilon$ | Dimensionless (m/m) | Note the distinction between engineering strain and logarithmic strain (for large deformation) |
| Elastic Modulus $E$ | Pa | Steel: ~210 GPa, Aluminum: ~70 GPa. Note temperature dependence |
| Density $\rho$ | kg/m³ | In mm system: tonne/mm³ (= 10⁻⁹ tonne/mm³ for steel) |
| Force $F$ | N (Newton) | Unify as N in mm system, N in m system |
Numerical Methods and Implementation
Discretization of Acoustic BEM
How do you solve the boundary integral equation numerically?
Discretize the surface into triangular or quadrilateral elements. Approximate sound pressure $p$ and normal velocity $v_n$ with nodal values.
$[H]$, $[G]$: Influence matrices. Integral values of the Green's function between each element pair.
Handling Singular Integrals
The most technically challenging part of BEM implementation is handling singular integrals (where the Green's function diverges as $r \to 0$):
- Weak singularity ($1/r$): Regularize using polar coordinate transformation
- Strong singularity ($1/r^2$): Define via Cauchy principal value
- Hypersingular integral: Hadamard finite-part integral
Non-Uniqueness Problem
In external BEM, the solution becomes non-unique at the internal eigenfrequencies of the closed boundary. Countermeasures:
- CHIEF method: Place additional equations at interior points (overdetermined system)
- Burton-Miller method: Linear combination of the standard BIE and its normal derivative BIE. Theoretically complete.
The Burton-Miller method is standard. It requires handling hypersingular integrals but reliably eliminates non-uniqueness.
FMM (Fast Multipole Method)
You mentioned that the dense matrix makes large-scale problems unsolvable...
Solved by FMM (Fast Multipole Method). Approximates calculations for groups of distant elements collectively:
- Memory: $O(N^2) \to O(N)$
- Computational cost: $O(N^2) \to O(N\log N)$
You mentioned that the dense matrix makes large-scale problems unsolvable...
Solved by FMM (Fast Multipole Method). Approximates calculations for groups of distant elements collectively:
BEM with 1 million elements has become practical. Implemented in Actran, COMSOL, etc.
Summary
FMM Reduces BEM Computational Cost from O(N²) to O(N log N)
Traditional BEM required O(N²) memory and computational cost relative to the number of nodes N, making large-scale models practically impossible. The Fast Multipole Method (FMM) published by Greengard and Rokhlin in 1987 brought about a revolution, making acoustic analysis of automotive body scale (hundreds of thousands of nodes) realistic. Nuances' VirtualLab Acoustics and LMS Sysnoise 5.x are known as products that implemented FMM-BEM early on.
Linear Elements (1st Order Elements)
Linear interpolation between nodes. Low computational cost but low stress accuracy. Beware of shear locking (mitigated by reduced integration or B-bar method).
Quadratic Elements (with Mid-Side Nodes)
Can represent curved deformation. Stress accuracy improves significantly, but degrees of freedom increase by about 2-3 times. Recommended: when stress evaluation is important.
Full Integration vs Reduced Integration
Full Integration: Risk of over-constraint (locking). Reduced Integration: Risk of hourglass modes (zero-energy modes). Choose appropriately for the situation.
Adaptive Mesh
Automatic refinement based on error indicators (e.g., ZZ estimator). Efficiently improves accuracy in stress concentration areas. Includes h-method (element subdivision) and p-method (order increase).
Newton-Raphson Method
Standard method for nonlinear analysis. Updates the tangent stiffness matrix each iteration. Achieves quadratic convergence within the convergence radius, but computational cost is high.
Modified Newton-Raphson Method
Updates the tangent stiffness matrix using the initial value or every few iterations. Cost per iteration is low, but convergence speed is linear.
Convergence Criteria
Force residual norm: $||R|| / ||F_{ext}|| < \epsilon$ (typically $\epsilon = 10^{-3}$ to $10^{-6}$). Displacement increment norm: $||\Delta u|| / ||u|| < \epsilon$. Energy norm: $\Delta u \cdot R < \epsilon$
Load Increment Method
Instead of applying the full load at once, apply it in small increments. The arc-length method (Riks method) can trace beyond limit points on the load-displacement curve.
Analogy: Direct Method vs Iterative Method
The direct method is like "solving simultaneous equations accurately with pen and paper"—reliable but takes too long for large-scale problems. The iterative method is like "repeatedly guessing to approach the correct answer"—starts with a rough answer but accuracy improves with each iteration. It's the same principle as looking up a word in a dictionary: it's more efficient to open it at an estimated location and adjust forward/backward (iterative method) than to search sequentially from the first page (direct method).
Relationship Between Mesh Order and Accuracy
1st order elements are like "approximating a curve with a ruler"—represented by straight line segments, so accuracy is limited. 2nd order elements are like a "flexible curve"—can represent curved changes, dramatically improving accuracy even at the same mesh density. However, computational cost per element increases, so judge based on total cost-effectiveness.
Practical Guide
Acoustic BEM in Practice
Typical applications include engine radiated noise, tire noise, transformer noise, and acoustic radiation from exhaust pipes.
Analysis Flow
1. Structural Vibration Analysis — Obtain surface vibration velocity $v_n$ via FEM
2. Create BEM Surface Mesh — Extract from FEM mesh or create independently
3. Set Boundary Conditions — Set $v_n$ (Neumann condition)
4. BEM Solution — Calculate surface sound pressure $p$
5. Sound Field Evaluation — Calculate sound pressure, radiated power at arbitrary observation points
Practical Checklist
BEM Mesh Guidelines
| Maximum Frequency [Hz] | Wavelength at 340m/s [m] | Maximum Element Size [mm] |
|---|---|---|
| 500 | 0.68 | 113 |
| 1000 | 0.34 | 57 |
| 2000 | 0.17 | 28 |
| 5000 | 0.068 | 11 |
| 10000 | 0.034 | 5.7 |
The mesh gets quite fine above 5kHz.