Numerical investigation of porous stack for a solar-powered thermoacoustic refrigerator

Abstract

This article proposes a numerical analysis for performance improvement of the stack, which represents a crucial element on the solar-powered thermoacoustic refrigerator. The stack is considered as a saturated parallelepiped homogeneous porous media. Numerical simulation of the flow and the heat transfer in the thermoacoustic refrigerator is carried out. The physical flow is governed by the modified Darcy–Brinkmann–Forchheimer model. The governing equations are solved numerically using the lattice Boltzmann method. Furthermore, the local thermal equilibrium assumption is applied to examine heat transfer. Particular attention is paid a new form of the lattice Boltzmann equation system describing the flow and the heat transfer in porous media and fluid regions. The effects of several parameters characterizing the thermal behaviour in the porous medium are studied. The parametric results lead to the optimization of the porous media form. The presence of the viscous dissipation term in the heat transfer formulation within the thermoacoustic system is particularly highlighted, due to its significant effects introduced by the Eckert number.

Introduction

Throughout the last decades, researchers have tried to turn all the systems towards green energy, in order to become environmentally friendly.^1–6 Because of high energy consumption of conventional refrigerator, the exploitation of renewable energy is an appropriate solution to increase the overall efficiency. Therefore, the thermoacoustic refrigerator (TAR) is proposed ⁷ (known as the thermoacoustic cooler or the thermoacoustic heat pump), as one of the recent advances in the field. The operating principle of the TAR is to convert acoustic energy into a temperature difference. ⁸ Generally, the operating gas used in thermoacoustic systems is either binary mixtures of rare gases that are harmless to the environment, ⁹ or sometimes, the ambient air itself thanks to its low cost and availability.

In 1777, Byron Higgins was the first scientist who observed the thermoacoustic effect. ⁹ In fact, he noticed that the heat may induce acoustic oscillation (sound). Later in 1896, Lord Rayleigh ¹⁰ gave his first qualitative interpretation for the Sondhauss tube operation. However, the inverse physical phenomenon, which is the creation of a temperature gradient from an intense acoustic wave, was demonstrated in the work of Merkli and Thomann ¹¹ in 1975 and Wheatley et al. ¹² in 1983. Since that discovery, the scientific research in thermoacoustic engines started to have concrete applications.^13–15 Nevertheless, industrialization and marketing were not possible because of low efficiency. ¹⁶ One of the solutions found to enhance the TAR efficiency is to couple it with a green energy source. Recent researches proved that photovoltaic and wind energy are the lowest cost technologies used to produce electricity from renewable resources.^17,18 Therefore, combining TAR and the photovoltaic system should be worthwhile. Coupling TAR and solar energy has already been studied in recent research papers^13–16,19 particularly in the thermo-driven TAR case. The solar-powered TAR is studied in detail by N Yassen. ²⁰ Moreover, he designs a stand-alone power system to feed a sinuous wave generator, which produces a standing wave in the thermoacoustic resonator. T Konaina and N Yassen ²¹ focused on the TAR for domestic usage powdered by solar energy. Consequently, they developed a computer code provides the system dimensions and the optimal design pressure as a function of the cooling capacity to ensure a convenient sound power. The code showed good reliability by comparing its results with some manufactured systems.

A detailed explanation of the operating mode of thermoacoustic systems was provided by Swift. ⁸ TAR has a simple architecture without moving parts, which makes it robust and reliable. In addition, it offers opportunities for production and miniaturization at low cost. It consists of four components: the resonator, the gas filling the resonator, the loudspeaker and the stack (called regenerator). Since the stack is the most crucial part of the TAR, particular attention will be paid to its improvement. The stack geometry is generally formed by thin plates, hence its name. Adeff et al. ²² experimentally proved the efficiency of using the porous medium (reticulated vitreous carbon (RVC)) as the core of thermoacoustic prime movers and refrigerator. Furthermore, the RVC is a low-cost solution, and its manufacture simplicity makes its usage suitable for all types of thermoacoustic devices. The RVC was simulated in previous works as a porous media.^16,22 Syeda H Tasnim et al. ²³ proposed an analytical study of the velocity and the temperature fields in a stack composed of parallel plates where the gaps between them are filled by porous media.

The discretization method used in this article is the lattice Boltzmann method (LBM) due to its low simulation time, its simplicity and accuracy. ²⁴ The LBM was used in previous work, on one hand, in the simulation of the TAR ²⁵ and, on the other hand, in the simulation of the flow and the heat transfer in porous media. Therefore, the LBM is chosen in this article to simulate the TAR having a porous media as a stack. Indeed, the LBM was adopted to analyse the heat transfer in the TAR. ²⁵ In the same work, the software ‘PowerFLOW’ was applied, which is based on the LBM to simulate the TAR with a stack formed of parallel plates. The LBM proved its performance in thermoacoustic systems. Besides, many research papers proved the performance of the LBM in the simulation of flow and heat transfer in porous media. For these reasons, the LBM was applied to simulate flow and heat transfer inside TAR having a porous media stack.

The main objective of this article is to investigate the viscous dissipation inside the porous stack. Conjugate problems, including the viscous fluid and the porous medium, are among important troubles for a whole bunch of technologies. For instance, air heaters, solid matrix, heat exchangers, electronic cooling, drying processes ²⁶ and pulse tube cryocooler ²⁷ are exposed to this kind of issues. In fact, enhancing heat transfer is the major advantage of porous media.^28–31

However, the detailed description of the flow and thermal fields of the acoustic waves, inside the thermoacoustic stack, modelled as a porous media, are not available in the literature concerning a rectangular resonator. In this work, we developed a mathematical model of the fluid motion and the heat transfer in three parts of the solar-powered TAR. Then, a new formulation of the Lattice Boltzmann equations is presented, describing the flow and the heat transfer in both porous media and fluid regions. Furthermore, a numerical study is carried out to investigate the effect of some relevant parameters on the thermal performance of the solar-powered TAR. These parameters include the Darcy number, the porous media porosity, the heat capacity ratio, the heat conductivity ratio and the Eckert number. The remainder of this article is organized as follows: section ‘Problem specification’ defines the problem specifications, section ‘Validation’ depicts the validation procedure, section ‘Results and discussion’ is devoted to results and discussion and section ‘Conclusion’ concludes the article.

Problem specification

