3.1.1 Finite differences scheme

In this numerical scheme we approximate pressure and velocity derivatives using finite differences, obtaining a linear, homogeneous set of equations. We begin by writing (2.58) as two first-order ordinary differential equations (ODEs) for the pressure and velocity oscillations

(3.3) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{d}p_{1}}{\text{d}x}=-\frac{\text{i}\unicode[STIX]{x1D70C}_{m}}{F_{\unicode[STIX]{x1D708}}}u_{1}, & \displaystyle\end{eqnarray}$$

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{d}u_{1}}{\text{d}x}=-\text{i}B^{2}p_{1}-\frac{A}{F_{\unicode[STIX]{x1D708}}}u_{1}, & \displaystyle\end{eqnarray}$$

where $A$ and $B$ are defined in (2.59)–(2.60). In order to solve the system of equations, the tube is divided into $n$ equally spaced intervals, each interval $\unicode[STIX]{x0394}x=L/n$ long. The derivatives $\text{d}p_{1}/\text{d}x$ and $\text{d}u_{1}/\text{d}x$ are expressed using forward differences, such that (3.3) and (3.4) are translated into $2n$ equations. The remaining 2 equations arise from boundary conditions, demanding $u_{1,1}=u_{1,n+1}=0$ , namely that the velocity is zero at the closed ends. Ultimately, the following system is obtained

(3.5) $$\begin{eqnarray}\unicode[STIX]{x1D63C}\left(\begin{array}{@{}c@{}}p_{1,1}\\ \vdots \\ p_{1,n+1}\\ u_{1,1}\\ \vdots \\ u_{1,n+1}\end{array}\right)=0,\end{eqnarray}$$

where $\unicode[STIX]{x1D63C}$ is the coefficients matrix, and the vector comprises pressure and velocity at all the discrete points in the system. We seek the non-trivial solution by setting $\text{det}(\unicode[STIX]{x1D63C})=0$ . Since $\unicode[STIX]{x1D63C}\in \mathbb{C}$ , guesses for two variables in $\unicode[STIX]{x1D63C}$ are required to allow the use of a gradient-based method, targeting $\text{Re}[\text{det}(\unicode[STIX]{x1D63C})]=\text{Im}[\text{det}(\unicode[STIX]{x1D63C})]=0$ . Our guesses are the resonant frequency, $\unicode[STIX]{x1D714}$ , and the temperature difference on the stack, $\unicode[STIX]{x0394}T$ . The value of $n$ , arbitrarily chosen to increase calculation accuracy, determines the number of eigenvalues in $\unicode[STIX]{x1D63C}$ , representing different harmonic modes. As our model assumes a monochromatic wave, neglecting effects of higher harmonics, we are only interested in the smallest eigenvalue of $\unicode[STIX]{x1D63C}$ , corresponding to the resonant frequency of the first natural mode,  $\unicode[STIX]{x1D714}$ .

3.1.2 Zero acoustic power – standing-wave approximation

In this approach, we follow Arnott et al. (Reference Arnott, Belcher, Raspet and Bass1994) and make the argument that the total acoustic power within the system – produced and dissipated – sums up to zero under onset conditions. In this stability limit, any increase in $\text{d}T_{m}/\text{d}x$ will enhance the oscillations.

In a low aspect-ratio straight tube, oscillations only slightly deviate from a standing-wave form, where pressure and velocity are $90^{\circ }$ out of phase. If the Helmholtz number, $\unicode[STIX]{x1D6FA}$ , satisfies

(3.6) $$\begin{eqnarray}\unicode[STIX]{x1D6FA}=\frac{\unicode[STIX]{x1D714}L_{s}}{a}\ll 1,\end{eqnarray}$$

the stack is very short compared with the wavelength, allowing us to assume the wave form is only weakly distorted by the presence of the stack. Accordingly, we may use the standing-wave approximation (Swift Reference Swift1988) for the pressure and velocity perturbations

(3.7) $$\begin{eqnarray}\displaystyle & \displaystyle p_{1}=p_{A}\text{cos}(x/\unicode[STIX]{x1D706}), & \displaystyle\end{eqnarray}$$

(3.8) $$\begin{eqnarray}\displaystyle & \displaystyle u_{1}=\frac{\text{i}p_{A}}{\unicode[STIX]{x1D70C}_{m}a}\text{sin}(x/\unicode[STIX]{x1D706}), & \displaystyle\end{eqnarray}$$

where $p_{A}$ is the (dimensional) maximum pressure amplitude. We use (2.15) to scale these equations, approximating $\unicode[STIX]{x1D70C}_{m}$ and $a$ as their reference values ( $\unicode[STIX]{x1D70C}_{0}$ and $\unicode[STIX]{x1D706}\unicode[STIX]{x1D714}$ , respectively) for simplicity, to obtain

(3.9) $$\begin{eqnarray}\displaystyle & p_{1}=p_{A}\text{cos}(x), & \displaystyle\end{eqnarray}$$

(3.10) $$\begin{eqnarray}\displaystyle & u_{1}=\text{i}p_{A}\text{sin}(x), & \displaystyle\end{eqnarray}$$

where $p_{A}$ is the scaled maximum pressure amplitude.

The time-averaged acoustic power flux is

(3.11) $$\begin{eqnarray}{\dot{W}}=\frac{1}{2\unicode[STIX]{x03C0}}\int _{0}^{2\unicode[STIX]{x03C0}}p_{1}(x,t)u_{1}(x,t)\,\text{d}t=\frac{1}{2}\text{Re}[\widetilde{p}_{1}u_{1}],\end{eqnarray}$$

where the tilde denotes a complex conjugate. Substituting (3.9) and (3.10) into (3.11) yields the expected result ${\dot{W}}=0$ , reproducing a well-known attribute of pure standing waves producing zero acoustic power. Seeking a non-trivial result, we take the derivative of (3.11) to obtain

(3.12) $$\begin{eqnarray}\frac{\text{d}{\dot{W}}}{\text{d}x}=\frac{1}{2}\text{Re}\left[\widetilde{u}_{1}\frac{\text{d}p_{1}}{\text{d}x}+\widetilde{p}_{1}\frac{\text{d}u_{1}}{\text{d}x}\right],\end{eqnarray}$$

into which (3.3) and (3.4) are substituted and only then the standing-wave pressure and velocity distributions may be used. This way the imperfect phasing between pressure and velocity is accounted for in the differential relations of (3.3) and (3.4). Ultimately we obtain

(3.13) $$\begin{eqnarray}\frac{\text{d}{\dot{W}}}{\text{d}x}=\frac{p_{A}^{2}}{2}\left(\unicode[STIX]{x1D70C}_{m}\text{Im}\left[\frac{1}{F_{\unicode[STIX]{x1D708}}}\right]\sin ^{2}(x)+\text{Im}[B^{2}]\text{cos}^{2}(x)-\frac{\unicode[STIX]{x1D70C}_{m}}{2}\text{Im}\left[\frac{A}{F_{\unicode[STIX]{x1D708}}}\right]\sin (2x)\right),\end{eqnarray}$$

which may be integrated over all the segments in the system and equated to zero, resulting in an implicit equation for $\unicode[STIX]{x0394}T$ on the stack. The numerical solution to this equation is the desired temperature difference $\unicode[STIX]{x0394}T_{onset}$ . Since $\unicode[STIX]{x0394}T_{onset}$ is not large, the gas properties do not vary significantly with temperature. For simplicity, all quantities besides $C_{m}$ are taken as constants, representative of the working fluid mixture and independent of  $T_{m}$ .

