Engineering Acoustics/Harmonic Generation

Nonlinear Generation of Harmonics

As described in the entry for the qualitative description of shocks, finite amplitude waves in any medium will undergo a steepening phenomena cumulating in the formation of shock wave. For strong flows and wave conditions where fluid velocity is similar in magnitude to sound speed, u/c_o ≥ O(1), where u is the particle velocity and c_o is the ambient sound speed, the transition to a shockwave occurs rapidly and can be termed a local effect. For weaker wave conditions where u/c_o << 1, but nonlinear effects are still observable, wave steepening occurs over many wavelengths and can be termed a cumulative effect. In this regime of wave strengths an important result of wave deformation is the accumulation of harmonic content in the propagating waveform.

Progressive Wave Deformation

To describe the harmonic content of the deformed wave profile consider some analysis can be carried out for the case of a plane wave propagating in the x⁺ direction driven by a boundary piston with velocity u_o = sin(ωt), where ω is the driving frequency, and t is the time variable. In this case the resulting sound field depends only on the x⁺ wave, thus falls under the simple wave assumption and can be defined using a reduced equation for progressive waves in an inviscid fluid:

${\frac {\partial u}{\partial t}}+\left(c+u\right){\frac {\partial u}{\partial x}}=0$

${\frac {\partial u}{\partial t}}+\left(c_{o}+\beta u\right){\frac {\partial u}{\partial x}}=0$

$\beta =1+{\frac {1}{2}}{\frac {B}{A}}$

The relation between β and B/A is given to highlight the relation to the acoustic parameter of nonlinearity. The equations given describe a propagating planar wave, in which the propagation velocity for any particular point is given by local value of (c_o + βu), as opposed to simply c_o in the case of an assumed linear wave. For an initially sinusoidal wave profile, the wave apexes propagate with the greatest velocity, while equilibrium points propagate only at the ambient sound speed. This progression is qualitatively depicted in Figure 1, where the trajectory (wave velocity) of the wave maximum, minimum, and equilibrium points are plotted. As the trajectories of different points on the wave are non-parallel, the wave will deform as it propagates. The rate of deformation depends on the magnitude of difference between the various trajectories in the wave profile, and those depend on both the induced particle velocity and the fluids value of B/A. As a result of this dependence, a fluid with a higher B/A value will exhibit more rapid wave deformation than a fluid with a lower B/A value if the same boundary velocity is applied to both.

Figure 1: Progressive deformation of an initially sinusoidal wave profile. Local trajectories are plotted on the x-t axis for the wave maximum (c + v), the local equilibrium (c_o), and the wave minimum (c + v).

According to the progressive wave equation given for an inviscid fluid - and the process depicted in Figure 1 - the wave maximum will eventually catch up and overtake the wave front to form a discontinuous shock. In real fluids this is not necessarily the only possible outcome, as all sound waves in real fluids will attenuate as they propagate to some extent. For many dissipative processes, the effect on the wave is proportional to ω², thus the generated higher harmonics are dissipated more severely than the fundamental frequency. In this regard the effects of dissipation hamper wave steepening, and for some wave amplitudes a quasi-steady wave-form can be reached where the nonlinear steeping effects are perfectly balanced by dissipative effects.

Provided sufficient amplitude, or a nearly inviscid fluid, shock waves are formed in a progressive wave. It is for this reason that cumulative deformation is more likely to result in shocks in water than air, as air can be highly dissipative at the required wave amplitudes. Although not discussed in the sections to follow, the generation of new harmonic content continues after shock formation. The analytical description of this process is not fundamentally different from the analysis given, as the weak shock assumption is employed. For more details on the post shock regime refer to the seminal paper by Blackstock or in a variety of reference texts on the nonlinear acoustics.

Frequency Analysis of Solutions Obtained Using the Method of Characteristics

The most intuitive analytical solution to the deformed wave profile, u(x,t), is obtained when using an approach method based on the method of characteristics. When applied to planar progressive waves this approach is described in the entry for the qualitative description of shock waves. Discussion focusing on the intermediate wave profile in addition to the shock formation properties are given in many reference texts on nonlinear acoustics, including those by Hamilton and Blackstock, Elfno, Beyer, or Pierce.