The system consists of a TAR powered by grid-connected photovoltaic panels as shown in Figure 1(a). The photovoltaic panel produces direct current (DC) that is converted into alternating current (AC) using the inverter, Indexed 2 in Figure 1(a). Two power metres are used to record the power generated by the photovoltaic panel and the power absorbed by the TAR. On a sunny day, the refrigerator consumes the electric energy generated by the photovoltaic panel. However, on cloudy days or at night, the TAR is powered by the grid. Moreover, the extra produced electricity by solar panels is sold to the power company. The loudspeaker converts electric power into acoustic power. The latter is then transformed into a temperature difference that can be explored in different ways such as refrigeration appliances or biomedical devices.

Figure 1.

Schematic diagram ((a) of the system and (b) of the thermoacoustic refrigerator).

In this investigation, the simulation only concerns the TAR. The photovoltaic panel sizing is devoted to future works. To this end, the TAR is designed according to a two-dimensional model as shown in Figure 1(b). The basic components of the system are the resonator, the loudspeaker and the stack. The resonator is a quarter wavelength tubes whose boundaries are exposed to the ambient temperature. The resonator has a closed extremity on the right side as shown in Figure 1(a). In the left resonator’s extremity is located a loudspeaker, ensuring acoustic power. The resonator is filled by air under atmospheric pressure where the Prandtl number is equal to 0.71, and the average temperature T_m is equal to 300°K. The fluid is animated by the wave emitted from the loudspeaker. The resonance frequency of the tube is 200 Hz. The crucial part of the TAR is the stack which is made from a homogeneous porous medium. The stack centre is located at the position x = x_c and the stack length is l_s, as shown in Figure 1(b). The values of l_s and x_c change throughout the article.

Governing equations

The governing equations, characterizing describing the flow and the heat transfer inside the TAR, are presented in the following section. It is worth noting that the local thermal equilibrium between the fluid and the solid phases assumption is taken into consideration in this article. The thermophysical properties of the fluid are constant, the fluid is assumed to be low compressible and the porous medium is homogeneous, isotropic and saturated.

A binary parameter is used to write the equations for both sections, fluid and porous, in a unified way. The binary parameter Θ is defined as^32,33

$$\mathrm{\Theta }=\{\begin{array}{cc}0,& \mathrm{where}\mathrm{the}\mathrm{air}\mathrm{reigns}\\ 1,& \mathrm{in}\mathrm{the}\mathrm{porous}\mathrm{media}\end{array}$$

(1)

Based on the previous assumptions and the binary parameter (equation (1)), the mass, the momentum and the energy equations are written as follows^23,33,34

• Continuity equation

$$\frac{\partial \mathrm{u}}{\partial \mathrm{x}}+\frac{\partial \mathrm{v}}{\partial \mathrm{x}}=0$$

(2)

• Momentum equation

$$\begin{array}{c}\frac{\partial \mathrm{u}}{\partial \mathrm{t}}+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1){\mathrm{u}}^2\right]}{\partial \mathrm{x}}+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1)\mathrm{uv}\right]}{\partial \mathrm{y}}\hfill \\ =-\frac{1}{\mathrm{\rho }}\frac{\partial \mathrm{p}}{\partial \mathrm{x}}+\nu (\frac{{\partial }^2\mathrm{u}}{\partial {\mathrm{x}}^2}+\frac{\partial \mathrm{u}}{\partial \mathrm{y}})-\mathrm{\Theta }(\frac{\mathrm{\epsilon }\nu }{\mathrm{K}}+\frac{\mathrm{\epsilon }{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{K}}}\left|\mathrm{u}\right|)\mathrm{u}\hfill \end{array}$$

(3)

$$\begin{array}{c}\frac{\partial \mathrm{v}}{\partial \mathrm{t}}+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1)\mathrm{uv}\right]}{\partial \mathrm{x}}+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1)\mathrm{v}\right]}{\partial \mathrm{y}}\hfill \\ =-\frac{1}{\mathrm{\rho }}\frac{\partial \mathrm{p}}{\partial \mathrm{y}}+\nu (\frac{{\partial }^2\mathrm{v}}{\partial {\mathrm{x}}^2}+\frac{\partial \mathrm{v}}{\partial \mathrm{y}})-\mathrm{\Theta }(\frac{\mathrm{\epsilon }\nu }{\mathrm{K}}+\frac{\mathrm{\epsilon }{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{K}}}\left|\mathrm{v}\right|)\mathrm{v}\hfill \end{array}$$

(4)

• Energy equation