4.1.1 Heat transfer mechanisms in thermoacoustic instability

Figure 2. Stability curves for conduction-driven (red) and phase-exchange (blue) thermoacoustic systems, in which the working fluid is air or an air–water mixture, respectively, shown as a function of the Womersley number, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . Solid and dashed curves are calculated using the approaches described in § 3.1. All calculations assume ideal sorption ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ ). The mean pressure is $p_{m}=1$ bar, the stack cold-side temperature is $T_{c}=23\,^{\circ }\text{C}$ , and scaled stack length and location are $L_{s}/L=0.026$ and $L_{0}/L=0.84$ , respectively. Experimental results, marked by filled dots, are taken from Tsuda & Ueda (Reference Tsuda and Ueda2017).

First, we compare the two mechanisms of heat transfer, namely conduction and phase exchange. In this comparison, shown in figure 2, we calculated $\unicode[STIX]{x0394}T_{onset}$ using the two methods presented in § 3.1 as a function of the Womersley number, here denoted by $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . The results of conduction-driven thermoacoustics are obtained by setting $C_{m}\equiv 0$ , reducing the equations to the known form of ‘classical’ thermoacoustics (Rott Reference Rott1969). The experimental results, reported by Tsuda & Ueda (Reference Tsuda and Ueda2017), show good agreement with the theoretical curves. In these experiments the working fluid comprises air and water vapour that evaporate/condense on the stack solid walls, matching the simple case of phase change derived by Raspet et al. (Reference Raspet, Slaton, Hickey and Hiller2002), which in our model corresponds to $\unicode[STIX]{x1D702}_{n}=1$ (see § 2.1). This agreement validates our model’s (as well as that of Raspet et al. Reference Raspet, Slaton, Hickey and Hiller2002) ability to accurately capture the physical mechanisms in both conduction-driven and phase-exchange thermoacoustics.

Compared with the conduction-driven system, the phase exchange substantially decreases the onset temperature difference. This distinction is best understood comparing the conduction-driven and phase-exchange heat-to-acoustic work terms (see (2.59)). In the limit of ideal sorption, using (3.2) to express $\text{d}C_{m}/\text{d}x$ as a function of $\text{d}T_{m}/\text{d}x$ and assuming $Pr\approx Sc$ , these terms differ only by the factor

(4.1) $$\begin{eqnarray}\unicode[STIX]{x1D713}\approx \frac{C_{m}}{1-C_{m}}\frac{l_{h}}{RT_{m}},\end{eqnarray}$$

which serves as a ‘booster’ for the temperature gradient, employing it to transfer heat through phase exchange of the reactive gas at the sorbent interface. The term $C_{m}/(1-C_{m})$ equals $1$ for $C_{m}=0.5$ and rapidly increases when the reactive gas constitutes more than half of the mixture. The non-dimensional group, $l_{h}/(RT_{m})$ , is typically much larger than $1$ (see table 1). If $\unicode[STIX]{x1D713}=1$ the two mechanisms contribute equally (approximately) to the conversion of heat to acoustic work, although choosing an appropriate reactive gas and setting the working conditions to increase $C_{m}$ , one may easily achieve $\unicode[STIX]{x1D713}\gg 1$ . As $\unicode[STIX]{x1D713}$ increases, more acoustic work is produced for the same temperature gradient or, in our case, onset is reached with a lower temperature difference along the stack.

Comparing between the two approaches taken for the stability analysis calculations, we notice the ‘zero acoustic power’ curves lie below their corresponding curves calculated using the numerical scheme. In the finite differences scheme $p_{1}(x)$ and $u_{1}(x)$ are the explicit solution to the boundary value problem in (3.3)–(3.4), while in the ‘zero acoustic power’ method these are assumed to obey a pure standing-wave behaviour, independent of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . This ‘ideal wave’ assumption leads to an artificial decrease in $\unicode[STIX]{x0394}T_{onset}$ since the method does not account for all the losses in the system. The better agreement with experimental results confirms the improved accuracy of the numerical scheme over the ‘zero acoustic power’ method.

Table 1. Characteristic values of parameters affecting the choice of reactive gas for phase-exchange thermoacoustic systems, shown for water, ethanol and methanol. Values are calculated for a mean pressure $p_{m}=1$ bar.

Another point worth mentioning is that, while the conduction-driven curves are unbounded, with $\unicode[STIX]{x0394}T_{onset}$ tending to infinity for $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\rightarrow 0$ and $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\rightarrow \infty$ , the temperature in phase-exchange systems is bounded by a ‘depletion temperature’, $T_{max}$ , namely the temperature at which the sorbent can no longer sorb gas (or, equivalently, the boiling temperature, $T_{b}$ , for a simple phase change process), such that $\unicode[STIX]{x0394}T_{onset}\leqslant T_{max}-T_{c}$ . The Clausius–Clapeyron equation, used in our model, relates the partial pressure (in our scaled notation, $C_{m}$ ) to the temperature for an equilibrium state of saturation. If the fluid temperature surpasses $T_{max}$ , the mean pressure must increase to allow the fluid to remain saturated. In our linear analysis the mean pressure is constant (see § 2.1); if we allow the gas mixture temperature to exceed $T_{max}$ the reactive gas will no longer be saturated and thus, it cannot condense into a liquid film. This limitation may, in principal, be corrected by extending the model to include $O(\unicode[STIX]{x1D700}^{2})$ terms, thus adding a second-order, $x$ -dependent correction to the mean pressure.

4.1.2 Effects of mixture composition

We now examine different gas mixture compositions and how their properties affect the temperature difference required for onset. The inert gas in the mixtures is either air or helium, two commonly used gases in thermoacoustic systems (Arnott, Bass & Raspet Reference Arnott, Bass and Raspet1992; Swift Reference Swift1992). The reactive gases selected for this comparison are water, ethanol and methanol, all used as condensable vapours that undergo evaporation/condensation. Accordingly, $l_{h}$ is the latent heat of evaporation; for water this value is obtained from steam tables, while for ethanol and methanol it is calculated according to Majer & Svoboda (Reference Majer and Svoboda1985). All other gas mixture properties are calculated using the kinetic theory of gases (Poling, Prausnitz & O’Connell Reference Poling, Prausnitz and O’Connell2001). Characteristic values of the parameters affecting the choice of the reactive gases are listed in table 1.

Figure 3. Stability curves for gas mixtures employing air (solid lines) and helium (dashed lines) as inert gases. The reactive gases are water (black), ethanol (red) and methanol (blue). The mean pressure is $p_{m}=1$ bar, scaled stack length and location are $L_{s}/L=0.04$ and $L_{0}/L=0.8$ , respectively, and ideal sorption is assumed ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ ). Curves are plotted against: (a) the Womersley number, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\equiv h(\unicode[STIX]{x1D714}/\unicode[STIX]{x1D708})^{1/2}$ , where the stack cold-side temperature is $T_{c}=20\,^{\circ }\text{C}$ , and (b) the averaged concentration of the reactive gas in the stack, $\overline{C_{m}}$ , for $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2.5$ .

