Asymptotic behavior of a frequency-domain nonlinearity indicator for solutions to the generalized Burgers equation Open Access

Nonlinearity indicators have been used to characterize and understand the cumulative nonlinear behavior of propagating, high-amplitude acoustic waves. Examples of nonlinearity indicators include the pressure waveform derivative skewness,^1,2 average and wave steepening factors,^3–5 bicoherence,⁶ and the Morfey and Howell-derived Q/S.^1,7–10 Indicators have been most useful for analyzing noise, e.g., from high-speed jets, but real-world system complexities can hinder connections between analytical and noise¹¹ treatments. To guide the characterization of nonlinear propagation in noise, this letter explores the asymptotic behavior of a frequency-domain, single-point nonlinearity indicator for solutions to the generalized Burgers equation (GBE), building upon prior work which has shown that high-harmonic decay proceeds more slowly than linear absorption predicts.^12,13 Similar trends have been predicted for noise propagation,^14,15 and the letter concludes by discussing the connections between the analytical behavior and that seen for jet noise.¹⁶ 

An ensemble-averaged, frequency-domain version of the GBE was recently used to derive three quantities that yield the change in sound pressure level (SPL) with distance caused by geometric spreading, absorption, and nonlinearity, respectively.¹⁷ The nonlinearity indicator is a single-point indicator, meaning it can be calculated directly from a single waveform measurement. It includes the Morfey and Howell-derived Q/S, defined as $Q/S=Qpp2/(prmsSpp)$⁠, where $Qpp2$ is the quadspectral density between the pressure and pressure-squared waveforms, p_rms is the waveform's root-mean-square pressure, and S_pp is the autospectral density. The calculated quadspectrum of the pressure and squared pressure reveals phase coupling between two different frequencies, which occurs from sum and difference-frequency nonlinear harmonic generation in steepening waves.^6,18

The indicator is derived from the time-domain version of the GBE for an arbitrarily diverging pressure waveform, p(t), in thermoviscous media, which may be written as

where r is the distance from the source; m is 0, 0.5, or 1 for planar, cylindrical, and spherical waves, respectively; δ is the diffusivity of sound; $c0$ is the equilibrium sound speed; τ is retarded time; β is the coefficient of nonlinearity; and $ρ0$ is the equilibrium density of air. The terms on the right-hand side of Eq. (1) show that the change in pressure over distance is related to the separate effects of divergence, absorption, and nonlinearity. If L_n is defined as the SPL for the nth harmonic of a periodic signal with fundamental frequency $ω1$⁠, the time-domain GBE in Eq. (1) can then be cast¹⁷ into a frequency-domain and ensemble-averaged⁷ form to give the total rate of change in L_n with distance as

where $η≡10log10(e) ≈ 4.34$⁠, $αn$ is the linear attenuation coefficient for the nth harmonic, and $νS$⁠, $να$⁠, and $νN$ are the rate of change in L_n specifically due to spreading, absorption, and nonlinearity, respectively.

The nonlinearity indicator, $νN$⁠, has been previously studied for two solutions to the GBE:¹⁷ the Blackstock Bridging Function (BBF),¹⁹ which ignores spreading and absorption, and the Mendousse solution,²⁰ which neglects spreading. For plane waves in a lossless medium (BBF), nonlinear losses at the shock cause the waveform to decay in amplitude, and after shock formation $νN$ is asymptotically negative for all frequencies [see Eqs. (25) and (26) in Ref. 17]. However, shocks in a lossy medium eventually thicken due to absorption, after which the waveform decay is slower than linearly predicted. This ongoing transfer of energy upward in the spectrum must correspond to a positive asymptotic value for $νN$⁠. Here, asymptotic expressions for planar and diverging waves in a thermoviscous medium¹³ are used to derive asymptotic values of $νN$ for such waves. The indicator is numerically calculated for both cases and compared against the analytical expressions. Connections are also made between the asymptotic $νN$ values for initial sinusoids and those for noise.^13,18 Knowing the asymptotic value of the $νN$ indicator, a single measurement can be used to help determine the long-range nonlinear progression of a waveform.