$$\begin{array}{c}[\mathrm{\Theta }((\mathrm{\epsilon }+(1-\mathrm{\epsilon })\frac{{\left(\mathrm{\rho }{\mathrm{C}}_{\mathrm{p}}\right)}_{\mathrm{s}}}{{\left(\mathrm{\rho }{\mathrm{C}}_{\mathrm{p}}\right)}_{\mathrm{f}}}-1)+1]\frac{\partial \mathrm{T}}{\partial \mathrm{t}}\hfill \\ +(\frac{\partial \left(\mathrm{uT}\right)}{\partial \mathrm{x}}+\frac{\partial \left(\mathrm{vT}\right)}{\partial \mathrm{y}})\hfill \\ =\mathrm{\alpha }[\mathrm{\Theta }((\mathrm{\epsilon }+(1-\mathrm{\epsilon })\frac{{\mathrm{k}}_{\mathrm{s}}}{{\mathrm{k}}_{\mathrm{f}}})-1)+1](\frac{{\partial }^2\mathrm{T}}{\partial {\mathrm{x}}^2}+\frac{{\partial }^2\mathrm{T}}{\partial {\mathrm{y}}^2})\hfill \\ +\mathrm{\mu }[2[{\left(\frac{\partial \mathrm{u}}{\partial \mathrm{x}}\right)}^2+{\left(\frac{\partial \mathrm{v}}{\partial \mathrm{y}}\right)}^2]+{[\frac{\partial \mathrm{u}}{\partial \mathrm{y}}+\frac{\partial \mathrm{v}}{\partial \mathrm{x}}]}^2]\hfill \\ +\mathrm{\Theta }(\frac{\mathrm{\epsilon }}{\mathrm{K}}+\frac{\mathrm{\epsilon }{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{K}}}\left|\mathrm{u}\right|)(\mathrm{u}+{\mathrm{v}}^2)\hfill \end{array}$$

(5)

where u, v, p and T are the axial and the transversal velocity components, the pressure and the temperature, respectively.

ρ, ν, α, C_p, k, ε, K and F_ε are the density, the kinetic viscosity, the thermal diffusivity coefficient, the heat capacity, the heat conductivity, the porous media porosity, the permeability and the Forchheimer number, respectively. Subscripts ‘s’ and ‘f’ design the porous media and the fluid part of the solid matrix, respectively. F_ε is defined as follows

$${\mathrm{F}}_{\mathrm{\epsilon }}=\frac{1,75}{\sqrt{150{\mathrm{\epsilon }}^3}}$$

(6)

Dimensionless governing equations

The dimensionless governing equations authorize the introduction and the study of similar problems with the same dimensionless parameter. The dimensionless reference quantities are defined as

$$\begin{array}{c}\mathrm{X}=\frac{\mathrm{x}}{\mathrm{H}},\mathrm{Y}=\frac{\mathrm{y}}{\mathrm{H}},\mathrm{U}=\frac{\mathrm{u}}{{\mathrm{u}}_0},\mathrm{V}=\frac{\mathrm{v}}{{\mathrm{u}}_0},\hfill \\ \mathrm{\tau }=\frac{\mathrm{t}{\mathrm{u}}_0}{\mathrm{H}},\theta =\frac{\mathrm{T}}{{\mathrm{T}}_{\mathrm{m}}},\mathrm{P}=\frac{\mathrm{p}}{\mathrm{\rho }{\mathrm{u}}_0^2}\hfill \end{array}$$

(7)

The dimensionless governing equations (2)–(5) are given by the following formulations

• Continuity equation

$$\frac{\partial \mathrm{U}}{\partial \mathrm{X}}+\frac{\partial \mathrm{V}}{\partial \mathrm{Y}}=0$$

(8)

• Momentum equation

$$\begin{array}{c}\frac{\partial \mathrm{U}}{\partial }+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1){\mathrm{U}}^2\right]}{\partial \mathrm{X}}+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1)\mathrm{UV}\right]}{\partial \mathrm{Y}}\hfill \\ =-\frac{\partial \mathrm{P}}{\partial \mathrm{X}}+\frac{1}{\mathrm{Re}}(\frac{{\partial }^2\mathrm{U}}{\partial {\mathrm{X}}^2}+\frac{{\partial }^2\mathrm{U}}{\partial {\mathrm{Y}}^2})+\mathrm{\Theta }{\mathrm{F}}_{\mathrm{X}}\hfill \end{array}$$

(9)

$$\begin{array}{c}\frac{\partial \mathrm{V}}{\partial }+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1)\mathrm{UV}\right]}{\partial \mathrm{X}}+\frac{\partial \left[\mathrm{\Theta }((\frac{1}{\mathrm{\epsilon }}-1)+1){\mathrm{V}}^2\right]}{\partial \mathrm{Y}}\hfill \\ =-\frac{\partial \mathrm{P}}{\partial \mathrm{Y}}+\frac{1}{\mathrm{Re}}(\frac{{\partial }^2\mathrm{V}}{\partial {\mathrm{X}}^2}+\frac{{\partial }^2\mathrm{V}}{\partial {\mathrm{Y}}^2})+\mathrm{\Theta }{\mathrm{F}}_{\mathrm{y}}\hfill \end{array}$$

(10)

where