Figure 3(a) presents stability curves of the selected gas mixtures as a function of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . Mixtures employing air as the inert gas outperform their matching helium-based mixtures in terms of $\unicode[STIX]{x0394}T_{onset}$ . A possible explanation is helium’s large heat capacity, compared with that of practically any gas, which requires a larger temperature difference for the same amount of heat to be transferred to the gas through conduction during an oscillation cycle, resulting in higher values of $\unicode[STIX]{x0394}T_{onset}$ . The closer a mixture is to the reactive component’s boiling temperature, the less dominant this effect is as the reactive gas becomes a larger fraction in the mixture, reflected through the concentration, $C_{m}$ . As $C_{m}$ increases, more of the gas in the mixture participates in the isothermal phase-exchange heat transfer mechanism, as clearly indicated by the parameter $\unicode[STIX]{x1D713}\propto C_{m}/(1-C_{m})$ (see table 1 for characteristic values). In these calculations, the cold side temperature, $T_{c}$ , is fixed at $20\,^{\circ }\text{C}$ for all mixtures, resulting in different values of $C_{m}$ , determined by the proximity of $T_{c}$ to the fluid freezing and boiling temperatures. The value of $T_{c}$ was chosen for practical reasons as it is close to the ambient temperature, making it easier to maintain. While water at $p_{m}=1$ bar coexists as liquid/vapour between $0<T<100\,^{\circ }\text{C}$ , EtOH and MeOH span a wider temperature range, and their boiling temperatures are significantly lower than that of $\text{H}_{2}\text{O}$ (see table 1). Accordingly, the values of $C_{m}$ in their mixtures for $T_{c}\leqslant T\leqslant T_{b}$ are substantially higher than for water (see table 1 for $C_{m}$ values at $T=T_{c}$ ), contributing to lower values of $\unicode[STIX]{x0394}T_{onset}$ .

The results shown in figure 3(b) compare stability curves of the considered gas mixtures as a function of the averaged concentration of the reactive component in the stack, $\overline{C_{m}}$ . The bar in $\overline{C_{m}}$ denotes axial averaging across the stack, and is an exception to its typical use for time averaging. All curves are calculated for $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2.5$ , a value representative of good performance for all mixtures (see figure 3 a). Here, different values of $\overline{C_{m}}$ are obtained by setting different temperatures as $T_{c}$ , such that the temperatures on the stack are in the range $T_{c}\leqslant T\leqslant T_{b}$ . As $T_{c}$ increases this range in narrowed, increasing $\overline{C_{m}}$ , thus decreasing $\unicode[STIX]{x0394}T_{onset}$ as the phase-exchange mechanism is more dominant. At a certain point the temperature range narrows to match $\unicode[STIX]{x0394}T_{onset}$ , setting the lowest possible temperature difference that would initiate onset, under the considered conditions. For all mixtures presented here the lowest temperature difference is $\unicode[STIX]{x0394}T_{min}\approx 5\,^{\circ }\text{C}$ , with $\overline{C_{m}}\approx 0.85$ .

At low values of $\overline{C_{m}}$ the curves may be separated into two clearly visible groups: air- and helium-based mixtures. This separation is a consequence of the results shown in figure 3(a), where mixtures employing helium as the inert gas exhibited higher values of $\unicode[STIX]{x0394}T_{onset}$ . As $\overline{C_{m}}$ increases, the reactive gas becomes a larger fraction in the mixture, gradually suppressing the differences between the two groups. At the highest concentration values for which onset is reached ( $\overline{C_{m}}\approx 0.85$ ), all curves appear to converge to the same value. The values of $\unicode[STIX]{x0394}T_{onset}$ in this case are so low that the temperature gradient imposed on the stack is lower than the adiabatic temperature gradient experienced by the gas as it compresses and expands. Accordingly, the conduction-driven mechanism cannot convert heat to acoustic power, leaving only the phase-exchange mechanism to produce the required acoustic power that triggers the instability. This heat-to-acoustic power conversion is governed by the parameter $\unicode[STIX]{x1D713}$ , dependent on $C_{m}$ , $T_{m}$ and the reactive gas latent heat, $l_{h}$ . Here $\overline{C_{m}}$ is practically equal for all mixtures, leaving the non-dimensional group $l_{h}/(RT_{m})$ to distinguish between the gas mixtures and their performance. As $C_{m}$ on the stack is close to $1$ , all temperatures are very close to the reactive component’s boiling temperature, $T_{b}$ , such that $\unicode[STIX]{x1D713}$ may be approximated as

(4.2) $$\begin{eqnarray}\unicode[STIX]{x1D713}\approx \frac{0.85}{1-0.85}\frac{l_{h}}{RT_{b}},\end{eqnarray}$$

where $\overline{C_{m}}\approx 0.85$ is the value for which $\unicode[STIX]{x0394}T_{min}$ is obtained. The non-dimensional group $l_{h}/(RT_{b})$ displays similar values for all three reactive gases chosen (see table 1). As a result, all mixtures converge to a very similar value of $\unicode[STIX]{x0394}T_{min}$ , although slight changes in $\unicode[STIX]{x0394}T_{min}$ between the different mixtures are perceptible. These changes may be explained by the factor $1/T_{m}$ multiplying $\unicode[STIX]{x1D713}$ in the phase-exchange heat-to-acoustic work term, which may also be approximated as $1/T_{m}\approx 1/T_{b}$ . The small variations in the values of $l_{h}/(RT_{b}^{2})$ for the three gases, listed in table 1, are reflected in the values of $\unicode[STIX]{x0394}T_{min}$ seen in figure 3(b).

4.1.3 Sorption kinetics

Figure 4. Stability curves for an air–water gas mixture at various values of the partition coefficient $K$ , as a function of the kinetics non-dimensional number, $R_{k}$ . The mean pressure is $p_{m}=1$ bar, the Womersley number is $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2.5$ , and scaled stack length and location are $L_{s}/L=0.03$ and $L_{0}/L=0.8$ , respectively.

The curves presented in figure 4 examine effects of reaction kinetics in a sorption process, where $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}},\unicode[STIX]{x1D702}_{D}\neq 1$ . These stability curves, calculated for different values of the partition coefficient $K$ , are plotted against $R_{k}\equiv \unicode[STIX]{x1D708}/(h^{2}k_{d})$ , a non-dimensional number linking the characteristic viscous and desorption time scales such that $\unicode[STIX]{x1D70F}_{k}\equiv \sqrt{R_{k}}\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . All curves are calculated for an air–water mixture. $R_{k}\rightarrow 0$ and $R_{k}\rightarrow \infty$ denote the limits of fast and slow kinetics relative to the oscillation time scale, respectively. The dashed lines represent the two extreme cases for the sorbent/sorbate pair: ‘ideal’ sorption ( $K\rightarrow \infty$ ), in which the sorbent has no limitation on its capacity for sorbing gas (also applicable for evaporation/condensation from/to a liquid film, unbounded in its capacity for condensed vapour – see § 2.2), and no sorption ( $K\rightarrow 0$ ), where the sorbent has zero capacity for sorbing the reactive gas such that (2.58) is reduced to the conduction-driven form (Rott Reference Rott1969). All the curves approach the limit of $K\rightarrow 0$ for vanishingly slow kinetics (very large values of $R_{k}$ ), where the oscillations are too fast to allow mass exchange at the sorbent interface. In the case of fast kinetics, sorption occurs instantly, resulting in a quasi-steady state of equilibrium between the sorbent and gas phases. Accordingly, at low values of $R_{k}$ all curves flatten towards a specific value, determined by the partition coefficient, $K$ , specifying the sorbent capacity for sorbing the reactive gas. If $K$ surpasses a certain value, here $K\approx 10^{-2}$ , all curves are monotonically increasing, with the ascent towards the horizontal asymptote of $K\rightarrow 0$ occurring at different values of $R_{k}$ . A curious limit is obtained for the case of ideal sorption occurring at an infinitesimally slow rate, for which the kinetics parameter does not recover $\unicode[STIX]{x1D702}_{n}\rightarrow 1$ nor $\unicode[STIX]{x1D702}_{n}\rightarrow \infty$ , but rather an intermediate result, unaffected by the reaction rate and strongly dependent on the parameters $\unicode[STIX]{x1D6EC}$ and $N^{s}$ , characterising the sorbent/sorbate pair.

