Skip to main content

Dynamic Modal Analysis

Pipelines are highly susceptible to flow-induced vibration, acoustic resonance, and seismic excitation. To ensure a design avoids resonant frequencies, the Pipe Stress FEA engine includes a dedicated Modal Analysis solver to extract the system's fundamental natural frequencies and mode shapes.

4.2.1 The Consistent Mass Matrix

To achieve high-fidelity dynamic results without requiring an excessively dense mesh, the engine utilizes a Consistent Mass Matrix [M][M].

Unlike a "lumped mass" approach—which simply divides element mass into point loads at the nodes—the consistent mass matrix distributes the mass continuously along the cubic shape functions of the Euler-Bernoulli beam element. This captures the true dynamic coupling between translational and rotational degrees of freedom.

The total mass per unit length (μ\mu) and the polar mass moment of inertia per unit length (JρJ_\rho) are calculated for each element:

μ=(ρsteelA)+mcoat+mfluid+mvalve\mu = (\rho_{steel} A) + m_{coat} + m_{fluid} + m_{valve}

Jρ=μ(ro2+ri22)J_\rho = \mu \left( \frac{r_o^2 + r_i^2}{2} \right)

The resulting 12×1212 \times 12 local mass matrix [Mlocal][M_{local}] contains the following distinct distribution coefficients (where m=μLm = \mu L and mtor=JρLm_{tor} = J_\rho L):

  • Axial & Torsional Mass: Distributed using linear shape functions, resulting in diagonal terms of 13m\frac{1}{3}m and off-diagonal coupling terms of 16m\frac{1}{6}m.
  • Bending Mass: Distributed using cubic Hermitian shape functions. The translational inertia terms are scaled by 156420m\frac{156}{420}m, while the rotational inertia terms (capturing the rotary momentum of the bending pipe) are scaled by 4420L2m\frac{4}{420}L^2 m.

Like the stiffness matrix, this local mass matrix is rotated using the element's direction cosine matrix [T][T] and assembled into the Global Mass Matrix [Mglobal][M_{global}].

4.2.2 Dynamic Boundary Conditions

Modal analysis is fundamentally a linear mathematical operation, meaning the non-linear iterative penalty solver used for standard static loads cannot be used. The engine must linearize the boundary conditions.

During the Modal state resolution, the engine applies the following boundary logic:

  1. Anchors: Fully stripped from the matrix to prevent singularities.
  2. Limit Stops & Gaps: Assumed to be Active (Rigid). Gaps are ignored, and the support is locked using the penalty stiffness (kp=1012k_p = 10^{12}).
  3. Snubbers: Dynamic snubbers instantly lock up and are treated as rigid supports along their active axis.
  4. Friction: Ignored.
Massless DOFs & Matrix Stability

In purely torsional degrees of freedom where no eccentric mass exists, the diagonal mass matrix entry can theoretically be zero. To prevent the eigenvalue solver from crashing on a singular matrix, the engine injects a mathematically insignificant stabilizer (106 kg10^{-6} \text{ kg}) to all diagonal free DOFs.

4.2.3 Eigenvalue Extraction

The natural frequencies of the undamped pipeline system are found by solving the generalized eigenvalue problem:

([K]ω2[M]){ϕ}=0([K] - \omega^2 [M])\{\phi\} = 0

Where:

  • [K][K] is the reduced global stiffness matrix.
  • [M][M] is the reduced global mass matrix.
  • ω2\omega^2 represents the eigenvalues.
  • {ϕ}\{\phi\} represents the eigenvectors (mode shapes).

The engine utilizes a highly optimized SciPy LAPACK solver (scipy.linalg.eigh) to extract the lowest 10 eigenvalues, as these represent the most critical fundamental frequencies governing pipeline design.

The angular velocity (ω\omega) is converted into the physical frequency (ff) in Hertz (Hz):

f=ω22πf = \frac{\sqrt{|\omega^2|}}{2\pi}

Mode Shape Normalization

The raw eigenvectors {ϕ}\{\phi\} representing the mode shapes indicate relative displacement, not actual physical deflection. To optimize the 3D visualizer in the frontend, the engine normalizes every mode shape so that the maximum nodal displacement vector is scaled exactly to 1.01.0.

4.2.4 Mass Participation

For each mode the engine reports the modal participation factor per global translational direction and the corresponding effective modal mass. Because eigh returns mass-orthonormal eigenvectors (ϕiT[M]ϕi=1\phi_i^T [M] \phi_i = 1), the participation factor against the rigid-body influence vector {rd}\{r_d\} (unit translation in direction dd) is:

Γd,i=ϕiT[M]{rd},meff,d,i=Γd,i2\Gamma_{d,i} = \phi_i^T [M] \{r_d\}, \qquad m_{eff,\,d,i} = \Gamma_{d,i}^2

The results table lists, per mode: frequency, period, Γx/Γy/Γz\Gamma_x/\Gamma_y/\Gamma_z, effective mass, and the cumulative participation percentage per direction.

Interpretation: a mode with high X-participation is the one an X-direction seismic or wind event excites most strongly. If the cumulative captured mass in a direction is well below ~90 % across the extracted modes, significant response lives in higher modes — stiffen the system, or treat the static-equivalent occasional results with extra conservatism.