Thermal conductivity of nanofluid in nanochannels

This paper concerns the behaviour of a copper–argon nanofluid confined in a nanochannel. Using molecular dynamics simulations, it is shown that in narrower channels, the thermal conductivity increases by approximately 20 % compared to macroscopic cases. The results suggest that the structured liquid layers surrounding the solid particles occupy a greater percentage of the system in narrower channels, thus enhancing the thermal conductivity of the nanofluid.

of the particle within the nanofluid. However, experimental studies on the heat transfer properties of nanofluids observed thermal conductivities far exceeding those predicted by classical approaches (Eastman et al. 2001;Choi et al. 2001;Patel et al. 2003).
A large body of theoretical work has been devoted to delineate the physical significance of this enhancement. Some studies initially speculated on the importance of the collisions of nanoparticles due to their Brownian motion (Keblinski et al. 2002;Jang and Choi 2004;Evans et al. 2006), a phenomenon not included in the classical models which assume a static system. It was, however, deemed unimportant as it is orders of magnitude slower than the thermal diffusion (Wang et al. 1999;Keblinski et al. 2002;Eapen et al. 2010;Evans et al. 2006). The possibility of a nanoscale convection induced by the Brownian motion of nanoparticles (Jang and Choi 2004) was also investigated with arguments both in favour of (Jang and Choi 2004;Koo and Kleinstreuer 2004;Sankar et al. 2008) and against (Evans et al. 2006;Eapen et al. 2007) its significance. The experimental (Yu et al. 2000) and computational (Li et al. 2008) observations of structured liquid atoms around the solid particles motivated research investigating the possibility of an enhanced thermal conductivity due to the solidlike nature of these nanolayers (Yu and Choi 2003;Li et al. 2010). The counter-argument is that the small thickness (≈1 nm) of these structures cannot possibly produce such large enhancements (Keblinski et al. 2002;Eapen et al. 2010). Finally, a number of studies suggest that the failure of theoretical models to predict experimental observations is ascribed to the assumption that nanofluids have uniformly dispersed particles (Eapen et al. 2010(Eapen et al. , 2007Keblinski et al. 2002). Allowing for the existence of long chains and aggregation of nanoparticles, a set of classical limits can be set (Hashin and Shtrikman 2004;Eapen et al. 2010),

Introduction
The heat generated by integrated circuits that increases in tandem with their constantly improved performance becomes a limiting factor in electronic devices, which may also lead to device failure. Nanofluids, suspensions of solid particles of nanoscale diameter in some fluid, are potential candidates for thermal management. Most classical theoretical heat conduction models Maxwell are based on an effective-medium theory according to which the thermal conductivity of a system is an average of the thermal conductivities of individual constituents (Maxwell 1881;Hamilton and Crosser 1962): where l and p are the thermal conductivities of the liquid and particle, respectively, and φ is the volume fraction which includes the values of all experimental data (Eapen et al. 2010). The evidence is compelling, as no experimental data show an anomalous enhancement in uniformly dispersed nanofluids (Eapen et al. 2010). However, theoretical models exist, which claim that aggregation diminishes the thermal conductivity, as it interferes with rapid motion of single nanoparticles (Xuan et al. 2003).
In addition to the unusual nature of nanofluids, liquid cooling usually manifests itself in micro-and nanochannel heat sinks (Tuckerman and Pease 1981;Kleiner et al. 1995). Past experimental (Doerr et al. 1998;Henderson and Swol 1984;Yu et al. 2000) and computational (Bitsanis et al. 1987;Schoen et al. 1987;Heinbuch and Fischer 1989;Priezjev et al. 2005;Drikakis 2010, 2011;Travis et al. 1997;Sofos et al. 2009a) studies have shown that under such spatial restrictions, fluids form dense, discrete layers close and parallel to the channel walls, which can even occupy the entire width of the channel, as it decreases to only a few atomic distances (Sofos et al. 2009a;Giannakopoulos et al. 2014). Solidliquid interactions on nanoscale cannot be modelled by the Navier-Stokes equations (Travis et al. 1997), as well as using no-slip condition assumptions commonly employed by continuum fluid dynamics Drikakis 2010, 2011;Sofos et al. 2009c;Kasiteropoulou et al. 2011;Li et al. 2010).
Previous studies have shown that confinement can affect the transport properties of liquids and that the thermal properties of nanofluidic systems deviate from the limits of continuum mechanics (Liu et al. 2005;Sofos et al. 2010). Nanoflows lead to quadratic temperature variations close to the centre of the channel, which induce a heat flux even in the absence of a temperature gradient Baranyai et al. 1992). Other studies found that the Kapitza (thermal interface) resistance is strongly correlated to the properties of both the channel walls and liquid (Barrat and Chiaruttini 2003;Kim et al. 2008). Investigations have also shown that the thermal conductivity is anisotropic within a nanochannel (Liu et al. 2005), and depends on the channel width (Liu et al. 2005;Sofos et al. 2009b).
The aim of the paper is to improve the understanding of thermal transport of an argon-copper nanofluid residing in an idealised channel. Using molecular dynamics (MD), the thermal conductivity of the suspension is calculated for different channel widths, and the physical mechanisms of liquid structure around a particle and their effects on the thermal conductivity are presented.