While large values of $K$ result in monotonically increasing curves, small values of $K$ clearly show a minimum at $R_{k}\approx 10^{-1}$ . Intuitively, one may assume the lowest temperature difference is always recovered for the limit of fast kinetics. For small $K$ values, however, the sorbent has very limited capacity for sorption, such that for $R_{k}\rightarrow 0$ mass exchange occurs in only a very small fraction of an oscillation cycle. During the remainder of the cycle, then, the sorbent is saturated and the phase-exchange mechanism is not active, resulting in an increase in $\unicode[STIX]{x0394}T_{onset}$ as $K$ decreases. An intentional, moderate increase in $R_{k}$ , namely a regulated delay in desorption rate, allows the sorbent to exchange mass with the gas mixture during a larger fraction of the oscillation cycle, thus lengthening the characteristic time in which the phase-exchange mechanism is active and consequently decreasing $\unicode[STIX]{x0394}T_{onset}$ .

4.2.1 Acoustic field characteristics

In the considered system configuration (shown in figure 1 a) the no-penetration boundary conditions at the two closed ends impose a near-standing-wave behaviour of pressure and velocity. Viscous dissipation at the solid boundaries, the presence of a stack and its asymmetric location within the tube all contribute to a slight distortion of the wave form, deviating from a standing wave, and thus allowing the production (or dissipation) of non-zero acoustic power. However, this deviation is normally very small in standing-wave systems (such as the one considered here, see figure 1 a), and the acoustic power flux

(4.3) $$\begin{eqnarray}{\dot{W}}={\textstyle \frac{1}{2}}\text{Re}[\widetilde{p_{1}}u_{1}]={\textstyle \frac{1}{2}}|p_{1}|\,|u_{1}|\cos (\unicode[STIX]{x1D719})\end{eqnarray}$$

is largely dependent on the phase angle between $p_{1}$ and $u_{1}$ , $\unicode[STIX]{x1D719}$ , which is typically in the range $|\unicode[STIX]{x1D719}|=90^{\circ }\pm 5^{\circ }$ (Swift Reference Swift1992).

Figure 5. Acoustic field characteristics for conduction-driven (red) and phase-exchange (blue) systems, compared to a pure standing wave (black, dotted): (a) real (solid line) and imaginary (dashed line) components of the pressure amplitude, scaled by the maximum pressure in the system, (b) real (solid line) and imaginary (dashed line) components of the velocity amplitude, scaled by the maximum velocity in the system, (c) phase angle between pressure and velocity, $\unicode[STIX]{x1D719}$ , throughout the system, and in the stack in particular (inset), for different mean concentrations of the reactive gas, $\overline{C_{m}}$ , and (d) acoustic power. The highest $\overline{C_{m}}$ value presented in (c) corresponds to the phase-exchange curves in (a,b,d). All properties are drawn against $\widehat{x}=x/L$ , the axial coordinate scaled by the system length. The shaded area marks the stack location in the system. The mean pressure is $p_{m}=1$ bar, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2$ , the scaled stack length and location are $L_{s}/L=0.04$ and $L_{0}/L=0.83$ , respectively, the heat flux input is $Q_{h}=25~\text{kW}~\text{m}^{-2}$ , and ideal sorption is assumed ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ ).

Figures 5(a) and 5(b) present the axial distribution of pressure and velocity, respectively, in conduction-driven (red) and phase-exchange (blue) systems. Each curve is scaled by its respective maximum value, allowing a clear, visual comparison of the wave forms. For reference, we added curves of a pure standing wave (black, dotted), where $p_{1}$ and $u_{1}$ simply obey cosine and sine functions, respectively (see § 3.1.2). In both figures, the conduction-driven curves only slightly deviate from the standing-wave curves; the largest deviation is unsurprisingly spotted in the stack, where acoustic oscillations are enhanced. The phase-exchange curves, calculated for a large value of the ‘reactive’ gas concentration ( $\overline{C_{m}}=0.91$ ) to demonstrate the magnitude of the effect, deviate strongly from the standing-wave curves. While in conduction-driven systems the stack may be thought of as an acoustic power point source, amplifying random displacements of gas parcels, in phase-exchange systems it also acts as a mass point source. The periodic mass exchange with the sorbent layer draws/releases mass from/to the gas mixture, consequently increasing the reactive gas mass flux in the stack, which translates directly to a sharp increase in velocity (both real and imaginary), seen clearly in figure 5(b) at $\widehat{x}\approx 0.83$ . The asymmetric location of the stack within the system shifts the imaginary component, quasi-sinusoidal velocity distribution towards the stack, where the velocity increases. As a result, the pressure distribution is also shifted from its quasi-cosine profile, and the absolute value of the pressure at the two closed ends clearly shows $|p_{1}(\widehat{x}=0)|>|p_{1}(\widehat{x}=1)|$ , corresponding to the increase in velocity near $\widehat{x}=1$ .

Figure 5(c) shows the phase angle between $p_{1}$ and $u_{1}$ along the system, for several values of $\overline{C_{m}}$ , the mean concentration of the reactive gas in the stack. As in figures 5(a) and 5(b), a pure standing-wave curve (black, dotted) and conduction-driven (red, dashed) curves are drawn for comparison. The inset displays a close-up view of the stack, where the phase angle changes rapidly. As reflected from the results of figure 5(a,b), the conduction-driven curve only slightly deviates from the standing-wave curve. In fact, the largest deviation, neglecting the sharp decrease near $\widehat{x}=0.5$ (where the pressure moves from leading the velocity to trailing it), is only $2.4^{\circ }$ . At the lowest concentration considered ( $\overline{C_{m}}=0.08$ ), the phase-exchange case is nearly identical to the conduction-driven curve. At such a low concentration, altering the working gas contributes merely to a decrease in the temperature difference along the stack (discussed in § 4.2.2), while the wave form is left unaffected. As $\overline{C_{m}}$ is increased, the phase angle distribution departs from the standing-wave behaviour, until at $\overline{C_{m}}=0.91$ a remarkable $35^{\circ }$ deviation from a pure standing wave ( $\unicode[STIX]{x1D719}=-125^{\circ }$ ) is obtained at the stack hot end. This difference in phase angles between the conduction-driven and phase-exchange mechanisms translates into an increase in acoustic power, as shown in figure 5(d) (calculated using (4.3)). Although the absolute pressure amplitudes decrease as $\overline{C_{m}}$ increases (not visible in figure 5 a due to scaling; discussed in § 4.2.2), the dramatic increase in phase angle leads to a higher acoustic power produced by the phase-exchange system.

4.2.2 Temperature difference and pressure amplitude

Figure 6. Limit cycle scaled temperature difference across the stack (a) and maximum pressure amplitude (b) as a function of the Womersley number, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , for different scaled cold side temperatures, $\widehat{T_{c}}=(T_{c}-T_{min})/(T_{max}-T_{min})$ . The mean pressure is $p_{m}=1$ bar, the scaled stack length and location are $L_{s}/L=0.04$ and $L_{0}/L=0.83$ , respectively, the heat flux input is $Q_{h}=25~\text{kW}~\text{m}^{-2}$ and ideal sorption is assumed ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ ).

