Theoretical and Numerical Study on Buongiorno’s Model with a Couette Flow of a Nanofluid in a Channel with an Embedded Cavity

Abstract

In the present paper, the fluid flow and heat transfer of a nanofluid are numerically investigated. More specifically, reference is made to a nanofluid, described by means of Buongiorno’s model, subjected to Couette flow. The considered domain consists of a channel that displays a cavity shortly after the inlet section. The transport model for the nanofluid, that is the mass conservation, momentum, and nanoparticles equation, is written in a dimensionless form and solved by employing the software package Comsol Multiphysics. Many ideas emerged from this work: the visualization of the velocity stream function, the dimensionless temperature, and nanoparticle concentration fields are provided, as a function of the governing parameters: Reynolds, Peclet, Lewis, Brownian diffusivity number, and thermophoretic diffusivity number. Concerning the nanofluid typical effects, the thermophoretic diffusion seems to affect the solution much more than the Brownian diffusion. The Nusselt number on the upper wall is calculated as well, and the results show that it proves to be, in most of the considered cases, an increasing function of the Reynolds number. Moreover, concerning the Nusselt number, the Brownian diffusion effects are shown to be negligible.

1. Introduction

The term “nanofluid” was first coined by Choi and Eastman [1,2] in 1995. Originally, the term described nanometer-sized copper particles dispersed in water to improve their thermal conductivity. Today, a nanofluid can be more accurately intended as an engineered colloidal suspension of nanoparticles into a base fluid. Nanofluids have shown a promising future as thermal fluids for various heat transfer applications [3]: in solar collectors [4,5], to enhance the heat-transfer coefficient of solar water heaters or to improve the capacity of thermal energy storage systems; and in the area of refrigeration, to enhance the performance of refrigeration systems. Later, many authors (see for example refs. [6,7,8,9,10]) studied the nanofluid flow and discussed their heat transfer performance. In the literature, a widely employed approach to describe the features of a nanofluid is the so-called Buongiorno model, which relies on the assumptions: (i) the flow is incompressible, (ii) there are no chemical reactions at nanofluids, (iii) external forces are negligible, (iv) the mixture is considered to be diluted, (v) viscous dissipation is negligible, (vi) radiation heat transfer is negligible, (vii) nanoparticles and base fluid are located in thermal equilibrium.

More specifically, Buongiorno presented the following seven mechanisms of slip between nanoparticle and fluid [11]: Inertia, Brownian diffusion, Thermophoresis, Diffusiophoresis, Magnus effect, Fluid Drainage, and Gravity settling. It is worth mentioning that the above-mentioned phenomena can produce a relative velocity between the base fluid and nanoparticles, that is, water. He also suggested a two-component four-equation nonhomogeneous equilibrium model for three basic equations: mass, momentum, and heat transfer. Furthermore, a comparison of length scales and the nanoparticle turbulent eddy time demonstrated that, in the presence of turbulent eddies, the nanoparticle moves with the base fluid homogeneously. As a consequence, the turbulence intensity effect is also double. The presence of turbulence in the flow, which prevents the agglomeration of particles, decreases the discrepancy between experimental and numerical results. Hence, the results obtained from an analysis of Buongiorno’s two-phase model are likely to be more accurate.

In this context, Azimikivi H. et al. [12] recount the effects of the rib shape with constant height, width and pitch, and nanoparticle diameter on natural convection heat transfer of low turbulence flow of AL2O3-water nanofluid inside a square cavity based on Buongiorno’s model. The steady Cu-water nanofluid flow in the presence of a magnetic field has been numerically investigated by Rawat S.K. et al. [13] under the effects of mixed convection, thermal radiation, and chemical reaction. Besides this, the modeled nanofluid problem is characterized by a single-phase nanofluid model (i.e., Tiwari and Das’s model) and a two-phase mixture nanofluid model (i.e., Buongiorno’s model) together with the Cattaneo–Christov Double diffusion model.