As depicted to in Figure 1, the essence of this approach is to project know values or ordinary differential eqations along characteristic paths defined by dx/dx = (c + u). For the case of a single propagating wave driven from a sinusoidal boundary, u_o = sin(ωt), the solution is particularly simple as each location in the domain corresponds to only one characteristic path carrying constant values of u, c, etc. At any position throughout the domain the characteristic trajectories can be identified and the corresponding value of u is used to construct the solution. In the following example, the implicit equations described in Elfno were used to calculate the wave profile:

$u=u_{o}\sin \left(\tau +{\frac {\beta u}{c_{o}^{2}}}x\right),$

$\tau =t-{\frac {x}{c_{o}}}$

For this particular example the wave amplitude was set to achieve shock formation after five wavelengths. For generality all magnitudes of the solution are given in non-dimensional form. For the fluid, the parameter of nonlinearity was set to B/A = 5.0, which corresponds to fresh water at 20^oC. The upper portion of Figure 2 gives the spatial wave profile in which the wave steepening is readily apparent. In the lower three panels of Figure 2, the frequency spectrum of the velocity wave is shown as while passing through the regions indicated by red bands. The frequency spectrum plots were obtained using Discrete Fourier transform (DFT) of the calculated time domain signal at each indicated location.

Figure 2: Harmonic contnet of a deforming progressive wave. β = 3.5, corresponding to B/A = 5.0 for distilled water at 20^oC.

In the frequency spectra, where f_o = 2πω, the progressive wave contains only the fundamental frequency at the driven boundary, x = 0. Importantly, it can also be seen the propagated and deformed wave profile contains only integer harmonics of the initial frequency, and the magnitude of these harmonics increases with propagation distance.

Direct Analytical Solution to Harmonic Profile

While the method of characteristics provides spatial and time domain wave profiles, from which frequency content can be analyzed, a more direct approach is available to directly solve for the harmonic content of the deformed progressive wave. Again consider the example of a planar wave driven by a sinusoidal boundary condition, u_o = sin(ωt). One common approach for obtaining a general solution when system input is periodic, is to assume a periodic system response expressed as a Fourier Series. As shown by Pierce, applying this approach to the progressive wave problem gives the series solution as:

$u\left(\omega \tau ,\sigma \right)=\sum _{n=1}^{\infty }B_{n}\left(\sigma \right)\sin \left(n\omega \tau \right)$

$\sigma ={\frac {x}{x_{\text{shock}}}}$

$\tau =t-{\frac {x}{c_{o}}}$

$B_{n}\left(\sigma \right)={\frac {2}{\pi }}\int _{0}^{\pi }u\left(\omega \tau ,\sigma \right)\sin \left(n\omega \tau \right)d\left(\omega \tau \right)$

From this form of the solution the coefficient magnitudes directly yield the harmonic amplitudes, while the complete Fourier series yields the solution in the time and space domains. The challenge in this approach is the evaluation of the Fourier coefficients. In the context of nonlinear acoustics this solution is first attributed to Fubini in 1935, where manipulation of the integral terms was used to achieve the integral form of Bessel’s function, thus the harmonic components can be directly defined according to:

$B_{n}(\sigma )={\frac {2u_{o}}{n\sigma }}J_{n}(n\sigma )$

Where J_n is the Bessel function of the first kind. Historically, the subsequent work by Blackstock in 1966 clarified the physical implications of this solution somewhat and also brought this solution to a wider audience in nonlinear acoustics. For additional details on this derivation refer to the texts by Pierce, Elfno, or Hamilton. For completeness, Figure 3 plots the continuous harmonic profile as a function of normalized propagation distance. The wave conditions are the same as those used in calculations for Figure 2. In this particular example the shock formation occurs at five wavelengths; however, the plotted harmonic profile is general to any value of x_shock.

Figure 3: Continuous harmonic profile of a deforming progressive wave. β = 3.5, corresponding to B/A = 5.0 for distilled water at 20^oC.