Here, we calculate the temperature difference across the stack and the maximum pressure amplitude ( $\unicode[STIX]{x0394}T$ and $|p_{1}|$ , respectively), so as to establish a more complete description of the system limit cycle. Figure 6(a,b) presents these two quantities as a function of the Womersley number, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . For each calculation, the boundary condition for the stack cold side temperature, $T_{c}$ , is varied. For the sake of even comparison, both the cold-side temperature and the temperature difference on the stack are scaled as follows

(4.4a,b ) $$\begin{eqnarray}\widehat{T_{c}}=\frac{T_{c}-T_{min}}{T_{max}-T_{min}},\quad \widehat{\unicode[STIX]{x0394}T}=\frac{\unicode[STIX]{x0394}T}{T_{max}-T_{c}},\end{eqnarray}$$

where $T_{max}$ is the reactive fluid ‘depletion temperature’ (see § 4.1.1), and $T_{min}$ is its triple-point temperature, below which the reactive gas component is completely in the sorbed phase and $C_{m}=0$ . For the case of simple phase change through evaporation/condensation $T_{min}$ is simply the freezing temperature. The shape of the curves in figure 6(a) resembles that of the stability curves in figures 2 and 3(a): each curve has a unique $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ value, for which the limit cycle temperature difference is minimised. Here $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\approx 2$ marks the minimum in all curves, a representative result for a standing-wave system (Swift Reference Swift2002). When $\unicode[STIX]{x1D70F}_{D}\equiv \unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}Sc^{1/2}$ , the ratio of the characteristic time scales for molecular diffusion and acoustic oscillations is either too small or too large, $\unicode[STIX]{x0394}T$ increases, so as to compensate for the poor match between these time scales.

The asymmetric increase in $\widehat{\unicode[STIX]{x0394}T}$ for $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}<2$ and $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}>2$ is a result of the match between $\unicode[STIX]{x1D70F}_{D}$ and the system travelling-wave component, manifested by the deviation of the wave form from that of a pure standing wave,

(4.5) $$\begin{eqnarray}\unicode[STIX]{x1D703}=||\unicode[STIX]{x1D719}|-90^{\circ }|.\end{eqnarray}$$

In a pure travelling-wave system, perfect thermal contact of the working gas with the boundary is required, i.e. $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}},\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D6FC}},\unicode[STIX]{x1D70F}_{D}\ll 1$ (Swift Reference Swift2002). For a given $\widehat{T_{c}}$ , as $\widehat{\unicode[STIX]{x0394}T}$ increases so does $\overline{C_{m}}$ , subsequently increasing $\unicode[STIX]{x1D703}$ . For $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}<2$ , the increased $\unicode[STIX]{x1D703}$ favours the decrease in $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ and the system sustains a limit cycle until $\widehat{\unicode[STIX]{x0394}T}$ nearly approaches $1$ , marking its highest possible temperature difference. For $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}>2$ each curve reaches a different $\widehat{\unicode[STIX]{x0394}T}$ before a limit cycle is no longer sustained at the given heat input. Each $\widehat{T_{c}}$ value corresponds to a characteristic $\overline{C_{m}}$ value, to which a specific $\unicode[STIX]{x1D703}$ matches. At $\widehat{T_{c}}=0.3{-}0.7$ , $\unicode[STIX]{x1D703}$ is relatively small and, while the curves depart from their optimal standing-wave point ( $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\approx 2$ ), the moderate increase in $\unicode[STIX]{x1D703}$ does not favour the increase in $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , and a limit cycle cannot be sustained before $\widehat{\unicode[STIX]{x0394}T}$ approaches 1. At $\widehat{T_{c}}>0.85$ , $\unicode[STIX]{x1D703}$ is initially large and the curves do not deviate substantially from $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2$ . As $\widehat{\unicode[STIX]{x0394}T}$ increases, the increase in $\unicode[STIX]{x1D703}$ dramatically decreases the resonant frequency (discussed in § 4.2.3), thus reducing the system losses – proportional to $\unicode[STIX]{x1D714}$ – and enabling the curves to approach $\widehat{\unicode[STIX]{x0394}T}=1$ (the $\widehat{T_{c}}=0.9$ in particular).

The maximum pressure amplitude in the system is shown in figure 6(b) as a function of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , for different values of $\widehat{T_{c}}$ . As the mean pressure is kept constant, $T_{max}$ is also a constant and as $\widehat{T_{c}}$ is increased, the system thermal potential decreases. As a result, the system can produce lower pressure amplitudes such that $|p_{1}|$ systematically decreases as $\widehat{T_{c}}$ is increased. For $\widehat{T_{c}}=0.3{-}0.7$ , the curves exhibit a local minimum that vanishes as $\widehat{T_{c}}$ is further increased. The $\overline{C_{m}}$ values near the minima in these curves are small and the system resembles a conduction-driven system, in which the curve for $|p_{1}|$ tends to infinity for $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\rightarrow 0,\infty$ and displays a minimum at $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\approx 2$ . However, this resemblance disappears as $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ is shifted from the minimum and the curves point downwards near their critical $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , where $\widehat{\unicode[STIX]{x0394}T}$ reaches its highest values. The increase in $\widehat{\unicode[STIX]{x0394}T}$ increases $\unicode[STIX]{x1D703}$ , and as the wave form deviates from a standing wave, the produced pressure amplitudes decrease. The two maximum points linking a curve’s ends with its minimum are noticeably different. At the left point ( $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\approx 1$ ), the increase in $\unicode[STIX]{x1D703}$ matches the decrease in $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ towards values favouring travelling-wave phasing, locally increasing $|p_{1}|$ . The right point marks the global maximum in each curve, and is the result of an increase in thermal potential, while $\unicode[STIX]{x1D703}$ only moderately increases, retaining the resemblance to a standing-wave system and thus increasing $|p_{1}|$ . As $\unicode[STIX]{x1D703}$ is further increased, $|p_{1}|$ again decreases until a limit cycle can no longer be sustained. For $\widehat{T_{c}}=0.85{-}0.9$ , any resemblance to a conduction-driven system disappears: the local minimum vanishes, leaving only a maximum point from which the curves decrease in either direction.

4.2.3 Frequency response

Figure 7. (a) Resonant frequency (blue, left axis) and the deviation from the frequency of an ‘ideal’ standing wave (red, right axis) as a function of the averaged reactive gas concentration in the stack, $\overline{C_{m}}$ , for $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2$ at mean pressures 1 and 10 bar. (b) Resonant frequency at $p_{m}=1$ bar as a function of the Womersley number, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , for different scaled stack cold side temperatures, $\widehat{T_{c}}=(T_{c}-T_{min})/(T_{max}-T_{min})$ . Curves are calculated for scaled stack length and location $L_{s}/L=0.04$ and $L_{0}/L=0.83$ , respectively, heat flux input $Q_{h}=25~\text{kW}~\text{m}^{-2}$ and for ideal sorption ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ ).

In general, the resonant frequency of a thermoacoustic system is determined by its geometry and the sound velocity through its working fluid. Most thermoacoustic systems employ a gaseous working fluid, in which the sound velocity is $a=\sqrt{\unicode[STIX]{x1D6FE}R_{g}T_{m}/M_{mix}}$ , such that the resonant frequency approximately scales as $\unicode[STIX]{x1D714}\propto M_{mix}^{-1/2}$ . The resonant frequency of a phase-exchange system is calculated and plotted in figure 7(a) (blue) as a function of $\overline{C_{m}}$ , for two representative mean pressures. At $\overline{C_{m}}\lesssim 0.65$ , the two curves seem to follow the rationale described above: water vapour is lighter than air and thus increasing its fraction in the gas mixture increases the resonant frequency. However, at $\overline{C_{m}}\approx 0.8$ a maximum appears, after which an increase in water concentration leads to a decrease in frequency until eventually, at very large $\overline{C_{m}}$ values, the system can no longer sustain a limit cycle, at the given heat input.