In addition, according to Buongiorno’s model, many papers have examined the convection of a nanofluid in a cavity or in a channel with an embedded cavity. In this regard, the forced convection heat transfer of a nanofluid in a channel has been widely investigated in refs. [14,15], while the mixed convection of a nanofluid (with various types of nanoparticles) in an open cavity has been pointed out by Mehrez et al. in ref. [16]. The single phase model has been employed to investigate the Couette forced convection of a nanofluid in ref. [17]: the governing PDEs are reduced to ODEs using a similarity transformation, and then analytically solved. The recent literature has paid attention, especially, to mixed and free convention, or, concerning forced convection, to the presence of MHD effects, as in refs. [18,19,20]. On the contrary, recently, the analysis of vortex formation in the case of laminar flow of a Newtonian fluid, flowing in a duct with an embedded cavity has been studied in ref. [21], where a brief analysis of the recent literature on the vortex formation in cavities is reported as well. In this context, the motivation in the current work is that flow of nanofluid over a channel/embedded cavity has been extensively discussed in the past by many authors, yet none of them made an attempt to combine the forced convection of a two-phase mixture of nanofluid with Couette flow.

More in detail, in the present paper, attention is paid to the Couette laminar flow of a nanofluid in a channel that encapsulates an embedded cavity. The nanofluid is described by employing the two-phase mixture model introduced by Buongiorno, and the presented equations are written in a dimensionless form and solved numerically by employing the commercial software package Comsol Multiphysics, based on a Galerkin finite element method. The obtained solution is shown to depend on the Reynolds number, the Peclet number, the Lewis number, the Brownian motion coefficient, and the thermophoresis coefficient. The dependence of the solution and of the Nusselt number on the parameters is then investigated and discussed.

2. Mathematical Model

Let us refer to the geometry sketched in Figure 1, that is, a channel with an embedded open cavity. The characteristic lengths of the geometry under consideration are evidenced in the figure (the geometry was built according to the size of the inlet channel H).

In accordance with Buongiorno’s model [11], let us consider the nanofluid as a two-component mixture (base fluid plus nanoparticles) with the flowing assumption: (i) incompressible flow, (ii) no chemical reactions, (iii) negligible external forces, (iv) dilute the mixture, (v) negligible viscous dissipation, (vi) negligible radiative heat transfer, (vii) nanoparticles and base fluid locally in thermal equilibrium. Definitely, in accordance also with Hussain and Ahmed [18], the transport model for the nanofluid is described by the following system of equations:

• Continuity equation

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

  (1)
• Momentum equations in X and Y direction

  $$\rho \left(U\frac{\partial U}{\partial X}+V\frac{\partial U}{\partial Y}\right)=-\frac{\partial P}{\partial X}+\mu \left(\frac{{\partial }^{2}U}{\partial X^2}+\frac{{\partial }^{2}U}{\partial Y^2}\right)$$

  (2)

  $$\rho \left(U\frac{\partial V}{\partial X}+V\frac{\partial V}{\partial Y}\right)=-\frac{\partial P}{\partial Y}+\mu \left(\frac{{\partial }^{2}V}{\partial X^2}+\frac{{\partial }^{2}V}{\partial Y^2}\right)$$

  (3)
• Energy equation

  $$U\frac{\partial T}{\partial X}+V\frac{\partial T}{\partial Y}=\alpha \left(\frac{{\partial }^{2}T}{\partial X^2}+\frac{{\partial }^{2}T}{\partial Y^2}\right)+\delta \left\{D_B\left(\frac{\partial \varphi }{\partial X}\frac{\partial T}{\partial X}+\frac{\partial \varphi }{\partial Y}\frac{\partial T}{\partial Y}\right)+\frac{D_T}{T_c}\left[{\left(\frac{\partial T}{\partial X}\right)}^2+{\left(\frac{\partial T}{\partial Y}\right)}^2\right]\right\}$$

  (4)