Nonlinear waveform distortion of sinusoids is based on the normalized distance $σ=x/x¯$⁠, where $x¯$ is the planar shock formation distance. The Mendousse solution is an infinite-series solution to the GBE for plane waves, but with thermoviscous absorption (α ≠ 0, m = 0).²⁰ The solution depends on the linear absorption coefficient calculated from the frequency of the initial sinusoid, $α1$⁠, and the corresponding Gol'dberg number, $Γ=1/x¯α1$⁠. As Γ → ∞, nonlinearity dominates linear absorption and the Mendousse solution becomes equal to the BBF. For strong waves (⁠$Γ≫1$⁠) in the old-age region (⁠$α1x=σ/Γ≫1$⁠) as defined by Blackstock et al.,²¹ the term in the series specifying the absorption decay becomes $e−nα1x$ instead of $e−n2α1x$⁠.¹³ In this letter, these two different decay rates are referred to as linear and quadratic exponential decays, respectively, because of their dependence on harmonic number, n. For the Mendousse solution, the spatial rate of change in L_n is given by the combined effects of absorption and nonlinearity, or $να+νN$⁠. Working from Eq. (2), the asymptotic value of $νN$ for the Mendousse solution can be expressed as

where $C(ωn)$ is a certain function independent of x. From Eq. (3), $νN$ is non-negative in the old-age region for all harmonics (⁠$νN,x→∞=0$ for n = 1). The expression shows an addition of energy due to nonlinearity and a slower old-age decay, unlike the loss of energy which occurs for the BBF (negative $νN$⁠). Note that the asymptotic behavior of $νN$ is independent of Γ, as is the asymptotic waveform shape.²¹ Since the waveform amplitude decays asymptotically, the sum of the effects of absorption and nonlinearity, $να+νN$⁠, remains negative.

For plane waves in a lossless medium, the first pressure discontinuity occurs at σ = 1, the shock reaches a maximum amplitude at σ = π/2, and a sawtooth wave forms at σ ≈ 3.²¹ Here, normalized waveforms out to σ = 20 are shown in Fig. 1(a) with Γ = 30. This moderate Γ results in a waveform that is close to, but not exactly, a sawtooth at σ = 3, and then begins to unsteepen due to strong absorption of the high frequencies. The waveform at σ = 20 is approaching the old-age region.

The change in waveform shape shown in Fig. 1(a) can also be observed in the harmonic spectral amplitudes, $Sn$⁠, shown in Fig. 1(b). The amplitude of the fundamental is always decreasing (negative $νN$⁠) as energy is given to higher harmonics, whose spectral amplitudes generally increase between σ = 1 and σ = π/2 as the waveform steepens. However, depending on the value of Γ, the low-frequency harmonic amplitudes (e.g., for n = 2) can also experience nonlinear losses as they drive higher-harmonic generation [see Fig. 6(c) in Ref. 17]. Shock thickening is manifested by the increasingly steeper rolloff for σ > π/2. The spatial rate of change of the spectral levels is described by $νN$⁠. The $νN$ values for several harmonics were shown in Fig. 6 of Ref. 17 as a function of σ for Γ = 30, but only up to σ = 3. The value of $νN$ as a function of frequency is calculated—by taking Fourier transforms to calculate Q/S in Eq. (2)—from the Mendousse solution waveforms at two distances and is shown in Fig. 1 (c). By the chain rule, Eq. (2) can be multiplied by $x¯$ to give $∂Ln/∂σ$ rather than $∂Ln/∂x$⁠. The calculated indicator, $x¯νN$⁠, is an accurate prediction of the spatial rate of change in L_n; the difference between $x¯νN$ and the numerical derivatives of L_n, $ΔL/Δσ$⁠, is much less than 1%. At σ = 3, the indicator values follow the analytical trend only for harmonics $n ≳ 10$⁠. At σ = 20, the $νN$ indicator has converged to the asymptotic expression in Eq. (3) with an error of less than 2% for n = 3 to n = 10.