Computational model
The computational model consists of argon atoms confined in a channel of nanometre characteristic dimensions as well as a copper particle suspended in the centre of the system. The walls of the channel are fixed perpendicular to the y direction (parallel to the xz plane). L y , the distance separating the walls (in the y direction), is a parameter of interest to this study and is therefore varied across different cases. Periodic boundary conditions are used along the x and z directions, emulating the perpetual continuation of the channel. In the y direction, the boundaries coincide with the walls, and fixed boundary conditions are used.
Each of the two walls consists of two perfect (111) fcc planes with density ρ * wall = 4ρσ 3 , forming closed-group surfaces and assuring enough density to prevent liquid particles from escaping the system. Figure 1 depicts the MD model used in this investigation.
The interatomic interactions between the argon atoms are modelled by the Lennard-Jones (LJ) potential: where i, j are the labels for two arbitrary particles in the system, r ij is their interatomic distance, and ǫ is the characteristic energy level. For computational efficiency, interatomic interactions beyond a cut-off distance r c = 2.2σ are disregarded. The LJ parameters for the potential between argon particles are ǫ ll = 1.0ǫ and σ ll = 1.0σ, while for the wall-argon interactions, the values ǫ wl = 0.6ǫ and σ wl = 0.75σ are used. The ǫ p and σ p take the following values: For the copper particle, ǫ p = 39.676ǫ and σ p = 0.738σ. For the particle-argon interaction, ǫ pl = 6.3ǫ and σ pl = 0.869σ. Finally, for the wall-particle interactions, the values ǫ wp = 39.0ǫ and σ wp = 0.85σ are employed. The number of argon atoms varies with the channel width, while retaining a constant density ρ * l = 0.84mσ −3 and temperature T * = 0.7ǫk −1 B , thus representing the liquid phase of argon (Sarkar and Selvam 2007;Li et al. 2010). Given the system potential, the equations of motion for the particle i are given by The wall particles are fixed onto their initial lattice sites by spring potentials urging them to return to their equilibrium positions r 0 via a restoring force given by where κ is the wall stiffness, a parameter vital to the realistic representation of the wall. The wall stiffness determines the strength of the bonds between the wall particles, as well as their motion. In the present study, the value κ = 500 ǫ σ −2 is used (Asproulis and Drikakis 2010).
In order to control the temperature of a system, each fcc plane of the walls is assigned a thermostat (Kim et al. 2008). The regulation of the temperature is achieved by explicitly rescaling the per-atom velocities averaging them up to the desired temperature. Both walls are set to the same temperature keeping the system in equilibrium. No thermostats or other artefacts are used on the liquid, portraying a behaviour as realistic as possible, and any excess viscous heat is transferred onto the walls through collisions, where it is further dissipated into the environment via the thermostats.
The two most widely used techniques for calculating the thermal conductivity are the 'direct' method and the Green-Kubo (GK) approach. The 'direct' method poses a onedimensional temperature gradient, and the resulting heat flux in conjunction with the Fourier's law of conduction is used to calculate the thermal conductivity. The use of the 'direct' method, however, is not reliable (McGaughey and Kaviany 2004) because at the nanoscale, the temperature difference in the order of 10 K may result in nonlinear temperature profiles.
The GK formalism follows a different approach for calculating the heat conduction of materials. A system in equilibrium has no net heat flux. However, perturbations in the system cause the heat flux vector to oscillate about the zero value. The thermal conductivity is then linked with the time taking for the energy fluctuations to dissipate. More formally, the thermal conductivity is calculated using where is the thermal conductivity; V is the volume of the system; κ B is the Boltzmann constant; d is the number of dimensions of the system; and the angle brackets indicate an autocorrelation function. The microscopic heat flux J is given by: where the indices i and j denote arbitrary liquid particles; v i is the velocity of atom i; E is its kinetic and potential energy; r ij is the interatomic distance between i and j; F ij denotes the two-body forces between i and j, and h is the microscopic enthalpy.
For a system in equilibrium, the heat flux autocorrelation function (HFACF) �J(0) · J(t)� should eventually decay to zero, so that its integral (and therefore the thermal conductivity) has a finite and well-defined value. The main issue with this approach is that the simulation must be performed for a sufficient number of time steps for the integral of the HFACF to reach a plateau. Moreover, for long autocorrelation timescales, the noise introduced by the HFACF may prevent the convergence of the integral. The calculations of the thermal conductivity in this investigation are carried out using the GK formalism, because it accommodates the small channel width (as small as 10σ).
The simulation time step is �t = 0.001τ ≈ 2 fs. An initial equilibration phase of 2 × 10 6 time steps is performed to allow the temperature and energy of the system to settle. The simulations are then performed for a further 2 × 10 7 time steps and the positions at each time step are taken from the microcanonical ensemble (NVE). For the calculation of HFACF, a correlation length of 2 × 10 3 time steps is used, which gives the autocorrelation function sufficient time to decay.