• Nanoparticles equation

  $$U\frac{\partial \varphi }{\partial X}+V\frac{\partial \varphi }{\partial Y}=D_B\left(\frac{{\partial }^2\varphi }{\partial X^2}+\frac{{\partial }^2\varphi }{\partial Y^2}\right)+\frac{D_T}{T_c}\left(\frac{{\partial }^{2}T}{\partial X^2}+\frac{{\partial }^{2}T}{\partial Y^2}\right)$$

  (5)

In Equations (1)–(5), (U, V) are the X and Y-component of the nanofluid velocity vector, P is the nanofluid pressure, T is the nanofluid temperature, $\varphi$ is the nanoparticles concentration. These dependent variables are given as a function of the position vector (X, Y). Furthermore, ρ is the nanofluid density, µ is the dynamic viscosity, T_c is the cold temperature, α = k/(ρ C_p) is the thermal diffusivity (k is the thermal conductivity of the fluid and C_p is the fluid-specific heat constant), δ = (ρ C_p)_p/(ρ C_p)_f (p and f denotes nanoparticle and the base fluid). Let us introduce the Brownian motion coefficient D_B, defined by the Einstein–Stokes equation, and the thermophoresis coefficient D_T, that is:

$$D_B=\frac{k_B T}{3 \pi \mu d_p}$$

(6)

$$D_T=0.26\frac{k}{2k+k_p}\frac{\mu }{\rho }\varphi$$

(7)

where k_B is Boltzmann’s constant equal to 1.38064852 × 10⁻²³ (m² kg)/(s² K), d_p is the nanoparticles diameter, k_p is the particle thermal conductivity. For the thermal energy Equation (4), non-linearities take place with both Brownian diffusion and thermophoretic diffusion.

The boundary conditions are prescribed as follows:

• Inlet

  $$U=U_0\frac{Y}{H}, V=0, T=\left(T_h-T_c\right)\frac{Y}{H}+T_c,\varphi ={\varphi }_0$$

  (8)
• Outlet

  $$P=0, \frac{\partial T}{\partial X}=0$$

  (9)
• Hot wall

  $$U=U_0, V=0, T=T_h, D_B\frac{\partial \varphi }{\partial Y}+\frac{D_T}{T_c} \frac{\partial T}{\partial Y}=0$$

  (10)
• Horizontal cold walls

  $$U=0, V=0, T=T_c, D_B\frac{\partial \varphi }{\partial Y}+\frac{D_T}{T_c} \frac{\partial T}{\partial Y}=0$$

  (11)
• Vertical cold walls

  $$U=0, v=0, T=T_c, D_B\frac{\partial \varphi }{\partial X}+\frac{D_T}{T_c} \frac{\partial T}{\partial X}=0$$

  (12)

In Equations (8)–(12), U₀ is the slip velocity of the hot wall and T_h is the temperature of the hot wall.

In detail, we investigate Couette flow: thus the boundary condition on the velocity field is the tangential motion of the upper surface with respect to the lower boundary. The effect of buoyancy is neglected. Moreover, the boundary conditions on the temperature foresee a uniform hot temperature T_h on the upper wall, a uniform cold temperature T_c on the lower wall and on the cavity walls, and a linearly varying temperature distribution both at the inlet and at the outlet sections. The hot and cold walls are marked with a red and blue line, respectively, in Figure 1.

On account of the dimensionless quantities

$$\left(x,y\right)=\frac{\left(X,Y\right)}{H}, \left(u,v\right)=\frac{\left(U,V\right)}{U_0}, p=\frac{P}{\rho U_0^2}, \mathsf{\Theta }=\frac{T-T_c}{T_h-T_c}, \mathsf{\Phi }=\frac{\varphi }{{\varphi }_0}$$

(13)

the system of Equations (1)–(5) and Equations (8)–(12) can be rewritten in the dimensionless form as follows:

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

(14)