The resonant frequency at the limit cycle inherently deviates from that of a pure standing wave since the latter cannot be practically sustained. This deviation, $\unicode[STIX]{x0394}f_{sw}=f-f_{sw}$ , is drawn in the right axis of figure 7(a) (red), and assists in understanding the competing phenomena leading to this unusual frequency response. In a diluted gas mixture ( $\overline{C_{m}}\ll 1$ ), the wave form resembles a standing wave, reflected through small values of $\unicode[STIX]{x0394}f_{sw}$ . Since a pure travelling-wave behaviour would substantially decrease the frequency, the minor travelling-wave component here slightly decreases it. As $\overline{C_{m}}$ increases, the travelling-wave component increases and the frequency deviates more and more from that characterising a standing wave. Towards higher vapour concentrations, $\unicode[STIX]{x0394}f_{sw}$ decreases sharply until, near $\overline{C_{m}}\approx 0.8$ , the increase in water vapour in the mixture, which tends to increase the resonant frequency, is balanced by a decrease emerging from the wave-form distortion. At the highest attainable concentrations, the wave is greatly distorted (see figure 5), leading to a considerable decrease in frequency despite the water-enriched mixture tending to increase it.

The same trend is obtained at higher mean pressures. As a representative example, curves calculated at a mean pressure $p_{m}=10$ bar (dashed) are drawn, displaying a similar behaviour but at slightly higher frequencies. Increasing the mean pressure inherently increases the depletion (or boiling) temperature $T_{max}$ ; for a given $\overline{C_{m}}$ , such increase in $p_{m}$ directly translates to an increase in $T_{m}$ (see (3.1)), subsequently increasing the sound velocity and thus increasing the resonant frequency. The largest $\overline{C_{m}}$ for which a limit cycle is attained moderately decreases with increasing $p_{m}$ . These terminal $\overline{C_{m}}$ values correspond to $\unicode[STIX]{x0394}T_{min}\approx 4.5,25\,^{\circ }\text{C}$ for $p_{m}=1,10$ bar, respectively. A good indicator for estimating these temperature differences is the parameter $\unicode[STIX]{x1D713}=(C_{m}/(1-C_{m}))(l_{h}/RT_{m})$ , representing the ratio of produced acoustic power by the conduction-driven and phase-exchange heat transfer mechanisms (see § 4.1.1). Since $\unicode[STIX]{x1D713}$ is a monotonically decreasing function of $p_{m}$ (see appendix A), a larger $\unicode[STIX]{x0394}T$ is required as $p_{m}$ increases, translating to a lower  $\overline{C_{m}}$ .

Studying the resonant frequency response to core variables affecting the system performance, we calculate $f$ as a function of the Womersley number, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , for different scaled cold-side temperatures, $\widehat{T_{c}}$ , and present the curves in figure 7(b). Increasing $\widehat{T_{c}}$ imposes a larger initial value of $\overline{C_{m}}$ on a system. As reflected in the results of figure 7(a), increasing $\widehat{T_{c}}$ in figure 7(b) initially increases $f$ until, for $\widehat{T_{c}}=0.9$ (corresponding to $\overline{C_{m}}\sim 0.9$ ), the resonant frequency decreases. The $\widehat{T_{c}}=0.3$ and $\widehat{T_{c}}=0.5$ curves display a similar trend of initial decrease and eventual flattening. The largest $\overline{C_{m}}$ in these curves is obtained at $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\approx 0.9$ (where $\widehat{\unicode[STIX]{x0394}T}$ approaches $1$ in figure 6 a); as $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ increases, $\unicode[STIX]{x0394}T$ decreases towards its minimum at $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\approx 2$ , consequently decreasing $\overline{C_{m}}$ . As a result, the resonant frequency decreases as the water vapour fraction in the mixture is reduced. Further increasing $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ shifts the system from its optimal standing-wave working point, and $\unicode[STIX]{x0394}T$ increases accordingly. This increases $\overline{C_{m}}$ , acting to further increase the frequency by increasing the mixture water vapour fraction while tending to decrease it by increasing the deviation from a standing wave. These two factors are balanced and the frequency remains approximately constant until a limit cycle can no longer be sustained at the given heat input. For $\widehat{T_{c}}>0.5$ the $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ range for which a limit cycle may be sustained narrows, and $\widehat{\unicode[STIX]{x0394}T}$ is confined within a smaller range. Consequently, the curves’ initial slope diminishes as the decrease in $\overline{C_{m}}$ subsides. For $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}>2$ , the increase in the travelling-wave components is augmented as $\widehat{T_{c}}$ increases, and tends to decrease $f$ . While for $\widehat{T_{c}}=0.7$ the increase in $\overline{C_{m}}$ is significant, resulting in a moderate decrease in $f$ (resembling the $\widehat{T_{c}}<0.7$ curves), the travelling-wave component growth at $\widehat{T_{c}}>0.7$ dominates the frequency response, leading to a sharp decrease in $f$ as $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ increases.

4.2.4 Varying heat flux input

Figure 8. (a) Scaled temperature difference across the stack, defined in (4.4), and (b) the heat flux input, scaled against the maximum pressure amplitude squared, as a function of the heat flux input, $Q_{h}$ , for different reactive gas concentrations. While in (b) only the limit cycle results (solid lines) are displayed, the dashed lines in (a) represent solutions to a one-dimensional conduction equation, obtained for $Q_{h}<Q_{crit}$ . The mean pressure is $p_{m}=1$ bar, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=1.4$ , the scaled stack length and location are $L_{s}/L=0.04$ and $L_{0}/L=0.83$ , respectively, and ideal sorption is assumed ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ ).

Part of the heat input to a thermoacoustic system is converted into acoustic power, sustaining the oscillations, and the rest is ejected as waste heat. For a given system with constant mean pressure, increasing the heat input to the system may increase the stack hot-side temperature, and/or the pressure and velocity amplitudes. The curves in figure 8(a) display $\widehat{\unicode[STIX]{x0394}T}$ as a function of the heat flux input to the system, $Q_{h}$ , for different $\overline{C_{m}}$ values. In all curves $Q_{h}=0$ trivially results in $\widehat{\unicode[STIX]{x0394}T}=0$ ; as $Q_{h}$ is increased $\widehat{\unicode[STIX]{x0394}T}$ initially increases (dashed), until the heat flux input reaches a critical value, $Q_{crit}$ , from which a limit cycle can be sustained. When $Q_{crit}$ is surpassed, $\widehat{\unicode[STIX]{x0394}T}$ moderately decreases before flattening towards a constant, representative value for each $\overline{C_{m}}$ . The weak dependency of $\widehat{\unicode[STIX]{x0394}T}$ on $Q_{h}$ is somewhat unintuitive, and may be explained by considering (2.65), namely the total power flux equation. The equation may easily be manipulated into the form

(4.6) $$\begin{eqnarray}\frac{\text{d}T_{m}}{\text{d}x}=\frac{Q_{h}-f_{1}(p_{1},u_{1},T_{m})}{f_{2}(u_{1},T_{m})+\text{diffusion terms}},\end{eqnarray}$$