$$\{\begin{array}{c}\multicolumn{1}{c}{{\mathrm{F}}_{\mathrm{x}}=-(\frac{\mathrm{\epsilon }}{\mathrm{DaRe}}+\frac{\mathrm{\epsilon }{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{Da}}}\left|\mathrm{U}\right|)\mathrm{U}}\\ \multicolumn{1}{c}{{\mathrm{F}}_{\mathrm{y}}=-(\frac{\mathrm{\epsilon }}{\mathrm{DaRe}}+\frac{\mathrm{\epsilon }{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{Da}}}\left|\mathrm{V}\right|)\mathrm{V}}\end{array}$$

(11)

• Energy equation

$$\begin{array}{c}[\mathrm{\Theta }((\mathrm{\epsilon }+(1-\mathrm{\epsilon }){\mathrm{R}}_{\mathrm{c}}-1)+1]\frac{\partial \theta }{\partial \mathrm{\tau }}+(\frac{\partial \left(\mathrm{U}\theta \right)}{\partial \mathrm{X}}+\frac{\partial \left(\mathrm{V}\theta \right)}{\partial \mathrm{Y}})\hfill \\ =[\mathrm{\Theta }((\mathrm{\epsilon }+(1-\mathrm{\epsilon }){\mathrm{R}}_{\mathrm{k}})-1)+1]\frac{1}{\mathrm{RePr}}(\frac{{\partial }^2\theta }{\partial {\mathrm{X}}^2}+\frac{{\partial }^2\theta }{\partial {\mathrm{Y}}^2})\hfill \\ +\mathrm{\Theta }\frac{\mathrm{Ec}}{\mathrm{Re}}[2[{\left(\frac{\partial \mathrm{U}}{\partial \mathrm{X}}\right)}^2+{\left(\frac{\partial \mathrm{V}}{\partial \mathrm{Y}}\right)}^2]+{[\frac{\partial \mathrm{U}}{\partial \mathrm{Y}}+\frac{\partial \mathrm{V}}{\partial \mathrm{X}}]}^2]\hfill \\ +\mathrm{\Theta }(\frac{\mathrm{Ec}}{\mathrm{ReDa}}+\frac{\mathrm{\epsilon }\mathrm{Ec}{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{Da}}}\left|\mathrm{U}\right|)({\mathrm{U}}^2+V^2)\hfill \end{array}$$

(12)

The dimensionless form often reveals the dimensionless parameters: the Reynolds number $\mathrm{Re}=({\mathrm{u}}_0\mathrm{H})/\mathrm{\nu }$ , the Darcy number $\mathrm{Da}=\mathrm{K}/\mathrm{H}$ , Prandtl number $\mathrm{Pr}=\nu /\mathrm{\alpha }$ , Eckert number $\mathrm{Ec}={\mathrm{u}}_0^2/({\mathrm{C}}_{\mathrm{p}}\mathrm{\Delta }{\mathrm{T}}_{\mathrm{ref}})$ , the heat capacity ratio ${\mathrm{R}}_{\mathrm{c}}=(\mathrm{\rho }{\mathrm{C}}_{\mathrm{p}})\mathrm{s}/(\mathrm{\rho }{\mathrm{C}}_{\mathrm{p}}{)}_{\mathrm{f}}$ and the heat conductivity ratio ${\mathrm{R}}_{\mathrm{k}}={\mathrm{k}}_{\mathrm{s}}/{\mathrm{k}}_{\mathrm{f}}$ .

Initial and boundary conditions

Initially, the fluid is motionless and the temperature inside the resonator is fixed at the ambient temperature ${\mathrm{T}}_{\mathrm{m}}$ :

• At the loudspeaker side, the velocity is expressed as follows ³⁵

$$\mathrm{u}\left(\mathrm{t}\right)={\mathrm{u}}_0\sin ({\mathrm{k}}_{\mathrm{w}}·{\mathrm{x}}_{\mathrm{B}})\cos \left(\mathrm{\omega }\mathrm{t}\right)$$

(13)

where ${\mathrm{u}}_0,{\mathrm{k}}_{\mathrm{w}},\mathrm{\omega }$ and t refer to the maximum amplitude of the acoustic velocity, the wavenumber ${\mathrm{k}}_{\mathrm{w}}=2\pi /\lambda$ , the wave pulsation $\mathrm{\omega }=2\pi \mathrm{f}$ and the simulation time, respectively, with λ and f are the wavelength and the sound wave frequency, respectively.

Here x_B is the axial position of the boundary measured from the resonator closed end. So, at the loudspeaker side, $\mathrm{x}=\mathrm{L}=\lambda /4$

$$\mathrm{u}={\mathrm{u}}_{\mathrm{o}}\cos \left(2\pi \mathrm{ft}\right),{\mathrm{u}}_{\mathrm{o}}=\frac{{\mathrm{P}}_{\mathrm{A}}}{\mathrm{\rho }{\mathrm{c}}_{\mathrm{s}}}$$

(14)

The temperature at x = 0 is written as ³⁶

$$\mathrm{T}={\mathrm{T}}_{\mathrm{m}}+\frac{{\mathrm{P}}_{\mathrm{A}}}{\mathrm{\rho }{\mathrm{C}}_{\mathrm{p}}}\cos ({\mathrm{k}}_{\mathrm{w}}·{\mathrm{x}}_{\mathrm{s}})\cos \left(\mathrm{\omega }\mathrm{t}\right)$$

(15)

where T_m is the ambient temperature, P_A is the pressure amplitude at the pressure antinode and ρ and C_p are the density and the specific heat at a constant pressure of the fluid, respectively. k_w and ω are the wavenumber and the pulsation as defined above.

• In the interfaces between the fluid region and the stack, the continuity of velocity and normal stress prevailing are adopted at the interfaces ³⁷

$$\{\begin{array}{c}{\mathrm{u}}_{|\mathrm{fld}}={\mathrm{u}}_{|\mathrm{prs}}={\mathrm{u}}_{|\mathrm{int}}\hfill \\ \frac{1}{\mathrm{\epsilon }}{\frac{\partial \mathrm{u}}{\partial \mathrm{x}}}_{|\mathrm{prs}}={\frac{\partial \mathrm{u}}{\partial \mathrm{x}}}_{|\mathrm{fld}}\hfill \\ \frac{1}{\mathrm{\epsilon }}{\frac{\partial \mathrm{v}}{\partial \mathrm{y}}}_{|\mathrm{prs}}={\frac{\partial \mathrm{v}}{\partial \mathrm{y}}}_{|\mathrm{fld}}\hfill \end{array}$$

(16)

where subscripts ‘int’, ‘fld’ and ‘prs’ refer to the interface layer, the fluid region and the porous region, respectively.

All the resonator solid walls are exposed to the ambient temperature and the no-slip boundary conditions are imposed

• At x = L

$$\{\begin{array}{c}\mathrm{u}=\mathrm{v}=0\hfill \\ \mathrm{T}(\mathrm{x}=\mathrm{L})={\mathrm{T}}_{\mathrm{m}}\hfill \end{array}$$

(17)

• At the horizontal sides of the resonator (y = 0 and y = H)

$$\{\begin{array}{c}\mathrm{u}=\mathrm{v}=0\hfill \\ \mathrm{T}(\mathrm{y}=0)={\mathrm{T}}_{\mathrm{m}}\hfill \\ \mathrm{T}(\mathrm{y}=\mathrm{H})={\mathrm{T}}_{\mathrm{m}}\hfill \end{array}$$

(18)

Numerical procedure

The numerical method adopted in this article is the LBM based on the Bhatnagar–Gross–Krook (BGK) approximation. ³⁸ The LBM is a new and promising numerical method applied to fluid flow simulation. It is used in many research fields, including fluid flows in porous media, ³⁹ oscillatory flow, ⁴⁰ multiphase flow, ⁴¹ magneto-hydrodynamics ⁴² and fluid fluxes at the nano and microscopic scales. ⁴³ Moreover, the LBM has proved that it is faster than the traditional computational fluid dynamics (CFD) methods. ⁴⁴ The double distribution function (DDF) is used to discretize the flow and the heat transfer inside the computational domain. In fact, its mean principal is to define the temperature and the density flow fields using two distribution functions. The first distribution function named ‘f’ is used to simulate the density and velocity fields, which are the particle density distribution function. The second one is dedicated to the thermal field and it is named ‘g’. Furthermore, the LBM is based on four steps, ⁴³ the collision step, the streaming step, the boundary conditions and the computing of macroscopic proprieties which are the density, the velocity and the temperature. ²⁴ Besides, a two-dimensional lattice was been chosen with nine discrete speeds on a square grid of steps Δx = Δy = 1. The model is called D2Q9 and it is illustrated in Figure 2.

Figure 2.

D2Q9 model.

The velocity field

In the LBM, for the D2Q9 model, the fluid velocity ‘u’ and the fluid density ‘ρ’ are calculated using the following equations. The system equation systems of fluid and porous media parts were combined.^34,43 The binary parameter Θ defined in equation (1) has been used

$$\mathrm{u}(\mathrm{M},\mathrm{t})=\frac{\mathrm{v}(\mathrm{M},\mathrm{t})}{\mathrm{\Theta }({\mathrm{c}}_0+\sqrt{{\mathrm{c}}_0^2+{\mathrm{c}}_1\left|\mathrm{v}(\mathrm{M},\mathrm{t})\right|}-1)+1}$$

(19)

$$\mathrm{\rho }(\mathrm{M},\mathrm{t})=\sum _{\mathrm{i}=0,8}{\mathrm{f}}_{\mathrm{i}}(\mathrm{M},\mathrm{t})$$

(20)

where

$$\mathrm{v}(\mathrm{M},\mathrm{t})=\frac{1}{\mathrm{\rho }(\mathrm{M},\mathrm{t})}\sum _{\mathrm{i}=0,8}{\mathrm{e}}_{\mathrm{i}}{\mathrm{f}}_{\mathrm{i}}(\mathrm{M},\mathrm{t})$$

(21)

$${\mathrm{c}}_0\frac{1}{2}(1+\mathrm{\epsilon }\frac{\mathrm{\Delta }\mathrm{t}}{2}\frac{\mathrm{\mu }}{\mathrm{\rho }(\mathrm{M},\mathrm{t})\mathrm{K}}),{\mathrm{c}}_1=\mathrm{\epsilon }\frac{\mathrm{\Delta }\mathrm{t}}{2}\frac{{\mathrm{F}}_{\mathrm{\epsilon }}}{\sqrt{\mathrm{K}}}$$

(22)

where v(M, t) is a temporal velocity;^26,34,45 ${\mathrm{c}}_0$ and ${\mathrm{c}}_1$ are two parameters, ²⁶ defined in (equation (22)) and e_i is the discrete velocity of a single particle. The value of ${\mathrm{e}}_{\mathrm{i}}$ ³⁴ is

$${\mathrm{e}}_{\mathrm{i}}=\{\begin{array}{cc}(0,0)\hfill & \mathrm{i}=0\hfill \\ \multicolumn{1}{c}{(\cos \left[\frac{(\mathrm{i}-1)\mathrm{\pi }}{2}\right],\sin \left[\frac{(\mathrm{i}-1)\mathrm{\pi }}{2}\right])}& \mathrm{i}=1,2,3,4\hfill \\ \multicolumn{1}{c}{(\sqrt{2}\cos \left[\frac{(2\mathrm{i}-9)\mathrm{\pi }}{4}\right],\sqrt{2}\sin \left[\frac{(2\mathrm{i}-9)\mathrm{\pi }}{4}\right])}& \mathrm{i}=5,6,7,8\hfill \end{array}$$

(23)

where f_i(M, t) is the density distribution function of particle located in the space position M(x, y) at the time step ‘t’.^34,43

The determination of the density distribution function over time and within the calculation domain is effectuated by the collision step. The BGK approximation is used^38,43

$${\mathrm{f}}_{\mathrm{i}}(\mathrm{M},\mathrm{t}+\mathrm{\Delta }\mathrm{t})={\mathrm{f}}_{\mathrm{i}}(\mathrm{M},\mathrm{t})[1-\frac{1}{{\mathrm{\tau }}_{\mathrm{f}}}]+\frac{{\mathrm{f}}_{\mathrm{i}}^{\mathrm{eq}}(\mathrm{M},\mathrm{t})}{{\mathrm{\tau }}_{\mathrm{f}}}+\mathrm{\Theta }\mathrm{\Delta }\mathrm{t}{\mathrm{F}}_{\mathrm{ei}}$$

(24)

where $\mathrm{\Delta }\mathrm{t}$ , ${\mathrm{\tau }}_{\mathrm{f}}$ and ${\mathrm{f}}_{\mathrm{i}}^{\mathrm{eq}}$ are the time increment, the relaxation time and the local equilibrium distribution function,^39,43 respectively, and these are defined as

$${\mathrm{\tau }}_{\mathrm{f}}=3-0.5$$

(25)

$$\begin{array}{c}{\mathrm{f}}_{\mathrm{i}}^{\mathrm{eq}}(\mathrm{M},\mathrm{t})={\mathrm{\omega }}_{\mathrm{i}}\mathrm{\rho }(\mathrm{M},\mathrm{t})\hfill \\ [1+3{\mathrm{e}}_{\mathrm{i}}·\mathrm{u}(\mathrm{M},\mathrm{t})+\frac{9{({\mathrm{e}}_{\mathrm{i}}·\mathrm{u}(\mathrm{M},\mathrm{t}))}^2}{2(\mathrm{\Theta }(\mathrm{\epsilon }-1)+1)}-\frac{3{\mathrm{u}}^2(\mathrm{M},\mathrm{t})}{2(\mathrm{\Theta }(\mathrm{\epsilon }-1)+1)}]\hfill \end{array}$$

(26)

where ν is the kinetic viscosity of the fluid and ${\mathrm{\omega }}_{\mathrm{i}}$ are the weighting factors for corresponding distribution functions ${\mathrm{f}}_{\mathrm{i}}$ . For the D2Q9 model, the weights are given as ³⁴

$${\mathrm{\omega }}_{\mathrm{i}}=\{\begin{array}{cc}\frac{4}{9}\hfill & \mathrm{i}=0\hfill \\ \frac{1}{9}\hfill & \mathrm{i}=1,2,3,4\hfill \\ \frac{1}{36}\hfill & \mathrm{i}=5,6,7,8\hfill \end{array}$$

(27)

${\mathrm{F}}_{\mathrm{ei}}$ is the forcing term due to the proprieties of the porous media ²⁴

$${\mathrm{F}}_{\mathrm{ei}}=\frac{\mathrm{F}({\mathrm{e}}_{\mathrm{i}}-\mathrm{u}(\mathrm{M},\mathrm{t})){\mathrm{f}}_{\mathrm{i}}^{\mathrm{eq}}}{\mathrm{R}{\mathrm{T}}_0}$$

(28)

where $\mathrm{F}={\mathrm{F}}_{\mathrm{x}}+{\mathrm{F}}_{\mathrm{y}}$ is the force term (equation (11)), $\mathrm{R}=1$ for the D2Q9 model and ${\mathrm{T}}_0$ is determined by the lattice Boltzmann model and equal to ‘1/3’. ²⁴

The next phase after the collision step is the streaming step ⁴³

$${\mathrm{f}}_{\mathrm{i}}(\mathrm{M}+{\mathrm{e}}_{\mathrm{i}}\mathrm{\Delta }\mathrm{t},\mathrm{t}+\mathrm{\Delta }\mathrm{t})={\mathrm{f}}_{\mathrm{i}}(\mathrm{M},\mathrm{t}+\mathrm{\Delta }\mathrm{t})$$

(29)

Temperature field

In the D₂Q₉ scheme, the temperature field (T(M, t)) in all the computational domain is given by equation (30). These equations are developed by combining the equations in a fluid region ⁴³ and those related to the porous media section ³⁴

$$\begin{array}{c}(\mathrm{\Theta }(\sigma -1)+1)\mathrm{T}(\mathrm{M},\mathrm{t})\hfill \\ =\sum _{\mathrm{i}=0,8}{\mathrm{g}}_{\mathrm{i}}(\mathrm{M},\mathrm{t}),\sigma =\mathrm{\epsilon }+(1-\mathrm{\epsilon }){\mathrm{R}}_{\mathrm{c}}\hfill \end{array}$$

(30)

where ${\mathrm{g}}_{\mathrm{i}}$ is the temperature distribution function. The temperature distribution function is calculated by the same steps as with the calculation of density distribution function. The collusion and the streaming steps are described by equations (31) and (36), respectively^24,34,43

$${\mathrm{g}}_{\mathrm{i}}(\mathrm{M},\mathrm{t}+\mathrm{\Delta }\mathrm{t})={\mathrm{g}}_{\mathrm{i}}(\mathrm{M},\mathrm{t})[1-\frac{1}{{\mathrm{\tau }}_{\mathrm{g}}}]+\frac{{\mathrm{g}}_{\mathrm{i}}^{\mathrm{eq}}(\mathrm{M},\mathrm{t})}{{\mathrm{\tau }}_{\mathrm{g}}}+\mathrm{\Theta }\mathrm{\Delta }\mathrm{t}{\mathrm{F}}_{\mathrm{gi}}$$

(31)

where ${\mathrm{R}}_{\mathrm{c}}$ is the heat capacity ratio (equation (11)), $\mathrm{\epsilon }$ is the porous media porosity and ${\mathrm{g}}_{\mathrm{ieq}}(\mathrm{M},\mathrm{t})$ is the equilibrium distribution function corresponding to ${\mathrm{g}}_{\mathrm{i}}$ ^43,46

$${\mathrm{g}}_{\mathrm{i}}^{\mathrm{eq}}(\mathrm{M},\mathrm{t})={\mathrm{\omega }}_{\mathrm{i}}\mathrm{T}(\mathrm{M},\mathrm{t})[1+\frac{3{\mathrm{e}}_{\mathrm{i}}·\mathrm{u}(\mathrm{M},\mathrm{t})}{(\mathrm{\Theta }(\mathrm{\epsilon }-1)+1)}]$$

(32)

where ${\mathrm{\omega }}_{\mathrm{i}}$ are the weighting factors presented in (equation (27)) and ${\mathrm{\tau }}_{\mathrm{g}}$ is the thermal relaxation time in the whole computational domain

$${\mathrm{\tau }}_{\mathrm{g}}=\mathrm{\Theta }({{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{mp}}-{{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{fld}})+{{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{fld}}$$

(33)

where ${{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{fld}}$ is the relaxation time in the fluid region ⁴³ (equation (34))

$${{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{fld}}=3\mathrm{\alpha }+0.5$$

(34)

where ${{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{mp}}$ is the relaxation time related to the temperature in the porous media ³⁴ (equation (35))

$${{\mathrm{\tau }}_{\mathrm{g}}}_{,\mathrm{mp}}=3\mathrm{Pr}({\mathrm{\tau }}_{\mathrm{f}}-0.5)\mathrm{C}+0.5\mathrm{with}\mathrm{C}=\frac{\mathrm{\epsilon }(1-\mathrm{\epsilon }){\mathrm{R}}_{\mathrm{k}}}{\mathrm{\epsilon }(1-\mathrm{\epsilon }){\mathrm{R}}_{\mathrm{c}}}$$

(35)

where α is the fluid thermal diffusivity.

The streaming step is given as

$${\mathrm{g}}_{\mathrm{i}}(\mathrm{M}+{\mathrm{e}}_{\mathrm{i}}\mathrm{\Delta }\mathrm{t},\mathrm{t}+\mathrm{\Delta }\mathrm{t})={\mathrm{g}}_{\mathrm{i}}(\mathrm{M},\mathrm{t}+\mathrm{\Delta }\mathrm{t})$$

(36)

Validation

The numerical code is developed using Fortran and is called porous media stack for a TAR (PMS-TAR). To validate the proposed model, two computing cases are considered. The former concerns a planar channel with a porous plug in the middle to test the efficiency of the LBM under this geometry. The fluid velocity vector is normal to the porous media interface plan. The simulation parameters are Pr = 0.71, Da = 10⁻² and Re = 1. The presence of porous media inside the channel affects strongly the velocity profile in the centreline (Figure 3(a)). It can be concluded that using the developed code based on the LBM, obtained results are in concordance with previously published papers. ⁴⁷

Figure 3.

Comparison of our results and the published ones ((a) velocity centreline validation, (b) velocity profile validation inside the porous media for oscillating flow and (c) thermal profile validation inside the porous media for oscillating flow at t = T_e/8).

The latter computing case is established to validate the velocity and the temperature profiles into a two-dimensional empty planar channel subject to oscillating flow. ²⁷ The resolution parameters used in this simulation are the Reynolds number Re = 5, the Womersley number α = H/2(ω/ν)^0.5 = 0.937, the angular frequency ω = 0.000785 and the wave amplitude A = 10⁻⁴. The velocity at t = T_e/8 becomes developed and reaches a maximum value along the Y-axis (around 0.5) (Figure 3(b)). For the same position on the Y-axis, the dimensionless temperature falls to a minimum value as shown in Figure 3(c). Once again, obtained numerical results prove the validity of the approach and show the concordance between the PMS-TAR software outcomes and similar research papers.

Results and discussion

In this article, the LBM is applied to determine the dynamic and thermal fields inside the TAR (Figure 1(b)). The two-dimensional closed tube is derived from Tasnim and Fraser’s ³⁶ geometry. The resonator dimensions are (H, 9.64H). The resonator is distributed as follows

$$\{\begin{array}{c}0\le \mathrm{X}<{\mathrm{X}}_1:\mathrm{Part}1\left(\mathrm{fluid}\mathrm{region}\right)\hfill \\ \multicolumn{1}{c}{{\mathrm{X}}_1\le \mathrm{X}<{\mathrm{X}}_2:\mathrm{Part}2\left(\mathrm{porous}\mathrm{media}\mathrm{region}\right)}\\ \multicolumn{1}{c}{{\mathrm{X}}_2\le \mathrm{X}<L:\mathrm{Part}3\left(\mathrm{fluid}\mathrm{region}\right)}\end{array}$$

where X₁ and X₂ are the stack end positions with X₁ = (x_c– l_c/2)/H and X₂ = (x_c+ l_c/2)/H; x_c and l_c are the stack centre position and the stack length, respectively. The oscillating flow is produced by the loudspeaker, which is located on the left side of the resonator (X = 0) and generates an acoustic wave whose velocity is sinusoidal. The emitted wave is reflected by the closed resonator end (x = L). Uniform temperature of 300°K prevails in the resonator and the stack. A frequency of 200 Hz is delivered by the loudspeaker. The Reynolds number and the Prandtl number are set to 25 and 0.71, respectively.

Dynamic field

In order to give a pictorial presentation of the flow inside the thermoacoustic resonator, the velocity field and U contours are recorded every quarter of period during the last simulation cycle. Figures 4(a)–(d), 5(a)–(c) and 6 illustrate the dynamic field in the three parts of the resonator. Figure 4 shows that the shape of the velocity profiles varies in the different parts. It can be depicted that the velocity reaches a maximum at the centreline of the geometry in the fluid parts (parts 1 and 3). Then, the motion of the fluid is accelerated in the first part and decelerated in the second one due to the presence of porous media and the closed end of the resonator. According to the same figure, the presence of annular effect can likewise be noticed at t = T_e/4 and t = 3T_e/4. The Richardson and Tyler ⁴⁸ effect occurs near the top-bottom walls. This effect is as a result of both the wall viscous phenomenon and the fluid flow inertial force that comes back from the last period in the middle of the channel. ²⁷

Figure 4.

Instantaneous U contours in the resonator at four instants: (a) t = T_e/4, (b) t = T_e/2, (c) t = 3T_e/4 and (d) t = T_e.

Figure 5.

Velocity fluctuation in a whole period at (a) X = 4.71, (b) X = 5.21 and (c) X = 7.57 for f = 200 Hz, Pr = 0.71, ε = 0.7, Da = 0.01, R_k = 30, R_c =2.5.

Figure 6.

Periodically variation of U at Y = 0.5, X = 4.71, 5,21, 7.57 for x_c = 5.21, f = 200 Hz, Pr = 0.71, ε = 0.7, Da = 0.01, R_k = 30, R_c = 2.5.

Figures 4 and 5(a) and (c) show that the maximum velocity is located in the channel centreline, within the fluid regions, parts 1 and 3, at t = T_e/2 and T_e. Then it decreases progressively to zero near the channel walls. However, at t = T_e/4 and 3T_e/4, the annular effect takes place. In part 2, the velocity profile is flattened, as shown in Figure 4(b). It is due to the presence of the solid matrix inside the porous media. The profile shows also an important dynamic boundary layer near the horizontal walls. ³⁴

In the three parts of the resonator, a periodic symmetry is observed in the temporal velocity profiles throughout the resonator area with a phase difference of half a period (Figure 5). It is interesting to note that the amplitude of velocity decreases as the X-position becomes higher and gets away from the loudspeaker. In the porous stack, the dynamic boundary approaches the top-bottom walls more than it does in fluid regions.

Figure 6 shows the periodic variations of velocity at three different locations along the centreline of the channel. The velocity profile has the same trend in different instances. This velocity has two local maxima at wt = 0° and wt = 360° and a local minimum at wt = 180°. Therefore, the amplitude of the velocity decreases gradually from the inlet to the resonator exit. The observed decrease can be explained by the presence of the porous media and the closed end of the resonator. ²⁷

Figure 7 demonstrates the effect of Darcy number on the velocity in the middle of the stack. It is clear that the presence of a porous medium affects deeply the flow. The amplitude of the velocity is proportional to the Darcy number as long as it has low values. In this case, the use of material with high permeability is recommended, to enhance the flow in the porous stack.

Figure 7.

Dimensionless velocity profiles at the middle of the stack at different values of the Darcy number for ε = 0.92, R_k = 10, R_c = 2.5 and Ec = 0.1.

Thermal field

The main focus of this section is to study the effect of different stack parameters and to specify the best combination for a TAR having a squared stack made of porous media. The performance of the TAR is evaluated using the average temperature difference between the stack extremities across a period at a steady state. It is calculated as follows (equation (37))

$$\mathrm{\Delta }\mathrm{T}=\frac{1}{{\mathrm{T}}_{\mathrm{w}}}\underset{{\mathrm{T}}_{\mathrm{f}}-{\mathrm{T}}_{\mathrm{e}}}{\overset{{\mathrm{T}}_{\mathrm{f}}}{\int }}({\mathrm{T}}_{\mathrm{m}}\left({\mathrm{n}}_2\right)-{\mathrm{T}}_{\mathrm{m}}\left({\mathrm{n}}_1\right))\mathrm{dt}$$

(37)

where

$${\mathrm{T}}_{\mathrm{m}}\left(\mathrm{x}\right)=\frac{1}{\mathrm{H}}\underset{\mathrm{o}}{\overset{\mathrm{H}}{\int }}\mathrm{T}(\mathrm{x},\mathrm{y})\mathrm{dy}$$

(38)

T_e is the time period of the wave, T_f is the simulation time, T_m(x) is the average of the temperature at the ‘x’ section, H is the width of channel, x = n₁ and x = n₂ are the stack extremities, and T(x, y) is the temperature in the (x, y) position.

In the present section, the effect of the system parameters is studied. Special attention is paid to the viscous dissipation term since it is often neglected in previous researches.

The effect of the Darcy number on the TAR performance is examined and presented in Figure 8. The lower is the value of the Darcy number, the larger the indicator is. This indicator is the performance called the time-averaged temperature difference $\mathrm{\Delta }\mathrm{T}$ , generated across the stack, over a cycle. The $\mathrm{\Delta }\mathrm{T}$ decreases from a high peak corresponding to low values of Darcy number to become quasi-asymptotic for those higher than 0.1. Therefore, the optimization of the TAR needs the use of porous material with low permeability. This type of material induces lower values of pressure drop; thus, it contributes more to the enhancement of heat transfer. In fact, the heat conductivity of the solid matrix is more important than that of the fluid, R_k = 10. As a consequence, the heat transfer inside the porous media is enhanced. ⁴⁹ For these reasons, the TAR performs more with stacks characterized by low values of Darcy number.

Figure 8.

Time-averaged temperature difference between the stacks extremities for different Darcy numbers for R_c = 2.5, R_k = 10, Ec = 0.1, ε = 0.7, x_c = 5.21 and l_s = 1.

Figures 9–12 show the time-averaged temperature difference generated across the stack, over a cycle, $\mathrm{\Delta }\mathrm{T}$ versus Eckert number (Ec) values ranging from 0 to 0.2. In these figures, the effect of Darcy number (Da), heat capacity ratio (R_c), porosity (ε) and heat conductivity ratio (R_k) are successively presented and interpreted. In all studied cases, the $\mathrm{\Delta }\mathrm{T}$ increases linearly with $\mathrm{Ec}$ . The maximum reached by $\mathrm{\Delta }\mathrm{T}$ for large values for $\mathrm{Ec}$ is strongly affected by each considered parameter: Da, R_k, R_c and ε. The heat transfer exchange between the stack extremities is enhanced by higher values of Ec. To expand and clarify, this is due to the contribution of the heat dissipation as a source of heating. These results are in concordance with those published previously.^34,50 Thus, the dissipation term is of great importance in the energy equation (12) and should from now on be considered.

Figure 9.

Dependence of the time-averaged temperature difference between the stacks extremities on the Eckert number for different Darcy numbers for R_c = 2.5, R_k = 10, ε = 0.7, x_c = 5 and l_s = 1.

Figure 10.

Dependence of the time-averaged temperature difference between the stack extremities on the Eckert number for different heat capacity ratio for Da = 0.001, R_k = 10, = 0.7, x_c = 5.33 and l_s = 1.

Figure 11.

Dependence of the time-averaged temperature difference between the stacks extremities on the Eckert number for different stack porosities Da = 0.001, R_k = 10, R_c = 20, x_c = 5 and l_s = 1.

Figure 12.

Dependence of the time-averaged temperature difference between the stacks extremities on the Eckert number for different heat conductivity ratio for Da = 0.001, ε = 0.3, R_c = 20, Ec = 0.2, x_c = 5 and l_s = 1.

Figure 10 shows the dependence of the $\mathrm{\Delta }\mathrm{T}$ on the Eckert number [0, 0.2] for three heat capacity ratios (2.5, 10, 20). The performance indicator increases with R_c. This can be explained by the contribution of the inertial term of the chosen material in the performance of the considered system. This figure shows again that the TAR performance increases with higher values of the heat capacity ratio. For a fixed Eckert number, Ec = 0.2, the $\mathrm{\Delta }\mathrm{T}$ across the thermoacoustic stack increases from 26.37°C to 40.13°C for a heat capacity ratio ranging from 2.5 to 20. These results are similar to previous researchers’ outcomes, studying stacks of parallel plates. Certainly, the efficiency of TAR increases using a stack material with high heat capacity and fluid with a low heat capacity. ⁹

Figure 11 illustrates $\mathrm{\Delta }\mathrm{T}$ across the stack as a function of the Eckert number (from 0 to 0.2) for three different porosities ε (0.3, 0.7 and 0.9). The squared stack dimensions are kept constant. For a fixed Eckert number (Ec = 0.2), as ε is decreased from 0.9 to 0.3, the temperature difference between the stack extremities increases a tiny bit from 38.3°C to 42.4°C (Figure 11). It can be concluded that TAR is more efficient with low stack porosities. The heat conductivity ratio is among the most crucial parameters in the choice of the stack material and the fluid to have a TAR performance. ⁹ For this reason, the effect of the heat conductivity ratio is presented in Figure 12. For large values of R_k, the $\mathrm{\Delta }\mathrm{T}$ decreases to acquire a value of 63.18°C, for Ec = 0.2. Therefore, low heat conductivity ratio values are required to improve TAR performance.

Figure 13 presents the temperature difference across the stack versus the stack centre positions varying from λ/16 to 3λ/16. The existence of a critical value of stack position is noticed, for which the considered system is more efficient. Therefore, the maximum temperature difference observed is 67.5°C at x_c =λ/8. This position coincides with the middle of the resonator and the velocity antinode. ⁵¹

Figure 13.

Time-averaged temperature difference generated across the stack extremities versus stack centre positions for Da = 0.001, R_c = 20, R_k = 1, x_c = 5.33 and l_s = 1.

Conclusion

A numerical investigation is performed to analyse the impacts of some parameters on the efficiency of a solar-powered TAR. The studied TAR is equipped with a rectangular section resonator and a homogeneous porous media stack. The considered parameters include the Darcy number, the stack porosity, the heat capacity ratio, the heat conductivity ratio, the Eckert number and the stack centre position. The Darcy–Brinkmann–Forchheimer extended model is employed for the momentum equation. The proposed approach follows the LBM using the BGK approximation. The time-averaged temperature difference between the stack ends over a cycle is used as an indicator to qualify the TAR performance. The used model consists of solving the Navier–Stokes and the ideal gas equations. The velocity and temperature fields are obtained and analysed. The results obtained in this article are as follows:

• The flow is accelerated in the first part and decelerated in the second one due to the presence of porous media and the closed end of the resonator. The velocity profiles are harmonic and present a phase angle for half a period.

• The amplitude of velocity decreases getting away from the loudspeaker. The dynamic boundary approaches the top-bottom walls in the porous stack; more than it does in fluid regions.

• The flow in the considered system is enormously affected by the Darcy number.

• The use of porous material with low permeability and high conductivity in the TAR enhances its efficiency.

• The viscous dissipation effect, considered internal source, affects the performance of the system. This impact becomes important using a porous material with higher values of the heat capacity ratio.

• The best stack centre location is the middle of the quarter wavelength resonator.