Boundary Layer Aeroelastic Stability Tool (BLAST)

A modular, object-oriented finite element solver for fluid-structure interaction stability problems.

Overview

This solver is designed to compute the global instability modes of fully coupled fluid-solid systems, handling arbitrary solid domain configurations.

• Monolithic coupling: Solves the fluid and solid domains simultaneously using exact linear algebra projections.
• Spatial and temporal analysis: Easily switch between temporal stability (tracking growth rates over time) and spatial stability (tracking spatial growth of perturbations).
• Dynamic polynomial eigenvalue solver: Features a custom Arnoldi-based shift-invert solver (neigs) that dynamically truncates polynomial degree based on the provided physics, automatically switching between generalized, quadratic, or higher-order polynomial eigenvalue problems.

---

Mathematical formulation

This solver linearizes the coupled FSI governing equations into a polynomial eigenvalue problem of the form:

[A_0 + lambda * A_1 + lambda^2 * A_2 + ...] * {x} = 0

The solver maps physical parameters to the non-dimensional eigenvalue (lambda) depending on the type of analysis:

• Temporal stability: * Target: Frequency (omega)

  • Mapping: lambda = i * omega_nondim

• Spatial stability: * Target: Wavenymber (kappa)

  • Mapping: lambda = i * kappa_nondim

The neigs solver automatically linearizes polynomials of degree N > 1 into a first-order companion matrix format and uses scipy.sparse.linalg.eigs (ARPACK) alongside scipy.sparse.linalg.factorized (UMFPACK) for highly efficient shift-invert spectral extraction.

---

Code architecture

The framework is modularized into several core components:

• fem.physics: This module contains the FluidProperties and SolidProperties objects, which are used to define flow and solid material properties, respectively.
• fem.mesh: It contains the objects FluidConfig and SolidConfig that allow the user to setup the problem (e.g., domain sizes, discretization parameters, boundary conditions, etc.). The module also contains the functions to generate meshes and an object Mesh to manage mesh data, including nodal coordinates, connectivities, and different node subsets (e.g., boundary parts).
• fem.dofs: It handles the degrees of freedom for the fluid domain (FluidDofHandler), the solid domains (SolidDofHandler) and the coupled problem (CoupledDofHandler).
• fem.bcs: This module handles the generation of Boolean projection matrices for varioius boundary conditions (e.g., precribed degrees of freedom or periodic boundary conditions). There are dedicated objects for the uncoupled fluid and solid problems (FluidBoundaryConditions and SolidBoundaryConditions, respectively) and a FSIBoundaryConditions that handles the coupling across the fluid and solid domains.
• fem.elements: It contains the objects to produce that interpolation data for different finite-element types (e.g., ElementL3 for quadratic line elements, or ElementQ9 for biquadratic rectangular elements).
• fem.assembly: It contains the objects to generate each system subdomain matrices (uncoupled fluid and solid matrices via the objects FluidMatrices and SolidMatrices).
• fem.solver: It contains the StabilitySolver object and neigs function that act as the main physics engine. It applies boundary projections, nondimensionalizes the system, assembles the temporal/spatial block matrices, and drives the eigenvalue solver.

---

Quick example

Minimal working example of how to initialize the system and find the unstable mode for a channel flow at a given Reynolds:

from fem.physics import FluidProperties
from fem.mesh import FluidConfig, generate_fluid_mesh
from fem.dofs import FluidDofHandler
from fem.elements import ElementL3
from fem.assembly import FluidMatrices
from fem.bcs import FluidBoundaryConditions
from fem.solver import StabilitySolver

# 1. Define physics
prop = FluidProperties(Re=10000, mu=0.0001, base_flow="poiseuille")

# 2. Define channel flow domain and compute mesh
config = FluidConfig(height=2.0, ny=96, spacing="geometric", growth_rate=(1.05, 0), boundaries={"bottom": "wall", "top": "wall"})
mesh = generate_fluid_mesh(config)

# 3. Define the DOFs handler
dofs = FluidDofHandler(mesh=mesh)

# 4. Define the element type
elem = ElementL3(n_gauss=3)

# 5. Assemble matrices
mats = FluidMatrices(element=elem, dofs=dofs, properties=prop).assemble()

# 6. Compute boolean projection matrices for applying boundary conditions
bcs = FluidBoundaryConditions(dofs=dofs)

# 7. Solve for temporal stability
input_k = 1.0
target_omega = 0.25
eigvals, eigvecs, _ = StabilitySolver(dofs=dofs, properties=prop).solve_temporal(input_k, bcs, mats, target_w=target_omega)
print(f"Unstable mode frequency: {eigvals[0]:.5f}")

---

Installation

We recommend using Conda to manage dependencies. To set up the environment, run:

git clone [https://github.com/david-roca/BLAST.git](https://github.com/david-roca/BLAST.git)
cd BLAST
conda env create -f environment.yml
conda activate fsi_env