where $f_{1}$ and $f_{2}$ are complex functions of the pressure, velocity and temperature (given explicitly in (2.65)), and the diffusion terms are a function of the temperature alone. Using (3.3) we may consider $f_{1}(p_{1},\text{d}p_{1}/\text{d}x,T_{m})$ and $f_{2}(\text{d}p_{1}/\text{d}x,T_{m})$ , explicitly dependent on $p_{1}$ and its derivative, and only implicitly dependent on $T_{m}$ through the gas mixture properties. Since this implicit relation is relatively weak, we assume $f_{1},f_{2}\neq f(T_{m})$ and turn to focus on the pressure dependency.

The pressure distributions in figure 5(a) display the two extremes of a phase-exchange system: conduction-driven system, where the phase-exchange mechanism is absent, and nearly isothermal system, in which only the phase-exchange mechanism converts heat to acoustic power. As opposed to velocity (figure 5 b), the pressure distribution is only weakly affected by the change in heat transfer mechanism. Its wave form can therefore be described, to a good degree of accuracy, by considering a pure standing wave, as given in (3.9). For $Q_{h}<Q_{crit}$ , a limit cycle cannot be sustained and (4.6) reduces to a one-dimensional heat equation in a binary gas mixture, in which $p_{1}=0$ . Increasing $Q_{h}$ in this region trivially increases $\widehat{\unicode[STIX]{x0394}T}$ . At $Q_{h}=Q_{crit}$ a jump in $p_{1}$ (followed by an increase as $Q_{h}$ increases) is obtained, leading to a local decrease in $\widehat{\unicode[STIX]{x0394}T}$ , indicated clearly in (4.6). As $Q_{h}$ is further increased, $p_{1}$ increases, swiftly overweighing the effect of the diffusion terms on $\widehat{\unicode[STIX]{x0394}T}$ , and the temperature gradient simplifies to

(4.7) $$\begin{eqnarray}\frac{\text{d}T_{m}}{\text{d}x}\approx G_{1}\frac{Q_{h}}{p_{A}^{2}}-G_{2},\quad Q_{h}\gg Q_{crit},\end{eqnarray}$$

where $p_{A}$ is the standing-wave pressure amplitude, and $G_{1}$ and $G_{2}$ contain the spatial distribution of $p_{1}$ and its derivative, as well as the losses multiplying $p_{A}^{2}$ , assumed independent of $T_{m}$ . If $Q_{h}\propto p_{A}^{2}$ , $\text{d}T_{m}/\text{d}x\approx \text{const.}$ and $\widehat{\unicode[STIX]{x0394}T}$ is nearly a constant. To illustrate this relation, we calculated $Q_{h}/|p_{1}|^{2}$ , namely the heat flux input divided by the limit cycle square maximum pressure amplitude, and plotted it against $Q_{h}$ in figure 8(b). Theoretically, such curves tend to infinity at $Q_{h}=Q_{crit}$ ; for practical reasons, only results above a threshold of $|p_{1}|=100~\text{Pa}$ are displayed. The curves quickly decay towards a constant value, in agreement with (4.7). This flat shape of the curves remains as $Q_{h}$ is further increased, and the approximation is valid as long as $p_{1}/p_{m}=O(\unicode[STIX]{x1D700})$ . It is worthwhile noting, that the same behaviour is observed for a conduction-driven system; the general form of (4.6) results in a nearly identical approximation, in which only the values of $G_{1}$ and $G_{2}$ in (4.7) vary.

Figure 9. (a) Scaled temperature difference across the stack, $\widehat{\unicode[STIX]{x0394}T}=\unicode[STIX]{x0394}T/(T_{max}-T_{c})$ , and (b) maximum pressure amplitude, $|p_{1}|$ , as a function of the non-dimensional group $R_{k}/K$ for different values of the partition coefficient $K$ , at $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2$ . (c) Value of $\widehat{\unicode[STIX]{x0394}T}$ for the limit of fast kinetics ( $R_{k}\rightarrow 0$ ) as a function of $K$ for different $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ values. The inset presents an analytic approximation for the critical partition coefficient, $K_{crit}$ , for which a maximum appears, as a function of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . (d) Value of $\widehat{\unicode[STIX]{x0394}T}$ for the limit of ideal sorption ( $K\rightarrow \infty$ ) as a function of $R_{k}/K$ for different $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ values. The inset presents an analytic approximation for the critical value $(R_{k}/K)_{crit}$ for which a minimum appears, as a function of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . The mean pressure is $p_{m}=1$ bar, the cold-side temperature is $T_{c}=20\,^{\circ }\text{C}$ , the scaled stack length and location are $L_{s}/L=0.02$ and $L_{0}/L=0.83$ , respectively, and the heat flux input is $Q_{h}=25~\text{kW}~\text{m}^{-2}$ .

4.2.5 Sorption kinetics

In order to study the effects of reaction kinetics on the system’s limit cycle, we relax the assumption of ideal sorption ( $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}},\unicode[STIX]{x1D702}_{D}\neq 1$ ) and calculate key elements that characterise the limit cycle – $\unicode[STIX]{x0394}T$ and $|p_{1}|$ , as affected by the sorption characteristics. Figure 9(a) presents the scaled temperature difference across the stack, $\widehat{\unicode[STIX]{x0394}T}$ , as a function of $R_{k}/K\equiv \unicode[STIX]{x1D708}/(h^{2}k_{a})$ , the ratio of the viscous time scale and the (forward) sorption reaction, for various values of $K$ . Although drawn against $R_{k}/K$ (instead of simply against $R_{k}$ ), the curves resemble their respective curves from the stability analysis (figure 4): the left end of the curves represents the limit of fast kinetics, in which each curve flattens towards a unique $\widehat{\unicode[STIX]{x0394}T}$ value, determined by the value of $K$ – the boundary’s sorptive affinity for a given gas component. The right end of the curves represents the limit of slow kinetics, a limit to which all curves approach. However, unlike figure 4, at the system limit cycle $K\rightarrow 0$ does not recover the known, conduction-driven thermoacoustics solution (dashed-dotted, red). The oscillatory flow, produced solely by the conduction-driven mechanism in this limit, in the presence of a temperature gradient, results in a hydrodynamic dispersion of heat down the gradient. Since $\text{d}C_{m}/\text{d}x\propto \text{d}T_{m}/\text{d}x$ , a concentration gradient is also imposed on the stack, giving rise to molecular diffusion as well as dispersion of mass down the concentration gradient. This undesired mass transfer requires more heating power to evaporate/desorb more of the reactive component, such that $\text{d}T_{m}/\text{d}x$ decreases, and eventually $\unicode[STIX]{x0394}T_{(K\rightarrow 0)}<\unicode[STIX]{x0394}T_{(\mathit{conduction-driven})}$ .

