Source notebook: examples/thermoacoustics/ita_and_cavity_modes.ipynb. This page is executed from it at build time.
A thermoacoustic network with an active flame carries two families of modes, and they are easy to confuse because a single eigenvalue solve returns them in one list with no label for the family.
Cavity modes (also called acoustic modes) are the resonances of the enclosure. They are set by the duct lengths, the sound speeds, and the reflection at each end. The flame shifts them and can feed or drain them, but they exist without it.
Intrinsic thermoacoustic (ITA) modes are set by the flame alone. The flame’s velocity-to-heat-release feedback closes a loop through the flame’s own acoustic jump: an unsteady heat release drives an acoustic wave, part of which travels back upstream to the flame’s reference station, where after the time lag \tau it modulates the heat release again. Nothing in that loop needs the ends of the tube. Remove the reflections entirely and the ITA modes remain; remove the flame and they disappear.
This notebook builds a Rijke tube that carries three of each below 1\ \mathrm{kHz}, then separates them two ways:
The companion notebook rijke_tube.ipynb treats the simpler tube whose fundamental is the only mode of interest. Here the point is that the spectrum has structure worth reading.
import warnings
import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from scipy.optimize import brentq
import nefes
from nefes.elements import catalog as cat
from nefes.elements import n_tau_flame
from nefes.perturbation import PerturbationBC, eigenvalue_trajectory
from nefes.plotting import COLORWAY, use_nefes_theme
use_nefes_theme()'nefes'Cold air enters through a mass-flow inlet, travels a cold duct of length L_1, crosses an acoustically compact flame, and leaves through a hot duct of length L_2.
Two choices shape everything that follows.
The mass flow is deliberately tiny, \dot m = 10^{-3}\ \mathrm{kg/s}, giving \overline{M}_1 \approx 2.5\times10^{-4}. The tube is near-quiescent, so the closed forms we compare against, which drop terms of order \overline{M}, are accurate to a fraction of a percent. This is a convenience for the verification, not a restriction of the solver.
The flame raises the total temperature by \Delta\overline{T} = \dot Q/(\dot m\,c_p) = 1000\ \mathrm{K}, so the hot duct is much longer than the cold one and carries a much faster sound speed. The resulting cavity resonances are irregularly spaced, which is what lets them interleave with the ITA ladder rather than sit on top of it.
The ends argument selects the acoustic closure at both terminals: a pressure-release open end (\widehat{p}=0, R=-1) for the physical tube, or a reflection-free termination (R=0) for the experiment in section 4. The mean flow is identical either way; only the linear acoustics change.
# Perfect gas
R, GAMMA = 287.0, 1.4
CP = GAMMA * R / (GAMMA - 1.0)
# Geometry: a uniform cross-section, a short cold duct and a long hot one
AREA = 0.01
L_COLD, L_HOT = 0.25, 0.75
# Mean operating point. The mass flow is small enough that the tube is acoustically quiescent.
MDOT = 1.0e-3
T_IN, P_OUT = 300.0, 1.0e5
# Flame heat release, sized by the total-temperature rise it must produce
DT_FLAME = 1000.0
QDOT = MDOT * CP * DT_FLAME
# Nominal flame transfer function and the frequency band
N, TAU = 1.0, 3.0e-3
F_MAX = 1000.0
def rijke(n, tau, *, ends="open"):
"""Build and solve the Rijke tube with a compact n-tau flame.
Parameters
----------
n : float
Interaction index (gain) of the n-tau flame transfer function.
tau : float
Flame time lag [s].
ends : {"open", "anechoic"}, optional
Acoustic closure applied at both terminals. ``"open"`` is the physical
pressure-release end (``R = -1``) that supports cavity modes; ``"anechoic"``
(``R = 0``) suppresses them, leaving only the intrinsic modes. The mean flow is
unaffected by this choice.
Returns
-------
Solution
The converged mean-flow solution; read its edge states by name and call the acoustic
analyses (``eigenmodes``, ...) as methods on it.
"""
bc = PerturbationBC.anechoic() if ends == "anechoic" else PerturbationBC.open_end()
net = nefes.Network(nefes.perfect_gas(R, GAMMA))
inlet = net.add(cat.mass_flow_inlet(MDOT, T_IN, perturbation_bc=bc))
cold = net.add(cat.duct(L_COLD, name="cold duct"))
flame = net.add(cat.heat_release_flame(QDOT))
hot = net.add(cat.duct(L_HOT, name="hot duct"))
outlet = net.add(cat.pressure_outlet(P_OUT, perturbation_bc=bc))
# cold -> flame edge id: reference station for the flame transfer function
net.connect(inlet, cold, area=AREA)
ref = net.connect(cold, flame, area=AREA)
net.connect(flame, hot, area=AREA)
net.connect(hot, outlet, area=AREA)
net.set_dynamic_source(flame, n_tau_flame(n, tau, ref_edge=ref))
sol = net.solve()
assert sol.converged
return sol
# The mean flow ignores the flame transfer function, so one solve fixes the mean state everywhere.
sol_passive = rijke(0.0, 0.0)
sol_passive.print_states()| edge | mdot [kg/s] | p [Pa] | h_t [J/kg] | rho [kg/m^3] | u [m/s] | T [K] | c [m/s] | M [-] | p_t [Pa] | area [m^2] |
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0.001 | 1e+05 | 3.0135e+05 | 1.1614 | 0.0861 | 300 | 347.19 | 0.00024799 | 1e+05 | 0.01 |
| 1 | 0.001 | 1e+05 | 3.0135e+05 | 1.1614 | 0.0861 | 300 | 347.19 | 0.00024799 | 1e+05 | 0.01 |
| 2 | 0.001 | 1e+05 | 1.3059e+06 | 0.26802 | 0.3731 | 1300 | 722.73 | 0.00051624 | 1e+05 | 0.01 |
| 3 | 0.001 | 1e+05 | 1.3058e+06 | 0.26802 | 0.3731 | 1300 | 722.73 | 0.00051624 | 1e+05 | 0.01 |
# Mean state either side of the compact flame: edge 1 is the cold approach, edge 2 the hot products.
cold_edge, hot_edge = sol_passive.edge(1), sol_passive.edge(2)
RHO1, C1, T1 = cold_edge["rho"], cold_edge["c"], cold_edge["T"]
RHO2, C2, T2 = hot_edge["rho"], hot_edge["c"], hot_edge["T"]
Z1, Z2 = RHO1 * C1, RHO2 * C2 # characteristic impedances
print(f"cold: T = {T1:7.2f} K c = {C1:6.2f} m/s Z = {Z1:6.2f} kg/(m^2 s)")
print(f"hot: T = {T2:7.2f} K c = {C2:6.2f} m/s Z = {Z2:6.2f} kg/(m^2 s)")
print(f"impedance ratio Z1/Z2 = {Z1 / Z2:.4f} (= sqrt(T2/T1) = {np.sqrt(T2 / T1):.4f} at uniform mean pressure)")cold: T = 300.00 K c = 347.19 m/s Z = 403.24 kg/(m^2 s)
hot: T = 1300.00 K c = 722.73 m/s Z = 193.71 kg/(m^2 s)
impedance ratio Z1/Z2 = 2.0817 (= sqrt(T2/T1) = 2.0817 at uniform mean pressure)With the flame inert (n=0) the network is acoustically passive. Its modes are the resonances of a two-region tube: plane waves in each uniform duct, pressure and velocity continuous across the compact flame, and \widehat{p}=0 at both open ends. Eliminating the four wave amplitudes leaves a single real dispersion relation,
\overline{Z}_1\,\tan\!\left(\frac{2\pi f L_1}{\overline{c}_1}\right) \;+\; \overline{Z}_2\,\tan\!\left(\frac{2\pi f L_2}{\overline{c}_2}\right) \;=\; 0, \qquad \overline{Z}_j=\overline{\varrho}_j\,\overline{c}_j ,
where \overline{Z}_j is the characteristic impedance of region j and \overline{c}_j its sound speed. The temperature jump enters only through \overline{c}_2/\overline{c}_1 and \overline{Z}_1/\overline{Z}_2; the roots are pure geometry and mean state.
Read what the relation does not contain: neither n nor \tau appears. That is the defining property of the cavity family.
We solve it by scanning for sign changes and polishing with Brent’s method. Multiplying through by \cos(\cdot)\cos(\cdot) first removes the tangent poles, so every sign change on the scan is a genuine root rather than a pole crossing.
def cavity_residual(f):
"""Pole-free form of the two-region cavity dispersion relation; its zeros are the resonances."""
k1, k2 = 2.0 * np.pi * f / C1, 2.0 * np.pi * f / C2
return Z1 * np.sin(k1 * L_COLD) * np.cos(k2 * L_HOT) + Z2 * np.cos(k1 * L_COLD) * np.sin(k2 * L_HOT)
scan = np.linspace(1.0, F_MAX, 20001)
vals = cavity_residual(scan)
cavity_analytic = np.array(
[brentq(cavity_residual, scan[i], scan[i + 1]) for i in range(scan.size - 1) if vals[i] * vals[i + 1] < 0.0]
)
# Nefes eigenmodes in the same frequency band (Beyn contour search).
eig_passive = sol_passive.eigenmodes(freq_band=(20.0, F_MAX), growth_band=(-600.0, 600.0), isentropic=True)
eig_passive| # | f [Hz] | growth [1/s] | damping ratio | residual | stability |
|---|---|---|---|---|---|
| 0 | 268.150 | -0.417 | +0.0002 | 8.3e-17 | stable |
| 1 | 599.015 | -0.586 | +0.0002 | 1.1e-16 | stable |
| 2 | 820.945 | -0.541 | +0.0001 | 2.2e-16 | stable |
order = np.argsort(eig_passive.freqs)
print(f"{'analytic [Hz]':>15}{'Nefes [Hz]':>14}{'growth [1/s]':>15}")
for f_an, i in zip(cavity_analytic, order):
print(f"{f_an:>15.3f}{eig_passive.freqs[i]:>14.3f}{eig_passive.growth_rates[i]:>15.3f}") analytic [Hz] Nefes [Hz] growth [1/s]
268.150 268.150 -0.417
599.015 599.015 -0.586
820.945 820.945 -0.541The frequencies agree to the printed precision. The growth rates are not exactly zero because the network keeps the O(\overline{M}) convective terms that the zero-Mach dispersion relation drops. At \overline{M}_1\approx 2.5\times10^{-4} that residual damping is a few tenths of \mathrm{s}^{-1}, which is the size of the discrepancy we should expect everywhere below.
Attaching the n-\tau flame transfer function \mathcal{F}(\omega)=n\,e^{-\mathrm{i}\omega\tau} adds the source block \mathbf{S}(\omega) to the perturbation operator \mathbf{A}(\omega)=\overline{\mathbf{J}}+\mathrm{i}\omega\mathbf{M}+\mathbf{P}(\omega)+\mathbf{S}(\omega).
The mode count goes from three to six. The three extra roots were not there before, and they are not the passive resonances moved around: the operator has genuinely acquired new solutions.
Both spectra are reported as certified complete: Nefes checks that the number of modes found matches the number predicted inside the search region (Beyn’s argument-principle count). That check matters here because the whole argument rests on a mode count. Without it we could not tell a missing mode from an incomplete search.
sol_active = rijke(N, TAU)
eig_active = sol_active.eigenmodes(freq_band=(20.0, F_MAX), growth_band=(-600.0, 600.0), isentropic=True)
eig_active| # | f [Hz] | growth [1/s] | damping ratio | residual | stability |
|---|---|---|---|---|---|
| 0 | 162.641 | +155.499 | -0.1522 | 1.1e-16 | unstable |
| 1 | 306.089 | -62.566 | +0.0325 | 6.2e-17 | stable |
| 2 | 464.186 | -355.322 | +0.1218 | 8.4e-17 | stable |
| 3 | 578.134 | +282.047 | -0.0776 | 1.3e-16 | unstable |
| 4 | 796.847 | +362.551 | -0.0724 | 1.4e-16 | unstable |
| 5 | 921.940 | -292.255 | +0.0505 | 8.3e-16 | stable |
# The two candidate reference sets: the passive cavity resonances (section 2) and the flame's own
# ladder f_k = (2k+1)/(2 tau), derived in section 4. Only the rungs below F_MAX are in band.
n_rungs = int(np.ceil((2.0 * TAU * F_MAX - 1.0) / 2.0))
ladder = np.array([(2 * k + 1) / (2.0 * TAU) for k in range(n_rungs)])
fig = go.Figure()
for f in cavity_analytic:
fig.add_vline(x=f, line=dict(color="#9aa5b1", width=1, dash="dot"))
for f in ladder:
fig.add_vline(x=f, line=dict(color="#cbd2d9", width=1, dash="dash"))
fig.add_scatter(
x=eig_passive.freqs,
y=eig_passive.growth_rates,
mode="markers",
name="passive flame (n = 0)",
marker=dict(size=13, symbol="circle-open", color="#52606d", line=dict(width=2)),
)
fig.add_scatter(
x=eig_active.freqs,
y=eig_active.growth_rates,
mode="markers",
name=f"active flame (n = {N})",
marker=dict(size=11, color=COLORWAY[0]),
)
fig.add_hline(y=0.0, line=dict(color="#9aa5b1", width=1.4, dash="dash"))
fig.add_annotation(
x=cavity_analytic[0],
y=520,
text="cavity resonances",
showarrow=False,
xanchor="left",
font=dict(color="#52606d", size=12),
)
fig.add_annotation(
x=ladder[0],
y=-520,
text="ITA ladder (2k+1)/(2τ)",
showarrow=False,
xanchor="left",
font=dict(color="#52606d", size=12),
)
fig.update_layout(
title=f"Six modes where there were three (n = {N}, τ = {TAU * 1e3:.0f} ms)",
xaxis_title="frequency [Hz]",
yaxis_title="growth rate −Im(ω) [1/s]",
width=950,
height=440,
)
fig.show()Three of the six modes are unstable. Each active mode sits near either a cavity resonance or a ladder rung, but proximity is only a visual hint at n=1: the flame has pulled every eigenvalue several tens of Hz, and the 922\ \mathrm{Hz} mode is nowhere near the 833\ \mathrm{Hz} rung it will eventually approach.
The next two sections replace proximity with tests that have a definite answer.
Set both terminals to R=0. No wave that leaves the domain ever comes back, so there is no enclosure and there can be no enclosure resonance. Whatever survives owes nothing to the ends.
For this configuration the problem reduces to a closed form. Anechoic ends mean the only wave entering the cold region is the one the flame sends upstream, and the only wave in the hot region is the one it sends downstream. Writing \widehat{p} continuous and the compact velocity jump \widehat{u}_2-\widehat{u}_1 = \big[(\gamma-1)/(\gamma\,\overline{p}\,A)\big]\,\widehat{q} with \widehat{q}=\overline{Q}\,\mathcal{F}(\omega)\,\widehat{u}_1/\overline{u}_1, the wave amplitudes cancel and leave a single scalar condition:
K\,e^{-\mathrm{i}\omega\tau} \;=\; -B, \qquad K \;=\; n\,\frac{\Delta\overline{T}}{\overline{T}_1}, \qquad B \;=\; 1+\frac{\overline{Z}_1}{\overline{Z}_2},
where K is the flame’s steady gain, the interaction index times the fractional temperature rise, and B is the acoustic load the two anechoic duct sides place on the compact flame. Splitting \omega=\omega_r+\mathrm{i}\omega_i into phase and magnitude gives the two families of solution directly:
f_k \;=\; \frac{\omega_r}{2\pi} \;=\; \frac{2k+1}{2\tau}, \qquad \sigma \;=\; -\omega_i \;=\; \frac{1}{\tau}\,\ln\!\frac{K}{B}, \qquad k=0,1,2,\dots
Three things follow, and all three are checkable.
The frequencies contain no L_1, no L_2, no \overline{c}. They are an odd-harmonic ladder whose only clock is the flame lag \tau; the phase condition \omega_r\tau=(2k+1)\pi is the loop closing on itself with a half-cycle of delay. They also do not depend on n.
The growth rate is the same on every rung, because the loop gain K/B is frequency-independent.
Instability requires K>B, that is n > n_{\mathrm{crit}} = \big(1+\overline{Z}_1/\overline{Z}_2\big)\,\overline{T}_1/\Delta\overline{T}, independent of \tau. As n\to 0 we get \sigma\to-\infty: the ITA modes leave the spectrum when the flame is switched off, which is exactly the fingerprint section 5 exploits.
# Passive and active tube, both with reflection-free ends.
sol_an_passive = rijke(0.0, TAU, ends="anechoic")
eig_an_passive = sol_an_passive.eigenmodes(freq_band=(20.0, F_MAX), growth_band=(-1200.0, 600.0), isentropic=True)
print(f"anechoic ends, flame off: {eig_an_passive.n_modes} modes <-- no enclosure, no cavity resonance\n")
sol_an = rijke(N, TAU, ends="anechoic")
eig_an = sol_an.eigenmodes(freq_band=(20.0, F_MAX), growth_band=(-1200.0, 600.0), isentropic=True)
eig_ananechoic ends, flame off: 0 modes <-- no enclosure, no cavity resonance
| # | f [Hz] | growth [1/s] | damping ratio | residual | stability |
|---|---|---|---|---|---|
| 0 | 166.667 | +26.028 | -0.0249 | 1.6e-16 | unstable |
| 1 | 500.000 | +26.028 | -0.0083 | 1.6e-16 | unstable |
| 2 | 833.333 | +26.028 | -0.0050 | 3.2e-16 | unstable |
K = N * (T2 - T1) / T1 # flame steady gain
B = 1.0 + Z1 / Z2 # anechoic acoustic load on the compact flame
sigma_analytic = np.log(K / B) / TAU
n_crit = B * T1 / (T2 - T1)
print(f"K = {K:.4f} B = {B:.4f} K/B = {K / B:.4f} n_crit = {n_crit:.4f}\n")
print(f"{'k':>3}{'f Nefes [Hz]':>15}{'(2k+1)/(2τ) [Hz]':>19}{'σ Nefes [1/s]':>16}{'σ analytic [1/s]':>19}")
for k, i in enumerate(np.argsort(eig_an.freqs)):
print(
f"{k:>3}{eig_an.freqs[i]:>15.4f}{(2 * k + 1) / (2.0 * TAU):>19.4f}"
f"{eig_an.growth_rates[i]:>16.4f}{sigma_analytic:>19.4f}"
)K = 3.3333 B = 3.0817 K/B = 1.0817 n_crit = 0.9245
k f Nefes [Hz] (2k+1)/(2τ) [Hz] σ Nefes [1/s] σ analytic [1/s]
0 166.6667 166.6667 26.0283 26.1674
1 500.0000 500.0000 26.0283 26.1674
2 833.3333 833.3333 26.0283 26.1674The frequencies land on the ladder to four decimal places and every rung shares one growth rate, as the closed form demands. The growth rates sit 0.14\ \mathrm{s^{-1}} below the analytic value, uniformly: that is the same O(\overline{M}) convective damping the passive cavity modes carried in section 2, and it is the only thing the zero-Mach algebra left out.
Sweeping \tau makes the ladder move. Every ITA frequency scales as 1/\tau and the growth rate as \ln(K/B)/\tau, while the cavity resonances, which know nothing of \tau, would not move at all.
taus = np.array([1.5e-3, 2.0e-3, 2.5e-3, 3.0e-3, 4.0e-3, 5.0e-3, 6.0e-3])
f_by_tau, g_by_tau = [], []
for tau in taus:
s = rijke(N, tau, ends="anechoic")
with warnings.catch_warnings():
warnings.simplefilter("ignore")
e = s.eigenmodes(freq_band=(20.0, F_MAX), growth_band=(-1200.0, 900.0), isentropic=True)
f_by_tau.append(np.sort(e.freqs))
g_by_tau.append(e.growth_rates.mean())
tau_fine = np.linspace(taus.min(), taus.max(), 400)
fig = make_subplots(
rows=2,
cols=1,
shared_xaxes=True,
vertical_spacing=0.08,
subplot_titles=("ITA frequencies ride the ladder (2k+1)/(2τ)", "one growth rate for every rung, σ = ln(K/B)/τ"),
)
for k in range(6):
fig.add_scatter(
x=tau_fine * 1e3,
y=(2 * k + 1) / (2.0 * tau_fine),
mode="lines",
showlegend=(k == 0),
name="analytic ladder",
line=dict(color="#9aa5b1", width=1.5, dash="dash"),
row=1,
col=1,
)
fig.add_scatter(
x=np.repeat(taus, [len(f) for f in f_by_tau]) * 1e3,
y=np.concatenate(f_by_tau),
mode="markers",
name="Nefes eigenmodes",
marker=dict(size=9, color=COLORWAY[0]),
row=1,
col=1,
)
fig.add_scatter(
x=tau_fine * 1e3,
y=np.log(K / B) / tau_fine,
mode="lines",
name="analytic σ",
line=dict(color="#9aa5b1", width=1.5, dash="dash"),
row=2,
col=1,
)
fig.add_scatter(
x=taus * 1e3,
y=g_by_tau,
mode="markers",
name="Nefes σ",
showlegend=False,
marker=dict(size=9, color=COLORWAY[0]),
row=2,
col=1,
)
fig.update_yaxes(title_text="frequency [Hz]", range=[0, F_MAX], row=1, col=1)
fig.update_yaxes(title_text="growth rate [1/s]", row=2, col=1)
fig.update_xaxes(title_text="flame time lag τ [ms]", row=2, col=1)
fig.update_layout(title="Anechoic tube: the ITA modes are the flame's clock", width=950, height=620)
fig.show()The number of ITA modes below 1\ \mathrm{kHz} is just the number of rungs that fit, so it grows with \tau. At \tau = 3\ \mathrm{ms} there are three, which is why the mode count in section 3 happened to double.
The anechoic experiment identifies the ITA family but changes the tube. On the physical tube, with its reflecting ends, we instead dial the flame gain n from its nominal value toward zero and watch where each eigenvalue goes.
The two families behave oppositely, and the closed forms of sections 2 and 4 say why. A cavity mode’s dispersion relation does not contain n, so as the flame weakens the mode relaxes onto its passive resonance with a bounded growth rate. An ITA mode’s growth rate is \ln(K/B)/\tau with K\propto n, so it dives toward -\infty and the mode leaves the spectrum.
eigenvalue_trajectory computes the full spectrum once with eigenmodes, then follows each mode step by step as n changes, using the eigenvector from the previous step. That continuation tracks individual modes through near-degeneracies, where comparing two independently computed spectra would mis-match them. The first call must return a complete certified spectrum; otherwise spurious modes may be followed by mistake.
We stop just short of n=0: a diverging growth rate has no meaningful endpoint, and it is the trend that classifies.
gains = np.linspace(1.0, 0.05, 73)
with warnings.catch_warnings():
warnings.simplefilter("ignore") # branches are expected to be lost as the ITA growth diverges
traj = eigenvalue_trajectory(
build=lambda n: rijke(n, TAU),
params=gains,
freq_band=(20.0, F_MAX),
growth_band=(-600.0, 600.0), # wide enough for a complete six-mode initial spectrum
isentropic=True,
param_name="FTF gain n",
)
traj| branch | f start → end [Hz] | growth start → end [1/s] | status |
|---|---|---|---|
| 0 | 162.6 → 176.6 | +155.5 → -929.7 | alive |
| 1 | 306.1 → 269.5 | -62.57 → -19.9 | alive |
| 2 | 464.2 → 483.4 | -355.3 → -1054 | alive |
| 3 | 578.1 → 597.2 | +282 → +34.61 | alive |
| 4 | 796.8 → 814.3 | +362.6 → +16.39 | alive |
| 5 | 921.9 → 861.1 | -292.3 → -958.7 | alive |
traj.plot(
show_markers=False,
width=950,
height=460,
title="Eigenvalue paths as the flame gain is dialed from n = 1 to n = 0.05",
).show()Three branches plunge; three do not. To turn that picture into labels we anchor on the passive spectrum, which is by definition the cavity family. A branch whose endpoint has relaxed onto a passive resonance is a cavity mode; one that has dived far away is ITA.
# Passive eigenvalues as complex omega, the anchor for the classification.
cavity_omega = 2.0 * np.pi * eig_passive.freqs - 1j * eig_passive.growth_rates
print(f"{'f @ n=1':>10}{'σ @ n=1':>11}{'f @ n=0.05':>13}{'σ @ n=0.05':>13} family")
for b in sorted(traj.branches, key=lambda br: br.freqs[0]):
# Distance from the branch endpoint to the nearest passive resonance, expressed in Hz.
gap = np.min(np.abs(cavity_omega - b.end)) / (2.0 * np.pi)
family = "cavity" if gap < 25.0 else "ITA"
print(f"{b.freqs[0]:>10.1f}{b.growth[0]:>+11.0f}{b.freqs[-1]:>13.1f}{b.growth[-1]:>+13.0f} {family}")
print(f"\ncavity resonances (flame off): {', '.join(f'{f:.0f}' for f in cavity_analytic)} Hz")
print(f"ITA ladder (2k+1)/(2τ) : {', '.join(f'{f:.0f}' for f in ladder)} Hz") f @ n=1 σ @ n=1 f @ n=0.05 σ @ n=0.05 family
162.6 +155 176.6 -930 ITA
306.1 -63 269.5 -20 cavity
464.2 -355 483.4 -1054 ITA
578.1 +282 597.2 +35 cavity
796.8 +363 814.3 +16 cavity
921.9 -292 861.1 -959 ITA
cavity resonances (flame off): 268, 599, 821 Hz
ITA ladder (2k+1)/(2τ) : 167, 500, 833 HzThe split is clean, and it agrees with the anechoic experiment: three ITA modes and three cavity modes.
Two things this exercise makes concrete.
The labels are reliable only when the initial spectrum was complete and the branches stay well separated along the sweep. Where two branches approach each other the families hybridize, the binary label is genuinely ill-defined, and the routine marks the branch as ambiguous rather than guessing.
The endpoint of a diving branch is meaningless. Push the sweep to exactly n=0 and the ITA growth rates run off to -\infty, and those branches are dropped. Read the trend, not the last point.
The same continuation works on any parameter you pass to the network builder, such as a duct length, an area, or a boundary reflection, not only the gain.
An active thermoacoustic network mixes two families of modes, and their dispersion relations say what separates them. The cavity relation contains the geometry and no flame parameter; the intrinsic relation contains the flame lag and no geometry.
eigenvalue_trajectory separates them on the physical tube: continued toward n=0, ITA branches dive in growth rate and cavity branches park on the passive resonances.The practical lesson is that an unstable mode’s family determines what influences it. Retuning the duct will not touch an ITA mode; retuning the flame’s time lag will not move a cavity resonance.