A complete analytical solution to the GBE with spherical spreading (α ≠ 0, m = 1) does not exist, and difficulties with solving the equation have been described in Refs. 22–24. However, an asymptotic solution (⁠$α1r ≫1$⁠) up through the first four harmonics shows a nonlinear decay of $r−ne−nα1x$ rather than a linear decay of $r−1e−n2α1x$⁠.¹³ The asymptotic value of $νN$ for spherical spreading and thermoviscous absorption can once again be determined from Eq. (2) to be

At least to n = 4, the $–(n – 1)/r$ term represents a larger effective geometric spreading decay with increasing n, and the $(n2 – n)/α1$ term represents a reduction in absorption from a quadratic exponential to a linear exponential decay. In the limit of large r, the absorption term dominates the spreading term and $νN$ is positive for n > 1. Once again, $να+νN$⁠, remains negative asymptotically.

Since there is no known analytical solution to the GBE with thermoviscous absorption and spherical spreading, a numerical solution is used.²⁵ To compare with the analysis from Sec. 3.2, a Gol'dberg number of 30 is desired. However, due to divergence there is much less nonlinear steepening for a spherical wave than for a plane wave of the same initial amplitude. For this reason, an effective spherical Gol'dberg number of Λ = 30 is used.²⁶ In addition, a spherically normalized distance is used to describe the waveform propagation, defined as $ζ≡r0/x¯ ln(r/r0)$⁠, where $r0$ is a known distance at which the waveform is sinusoidal and $x¯$ is the shock formation distance for plane waves. With these parameters, the waveform can be characterized similarly to the Mendousse solution waveform so long as $kr≫1$⁠. For this example, the fundamental frequency was 10 kHz, $r0$ was 1 m, the initial waveform amplitude was 415 Pa, the temperature was taken to be 20  °C, and the atmospheric pressure was 1 atm, giving Λ = 30 and Γ = 348.

The waveform shock formation is shown in Fig. 2(a), showing some trends similar to those of the Mendousse waveforms in Fig. 1(a). The shock appears to be steepest for ζ ≈ π/2, after which the waveform begins to unsteepen. However, unlike the Mendousse waveforms, at ζ = 3.75 the waveform has nearly reached sinusoidal shape again. The rapid unsteepening occurs because ζ is a function of ln(r) rather than r, meaning that the physical distance propagated to get the same scaled values of σ and ζ is much larger for the diverging case. Absorption therefore causes a large decrease in the higher harmonic amplitudes, contributing to a thicker shock. Compared to the Mendousse spectral amplitudes in Fig. 1(b), the $Sn$ in Fig. 2(b) show a much larger overall decrease in amplitude due to divergence from ζ ≈ π/2 to ζ = 3, as well as significantly more high-frequency absorption at ζ = 3 and ζ = 3.75.

Given the ln(r) dependence of ζ, an asymptotic analysis can be performed at a relatively small scaled distance. Similar to the scaling by $x¯$ for the Mendousse case in Fig. 1(c), by the chain rule Eq. (2) can be multiplied by $ε ≡ x¯ex¯ζ/r0=x¯r/r0$ to give $∂Ln/∂ζ$ rather than $∂Ln/∂r$ in Fig. 2(c). The calculated $νN$ indicator in Fig. 2(c) converges to the analytical value from Eq. (4)—labeled ζ → ∞—at a much smaller normalized distance than does the Mendousse solution in Fig. 1(c). The waveform is approaching the old-age region at ζ = 3.75, with $α1r ≈ 0.64$⁠. The numerical derivative due to nonlinearity, calculated from $ΔLn/Δζ – ενS – ενα$⁠, is also shown because the $νN$ calculation deviates slightly from the actual spatial change. The slight discrepancy between the two is at least partially due to the extremely low spectral amplitudes of the harmonics, as seen at ζ = 3.75 in Fig. 2(b). Note that the nonlinear numerical derivative is calculated assuming complete accuracy in ν_S and ν_α, which both have well-understood mechanisms.²¹ Finally, to emphasize the distinction between the planar and diverging cases, the Mendousse solution asymptotic indicator value, $ενN,σ→∞$⁠, is shown in Fig. 2(c) and differs substantially from all other curves.