As in all phase-exchange results, $\widehat{\unicode[STIX]{x0394}T}\leqslant 1$ is bounded by the reactive fluid depletion temperature, $T_{max}$ (see § 4.1.1). Additionally, the stability analysis showed that all curves are further bounded by the limit $K\rightarrow 0$ , recovering the conduction-driven result. Here $K\rightarrow 0$ does not recover the conduction-driven result, and, surprisingly, the curves are not physically bounded by the $K\rightarrow 0$ nor by the conduction-driven lines. Curves of small, finite $K$ values occasionally exceed these lines for sufficiently fast kinetics ( $R_{k}/K\lesssim 10$ ). In order to closely examine this phenomenon, we consider the limit of fast kinetics ( $R_{k}\rightarrow 0$ ), neutralising the effect of reaction kinetics and focus attention to the effect of $K$ and $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ on $\widehat{\unicode[STIX]{x0394}T}$ . We calculate $\widehat{\unicode[STIX]{x0394}T}$ as a function of $K$ , and plot different curves in figure 9(c), each representing a unique $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ value. The $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}=2$ curve matches the $R_{k}/K\rightarrow 0$ limit results in figure 9(a). All the curves in figure 9(c) approach a constant value as $K\rightarrow \infty$ (i.e. ideal sorption), determined by the value of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . Interestingly, some curves display a maximum point while others are cut off as the curve approach $\widehat{\unicode[STIX]{x0394}T}=1$ . Theoretically, all curves exhibit a maximum point, although for certain $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ this point may be obtained at $\widehat{\unicode[STIX]{x0394}T}>1$ or $K<0$ , both representing non-physical cases. These $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , either too small or too large, yield a relatively large $\widehat{\unicode[STIX]{x0394}T}$ even when $K\rightarrow \infty$ (see figure 6 a, where $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}=\unicode[STIX]{x1D702}_{D}=1$ is assumed), and as the boundary’s affinity for sorption is decreased, the short distance to the $\widehat{\unicode[STIX]{x0394}T}=1$ line is insufficient for a maximum to appear. This result is confirmed by an analytic approximation detailed in appendix B. The results of this approximation, presented in the inset of figure 9(c), show the qualitative trend of $K_{crit}$ , the $K$ value leading to a maximum in the curves of figure 9(c), as a function of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ for $Sc=0.6$ . As clearly seen, a positive (physical) value of $K_{crit}$ (shaded) is obtained only for a finite, intermediate range of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . This result is reproduced for any $Sc>0$ .

The curves in figure 9(a) generally display two minima: the first appears at $R_{k}/K\approx 50$ in all curves, and the second appears at $R_{k}\approx 10^{-1}$ , marked by a dotted line connecting the points. While the two minima are suppressed in the limit of $K\rightarrow 0$ , only the second minimum is flattened for $K\rightarrow \infty$ . For $K\rightarrow 0$ , $\unicode[STIX]{x1D702}_{n}\rightarrow \infty$ and $\unicode[STIX]{x1D702}_{\unicode[STIX]{x1D708}}/\unicode[STIX]{x1D702}_{D}\rightarrow (1-F_{\unicode[STIX]{x1D708}})/(1-F_{D})$ are recovered, regardless of the value of $R_{k}$ . This limit truly suppresses the effect of reaction kinetics since the material now has no sorption affinity, and the minima vanish. If $K$ is small but finite, the value of $R_{k}$ dramatically affects $\widehat{\unicode[STIX]{x0394}T}$ ; the second minimum ( $R_{k}\approx 10^{-1}$ ) appears also in the stability analysis calculations (figure 4), and is a result of the boundary’s limited sorption capacity, that may be compensated by a regulated delay in desorption rate (see detailed explanation in § 4.1.3). When $K>10^{-1}$ , this minima disappears, leaving a short flat region before increasing towards the limit of $K\rightarrow 0$ . The first minimum ( $R_{k}/K\approx 50$ ), absent in the stability analysis calculations, remains for all $K$ . Investigating the nature of this minimum, we calculate additional curves of $\widehat{\unicode[STIX]{x0394}T}$ versus $R_{k}/K$ for the limit of ideal sorption ( $K\rightarrow \infty$ ) at various $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ values. All these curves, presented in figure 9(d), show at least one minimum. For $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}>2$ , the minimum shifts towards small $R_{k}/K$ . Here the boundary has a high affinity for sorption, and the value of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ sets the oscillation rate, for which an optimum of the reaction kinetics exists, such that $\widehat{\unicode[STIX]{x0394}T}$ is minimised. As $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ is increased, the optimum is obtained for faster kinetics, or smaller $R_{k}/K$ values. Theoretically, there exists a critical value $(\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}})_{crit}$ , for which $R_{k}/K=0$ ; if $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}>(\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}})_{crit}$ , the minimum point lies at $R_{k}/K<0$ – a non-physical result. Practically, since $\widehat{\unicode[STIX]{x0394}T}\leqslant 1$ sets a physical constraint, $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}>(\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}})_{crit}$ curves do not yield any valid limit cycle results, and therefore curves showing no minimum cannot be drawn. For $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}<2$ , the second minimum point, suppressed for $K\rightarrow \infty$ at $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}\geqslant 2$ , reappears and as $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ is further decreased the first minimum vanishes, leaving a single minimum that shifts towards large $R_{k}/K$ . If $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ is dramatically decreased, the minimum is pushed towards very large values of $R_{k}/K$ , in which the curves approach the $K\rightarrow 0$ limit, where the minimum is, again, suppressed. As before, the $\widehat{\unicode[STIX]{x0394}T}\leqslant 1$ constraint makes it impossible to draw such curves, as their results exceed $\widehat{\unicode[STIX]{x0394}T}=1$ well before the minimum disappears. In the absence of numerical evidence to support our claim, we seek an analytic approximation proving that a minimum is obtained for only a finite range of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ . The details of this approximation appear in appendix B, and its results, presented in the inset of figure 9(d), show $(R_{k}/K)_{crit}$ , namely the $R_{k}/K$ value leading to a minimum in the curves of figure 9(d), as a function of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ for $Sc=0.6$ . These results indicate that a minimum in $\widehat{\unicode[STIX]{x0394}T}$ should only be obtained for a finite, intermediate range of $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D708}}$ , since $(R_{k}/K)_{crit}<0$ does not represent a physically valid solution.

The effect of reaction kinetics on a thermoacoustic system vary according to the system’s objective: for minimising $\unicode[STIX]{x0394}T$ , controlling $R_{k}$ proves to be beneficial even if $K\rightarrow \infty$ , but works to the contrary if one wishes to maximise ${\dot{W}}$ , since the minimum in $\widehat{\unicode[STIX]{x0394}T}$ imposes an undesired minimum on $\unicode[STIX]{x1D703}\equiv ||\unicode[STIX]{x1D719}|-90^{\circ }|$ such that ${\dot{W}}$ is also lowered. This effect is clearly seen in figure 9(b), showing the maximum pressure amplitude, $|p_{1}|$ , as a function of $R_{k}/K$ . The nearly monotonically increasing nature of the curves, linking the limits of fast and slow kinetics, is disturbed by a local maximum point, exactly matching the minimum points at $R_{k}/K\approx 50$ in figure 9(a). The increase in $|p_{1}|$ cannot be attributed to the decrease in $\widehat{\unicode[STIX]{x0394}T}$ directly, since the latter decreases the system’s thermal potential. However, the subsequent effect of the minimum in $\widehat{\unicode[STIX]{x0394}T}$ on $\unicode[STIX]{x1D719}$ alters the wave form such that a nearly pure standing wave is obtained, for which $|p_{1}|$ increases in spite of the decrease in ${\dot{W}}$ . All values of $|p_{1}|$ here are smaller than the conduction-driven reference values. Although $\widehat{\unicode[STIX]{x0394}T}$ may exceed the conduction-driven value for $R_{k}\rightarrow 0$ , thus maximising the system thermal potential, the wave form deviates from a standing wave appreciably, resulting in a small, minimal value of $|p_{1}|$ . At the local maximum ( $R_{k}/K\approx 50$ ), the resemblance of wave form to a pure standing wave increases $|p_{1}|$ , such that for finite $K$ values the $K\rightarrow 0$ limit is breached. The local minimum in $\widehat{\unicode[STIX]{x0394}T}$ decreases the system thermal potential, such that the maximal phase-exchange pressure amplitude ( $|p_{1}|\approx 2000~\text{Pa}$ ) is considerably smaller than the conduction-driven value ( $|p_{1}|=2488~\text{Pa}$ ).