Results and discussion
The thermal conductivity has been calculated for channels of various channel widths. The radius of the particle remains fixed resulting in a diminishing particle volume fraction with increasing channel width. Figure 2 shows the thermal conductivity as a function of the channel width for the nanofluid, including the copper particle, and for the pure argon case. The results clearly show that the addition of the particle enhances the thermal conductivity for all channel widths used here. The narrower the channel is, the more significant the effect of nanoparticle becomes. This behaviour of the nanofluid is expected: the system size increases, while the size of the particle remains the same resulting in a reduction in the volume fraction of the particle. For wider channels, the thermal conductivity of the nanofluid is reduced, tending to approach its value obtained in the pure argon case (red curve in Fig. 2). With the exception of the first point corresponding to the narrowest channel, the thermal conductivity against the channel width is governed by a linear relationship (black curve in Fig. 2): The increase in the thermal conductivity when increasing the volume fraction is also predicted by the conventional Maxwell model (Eq. 1). In Fig. 3, the thermal conductivity of the present MD calculations for a confined nanofluid (blue) and the predictions of the analytical model (black) as a function of the volume fraction are shown. As the volume fraction approaches zero (i.e. pure argon), both the MD and the Maxwell model converge to ≈ 0.132 W /mK, which is the experimentally observed value for liquid argon at temperature 0.71 k B /ǫ (Müller-Plathe 1997). However, the thermal conductivity of the confined nanofluid increases significantly faster. This is contrary to the linear relationship of the Maxwell model as well as previous MD studies on unconfined nanofluids, which show an initial sudden increase followed by a more gradual enhancement (Sarkar and Selvam 2007).
(7) = −0.0009L y + 0.167 To shed light on these discrepancies, we consider the effect of the particle on the liquid structure. The radial distribution function (RDF) of the liquid around a nanoparticle, for all channel heights, shows how the liquid density is distributed around the solid sphere (Fig. 4) (Sarkar and Selvam 2007).
As the channel width increases, the density of the nanolayer around the particle decreases, with the density differences being larger between the layers 10σ and 17σ. This liquid layering has also been found to enhance the thermal conductivity of unconfined liquids (Yu and Choi 2003;Li et al. 2010). In confined liquids, however, as the channel width decreases, the particle is surrounded by wider and more dense structured layers, which can be the reason of the thermal conductivity enhancement. Furthermore, as the channel width decreases, the ratio of the surface area of the sphere to the liquid volume increases exponentially. Hence, in narrower channels, the nanolayers occupy a significantly larger percentage of the overall system, thus contributing further the thermal conductivity enhancement.
In addition to the nanolayers around the particle, the interaction between the wall and liquid atoms forms structured layers near the walls (Fig. 5). In contrast to the spherical nanolayers, which decrease in density with increasing channel size, the local maxima of the layers close to the wall seem to increase. The two phenomena are not independent. For small channel widths, the influence of the particle extends across the entire channel width, pulling the argon atoms located near the walls towards the centre of the channel, thus reducing the local maxima of the density profiles. However, this does not seem to have a negative effect on the thermal conductivity of the nanofluid. The results suggest that in narrower channels, the presence of the nanoparticle and the formed liquid nanolayers around it have a more significant effect on the thermal conductivity (see Figs. 2, 3).
The anisotropic behaviour of thermal conductivity with respect to the x−, y− (normal to the walls) and z-direction has also been investigated. The thermal conductivity in the y direction increases when increasing the channel width and gradually reaches a plateau (Fig. 6). This is the opposite to the behaviour of the total thermal conductivity (Fig. 2). Therefore, we suggest that the increased thermal conductivity when decreasing the channel width is due to the behaviour of the nanofluid in the parallel directions (x and z). The anisotropy in the y-direction is attributed to the smaller collision frequencies between the argon atoms, near the wall, which results in a less efficient transfer of energy. Past papers have confirmed a similar behaviour in pure liquids (with no nanoparticles) (Murad et al. 1993;Liu et al. 2005). The slight deviation of the results from the plateau for the largest channel width (Fig. 2) is due to fluctuations arising from the Green-Kubo approach.
The anisotropic behaviour in the direction normal to the walls is also reflected on the mean square displacement (MSD) for different channel widths (Fig. 7b). Although the MSD is considered to provide uncertain results for the diffusion coefficient in confined fluids, it qualitatively shows here the anisotropy in the y-direction (Fig. 7b), where the MSD significantly differs from the results in the other directions, gradually reaching a plateau for different channel widths.

Conclusions
The thermal behaviour of an argon-copper nanofluid confined in a nanochannel has been studied by molecular dynamics simulations. The results reveal that in addition to the increasing volume fraction, confinement can also enhance the thermal conductivity of the nanofluid. Furthermore, an investigation of the liquid distribution shows that this is due to the organised liquid structures forming around the nanoparticle. In narrow channels, the density and thickness of these nanolayers increases significantly. Furthermore, the high ratio of the surface area of the nanoparticle to the liquid volume at narrower channels implies that the structured liquid layers occupy a greater percentage of the overall system. We believe that this is responsible for the enhancement of the thermal conductivity. Finally, the anisotropy of the nanofluid was investigated showing that the thermal conductivity in the direction normal to the walls increases when increasing the channel width. This behaviour deviates from that of the total thermal conductivity, and it is a result of the impaired atomic motion under severe spatial restriction, which also leads to less efficient energy transfer.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.