Various values of ζ were tested for the asymptotic behavior, with ζ = 3.75 selected as the smallest distance for which the convergence was good. The error between the $νN$ curve, numerical derivative, and asymptotic value is less than 8% for harmonics 3–6. For n > 6, the round-off error in the spectrum is too large for the $νN$ calculation to be accurate, so larger harmonic numbers are not shown. A substantial difference is seen between the Mendousse and the spherical spreading asymptotic values for the indicator, and $νN$ clearly follows the expression in Eq. (4), which was derived for n $≤$ 4. The agreement to n = 6 suggests the validity of the expansion in Ref. 13 to harmonics beyond n = 4.

This work has extended the prior analysis of Reichman et al.,¹⁷ where a quadspectral, frequency-domain nonlinearity indicator, $νN$⁠, was used to quantitatively determine the nonlinear rate of change in sound pressure level over distance for acoustic waves. Here, $νN$ for planar and diverging waves in a thermoviscous medium was shown to be positive and increasing with frequency in the old-age region, signifying a slower decay than predicted by linear theory. The derived asymptotic expressions in Eqs. (3) and (4) are found to agree with $νN$ calculated from analytical and computational solutions to the GBE.

The analytical work for sinusoids and thermoviscous media guides understanding of broadband noise evolution in air. In Ref. 16, experimental jet noise was numerically propagated to beyond 3 km—for atmospheric absorption and dispersion—and the nonlinear gain (NG), i.e., the change in sound pressure level due to nonlinearity, was calculated. The spatial derivative of the NG, $∂NG/∂r$⁠, is essentially equivalent to $νN$ and was compared for various frequencies as a function of distance. Figures 2 and 3 of Ref. 16 show $∂NG/∂r$ decreasing asymptotically to zero in the peak-frequency region, and for higher frequencies $∂NG/∂r$ asymptotically converges to a successively larger value that is constant with distance. These trends seen for jet noise propagation are the same trends shown here for both planar and diverging waves in a lossy medium.

A quantitative connection can be found between the jet noise example and the sinusoidal results shown in Sec. 4. As r becomes large, the spreading term in Eq. (4) becomes much smaller than the absorption term. Neglecting the spreading term and taking the ratio of $νN,r→∞$ for two different harmonics, n_p and n_q, yields $(nq2−nq)/(np2−np)$⁠. If the jet noise peak frequency in Ref. 16 is taken to be 113 Hz, the ratio of the asymptotic $νN$ values for 8, 10, and 12.5 kHz for an initial sinusoid in a thermoviscous medium should be 1.62, 2.54, and 3.97, respectively. Using the asymptotic values for numerically propagated jet noise in air from Fig. 3 of Ref. 16, the actual ratios are 1.62, 2.50, and 3.71. The ratios from the data are very close to, but slightly less than the predicted ratios. The slight difference could be due to neglecting the spreading term in Eq. (4), or to the differences in absorption between air in Ref. 16 and the thermoviscous medium assumed in the theoretical analysis. Regardless, the close agreement suggests that $νN$ calculations help to understand the asymptotic behavior of noise waveforms as well.

The results in this letter demonstrate that a single-point waveform measurement can be used to determine whether a nonlinear wave—or noise waveform¹⁶—has progressed into the old-age region. The expressions in Eqs. (3) and (4) could be extended to find asymptotic expressions for other types of spreading (e.g., cylindrical), as well as calculating asymptotic trends in $νN$ for experimental noise that includes other atmospheric effects (e.g., turbulence).

This work was carried out in part under an Air Force Research Laboratory Small Business Innovation Research (SBIR) program with Blue Ridge Research and Consulting. Support from the Utah NASA Space Grant Consortium is also gratefully acknowledged.