$$u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}=-\frac{\partial p}{\partial x}+\frac{1}{\mathrm{Re}}\left(\frac{{\partial }^{2}u}{\partial x^2}+\frac{{\partial }^{2}u}{\partial y^2}\right)$$

(15)

$$u\frac{\partial v}{\partial x}+V\frac{\partial v}{\partial y}=-\frac{\partial p}{\partial y}+\frac{1}{\mathrm{Re}}\left(\frac{{\partial }^{2}v}{\partial x^2}+\frac{{\partial }^{2}v}{\partial y^2}\right)$$

(16)

$$u\frac{\partial \mathsf{\Theta }}{\partial x}+v\frac{\partial \mathsf{\Theta }}{\partial y}=\frac{1}{\mathrm{Re} \mathrm{Pr}}\left\{\left(\frac{{\partial }^2\mathsf{\Theta }}{\partial x^2}+\frac{{\partial }^2\mathsf{\Theta }}{\partial y^2}\right)+N_B\left(\frac{\partial \mathsf{\Phi }}{\partial x}\frac{\partial \mathsf{\Theta }}{\partial x}+\frac{\partial \mathsf{\Phi }}{\partial y}\frac{\partial \mathsf{\Theta }}{\partial y}\right)+N_T\left[{\left(\frac{\partial \mathsf{\Theta }}{\partial x}\right)}^2+{\left(\frac{\partial \mathsf{\Theta }}{\partial y}\right)}^2\right]\right\}$$

(17)

$$u\frac{\partial \mathsf{\Phi }}{\partial x}+v\frac{\partial \mathsf{\Phi }}{\partial y}=\frac{1}{\mathrm{Re}\mathrm{Pr}\mathrm{Le}}\left[\left(\frac{{\partial }^2\mathsf{\Phi }}{\partial x^2}+\frac{{\partial }^2\mathsf{\Phi }}{\partial y^2}\right)+\frac{N_T}{N_B}\left(\frac{{\partial }^2\mathsf{\Theta }}{\partial x^2}+\frac{{\partial }^2\mathsf{\Theta }}{\partial y^2}\right)\right]$$

(18)

$$\mathit{Inlet}u=y, v=0, \mathsf{\Theta }=y, \mathsf{\Phi }=1$$

(19)

$$\mathit{Outlet}p=0, \frac{\partial \mathsf{\Theta }}{\partial x}=0, \frac{\partial \mathsf{\Phi }}{\partial x}=0$$

(20)

$$\mathit{Hot} \mathit{wall}u=1, v=0, \mathsf{\Theta }=1, N_B\frac{\partial \mathsf{\Phi }}{\partial y}+N_T\frac{\partial \mathsf{\Theta }}{\partial y}=0$$

(21)

$$\mathit{Horizontal} \mathit{cold} \mathit{wall}u=0, v=0, \mathsf{\Theta }=0, N_B\frac{\partial \mathsf{\Phi }}{\partial y}+N_T\frac{\partial \mathsf{\Theta }}{\partial y}=0$$

(22)

$$\mathit{Vertical} \mathit{cold} \mathit{wall} u=0, v=0, \mathsf{\Theta }=0, N_B\frac{\partial \mathsf{\Phi }}{\partial x}+N_T\frac{\partial \mathsf{\Theta }}{\partial x}=0$$

(23)

The dimensionless boundary conditions for the velocity, temperature, and concentration fields are reported in Figure 1. In Equations (14)–(23) the following dimensionless groups have been used:

$$\mathit{Reynolds} \mathit{number}\mathrm{Re}=\frac{\rho u_0 H}{\mu }$$

(24)

$$\mathit{Prandtl} \mathit{number}\mathrm{Pr}=\frac{\mu C_P}{k}$$

(25)

$$\mathit{Lewis} \mathit{number}\mathrm{Le}=\frac{k}{\rho C_P D_B}$$

(26)

$$\mathit{Brownian} \mathit{diffusivity} \mathit{coefficient}{N}_B=\frac{\delta D_B \rho C_P {\varphi }_0}{k}$$

(27)

$$\mathit{Thermophoretic} \mathit{diffusivity} \mathit{coefficient}{N}_T=\frac{\delta D_T \rho C_P \mathsf{\Delta }T}{k T_c},$$

(28)

where ΔT = T_h − T_c.

3. Mesh and Validation

First, let us introduce, for a fixed geometry configuration (Figure 1), a sensitivity analysis of the model described through Equations (14)–(23). The constants range values used for investigating the thermophoresis and Brownian diffusion are: Pr = 6.2, 50 ≤ Re ≤ 500, 0.1 ≤ Le ≤ 5, 0.1 ≤ N_B ≤ 1, and 0.1 ≤ N_T ≤ 1.

The mesh used for the domain discretization is made up of quadrilateral elements. In Figure 2 some of the considered meshes are shown. In particular, in order to identify the optimal number of computational elements, the following quantities are evaluated as average in the domain of interest: velocity, temperature, particle concentration, and Nusselt number on the hot wall.

For the sake of clarity, the mean volume—related to a unit depth—quantities were evaluated according to the expression

$$\overline{\xi }=\frac{1}{\mathsf{\Omega }} {\displaystyle \iint }\xi d\mathsf{\Omega }$$

(29)

where ξ is the general non-dimensional quantities and $\mathsf{\Omega }$ is the non-dimensional volume under consideration. Otherwise, the mean Nusselt number is evaluated along the hot wall by the expression

$$\overline{\mathrm{Nu}}=\frac{1}{l} {{\displaystyle \int }}_0^l\frac{d\mathsf{\Theta }}{dx} dx$$

(30)

where l = l₁ + l₂ + l₃ is the total length of the hot wall.

In

Table 1. Grid independence test. N is the mesh elements number, $\overline{u}$ is the mean velocity, $\mathsf{\Theta }$ is the mean temperature, $\mathsf{\Phi }$ is the mean concentration and $\overline{\mathrm{Nu}}$ is the mean Nusselt number. The mesh validation is done for Re = 100, Pr = 6.2, NB = NT = 0.1 and Le = 10.

#	N	$\overline{\mathit{u}}$	$\overline{\mathsf{\Theta }}$	$\overline{\mathsf{\Phi }}$	$\overline{\mathbf{Nu}}$
1	1246	0.60210	0.43391	1.04165	0.98558
2	30,449	0.60199	0.43461	1.03690	0.98569
3	59,520	0.60198	0.43463	1.03601	0.98571
4	77,826	0.60198	0.43463	1.03600	0.98572
5	120,360	0.60198	0.43464	1.03599	0.98573

the grid independence tests are reported. In the following, grid #3 will be employed.

4. Results and Discussion

Let us now discuss the main features of the solutions of Equations (14)–(23). In Figure 3, Figure 4, Figure 5 and Figure 6, the velocity field stream function, the isotherms, and iso-concentration lines are reported for Pr = 6.2 and by considering the variation of one parameter (Re, Le, N_B, or N_T) while the other ones, respectively, are (100, 0, 0.1, 0.1).

Figure 3 refers to different values assumed by the Reynolds number and shows that this parameter does not seem to affect sensibly the solution. Moreover, the figure shows that for increasing values of Re the stream function and isotherm gradients occur only in the channel region, not involving the open cavity, while the concentration gradients appear mainly in the cavity region.

Figure 4 refers to different values of the Lewis number and shows that for increasing values of Le the concentration vortex tends to occupy the whole cavity. Also, with reference to the channel region, for increasing values of the Lewis number, the concentration gradients tend to move close to the channel walls.

Figure 5 refers to different values of the Brownian diffusivity coefficient, showing that the solution is not strongly affected by this parameter: in fact, slight differences appear not only in the stream function distribution but also in the temperature and in the concentration distribution. Finally, Figure 6 refers to different values of the thermophoretic diffusivity coefficient N_T, showing that also this parameter does not strongly affect the obtained solution.

In Figure 7, the dimensionless concentration (lefthand frame) and the dimensionless temperature (righthand frame) are reported versus the y- coordinate at x = (l₁ + l₂/2)/l, and for Pr = 6.2, Le = 1, N_B = N_T = 0.1 and for different values assumed by the Reynolds number. The considered X-coordinate represents the vertical middle section of the cavity. The figure shows that while for Re = 10 the functions are monotonic, for higher values local minimum and maximum arise. Moreover, the concentration assumes value 1 for any considered Re at y = 0.75. This happens because the velocity value at that point is the same for all the conditions considered. Within the cavity, strong variations of the concentration occur. The dimensionless temperature distribution varies with the vertical coordinate only in the bottom cavity region, and for every mall value of the Reynolds number is a monotonic function.

In Figure 8, the y-velocity is reported versus y at the vertical middle section of the cavity, for Pr = 6.2, Le = 1, N_B = N_T = 0.1, and for different values assumed by the parameter Re. The figure shows that the formation of a clockwise rotating cell in the cavity occurs and that this cell occupies a bigger region of the cavity for increasing values of Re. Moreover, the maximum of the dimensionless velocity increases for increasing values of Re.

In Figure 9 the dimensionless concentration and temperature distribution are reported as a function of the Lewis number. The righthand frame of the figure shows that the dimensionless temperature does not significantly vary: in fact, very small variations occur and are visible in the evidenced region, where the dimensionless temperature is an increasing function of Le. Moreover, a local maximum in the region close to the center of the cavity arises. On the contrary, the concentration strongly depends on the value assumed by the Lewis number. As discussed concerning Figure 7, the boundary condition imposes a concentration equal to 1 for any Le at the upper corner of the cavity.

In Figure 10, the dependence of the concentration and temperature fields on the Brownian diffusivity coefficient N_B are reported, showing that smaller values of N_B yield higher concentrations and to less uniform behavior of this function. On the contrary, N_B does not affect the trend of the dimensionless temperature.

Similar considerations can be done with respect to the parameter N_T, as shown in Figure 11.

Finally, in Figure 12 and Figure 13 the averaged Nusselt number on the hot wall is reported versus Re. Figure 12 analyzes the behavior with respect to the Brownian diffusivity coefficient (left) and to the thermophoretic diffusivity coefficient (right) and Figure 13 to the Lewis number. Figure 12 shows that the parameter N_B does not affect the behavior of the Nusselt number, while N_T and Le do, as evident in the right frame of Figure 12 and in Figure 13. In particular, the average Nusselt number is a decreasing function on N_T, while it is an increasing function of Le. It is interesting to notice that for high values of N_T, that is, N_T = 0.5, the function is not monotonic, thus displaying a minimum for Re = 100.

5. Conclusions

In the present paper, the forced convection of a nanofluid in a channel with an embedded cavity is numerically investigated. Couette flow is considered and the Buongiorno’s two-phase mixture model is employed for the description of the nanofluid. The dimensionless velocity, temperature, and concentrations field are solved by employing the software package Comsol Multiphysics, as a function of the parameters Re, Pr, Le, N_B, and N_T and discussed. In particular, for increasing values of Re, the stream function and isotherm gradients occur only in the channel region, not involving the open cavity, while the concentration gradients appear mainly in the cavity region. For increasing values of the Lewis number, the concentration vortex tends to occupy the whole cavity. Moreover, with reference to the channel region, for increasing values of Le the concentration gradients tend to move close to the channel walls.

Concerning the nanofluid typical effects, the thermophoretic diffusion seems to affect the solution much more than the Brownian diffusion. Finally, the average Nusselt number on the hot wall is calculated as a function of the dimensionless parameters and discussed, showing that the Brownian diffusion effects are negligible.