The term nanofluid was used by Choi  for the first time. After that many researchers continued his works and focused on the modeling of the thermal conductivity of nanofluid    . Recently, the concentration is on the heat transfer and fluid flow behavior of nanofluid.
Most experimental studies for nanofluid are done on macro and micro-scales     . Wen and Ding  studied the heat transfer of Al2O3 -water nanofluid in a copper tube under the constant heat flux. Their measurements in the heat transfer of Al2O3 -water nanofluid showed an enhancement in the entrance region of the tube. They explained that the particle migration reduces the thickness of the thermal boundary layer and causes this behavior in the heat transfer of the nanofluid. Heris et al.  investigated CuO-water and Al2O3 -wa- ter nanofluids in an annular tube. They compared the results of the experimental study and the homogeneous model and realized the under-estimation of the homogeneous model in the heat transfer enhancement, especially in higher nanoparticle volume concentration. Homogeneous (single phase) and two-phase models are commonly used in the numerical study of the heat transfer and fluid flow of the nanofluid. In single phase model, the velocity and temperature of the base fluid and the particles are the same. Most studies in this field are done by using the single phase model    . In the two-phase model, the base fluid and particles are considered as two different phases in despite of the single phase model and both phases have different velocities and temperatures. Also, the interactions between the phases are so important in the governing equations of the two-phase models.
Behzadmehr et al.  studied the turbulent convection of the nanofluid in a circular tube by using a two phase mixture model. After comparing their results with an experimental study, they reported that the mixture model has more coincidence to the experimental study compared to the homogeneous model. Mirmasoumi and Behzadmehr  also used the two phase mixture model for studying the mixed convection of the nanofluid in a tube. Mirmasoumi and Behzadmehr  and Akbarinia and Laur  studied the effect of the nanoparticle size on the mixed convection of the nanofluid by the mixture model. An increase in the heat transfer of the nanofluid with decrease in the nanoparticle size was observed in both studies. Kurowski et al.  simulated the nanofluid flow by three different mixtures, homogeneous and Eulerian-Lagrangian models in a minichannel. The results of all models were almost the same. Fard et al.  used single and two-phase models for studying the heat transfer of the nanofluid inside a tube. They investigated a 0.2% CuO-water and compared the results with experimental ones and reported that the average relative error between the two-phase model and experimental data was 8% while it was 16% for the single phase model.
Mohammad Kalteh et al.  investigated the nanofluid forced convection experimentally and numerically inside a wide microchannel heat sink. They used two-phase Eulerian-Eulerian model in numerical method and reported that the average Nusselt number increases with increase in Reynolds number and nanoparticle volume concentration. Keshavarz and Mohammadi  investigated the thermal performance of Al2O3 -water nanofluid in minichannel heat sink by using single and two-phase models. In their study, it has been observed that the two-phase models have more coincidence with the experimental results in comparison to the single phase model, however, it is not sensible in low volume concentration but by increasing the particle volume concentration up to 1% or in high Reynolds number, the deviation between the single phase model and the experimental data increases. Shariat et al.  studied the impact of nanoparticle mean diameter and the buoyancy force on laminar mixed convection nanofluid flow in an elliptic duct by using two-phase mixture model. They reported that in certain Reynolds and Richardson numbers, the average Nusselt number increases with decrease in the nanoparticle diameter as well as with increase in Richardson number.
In this paper, the nanofluid flow and laminar forced convection in a sudden expansion microchannel with different expansion ratio and isotherm walls are studied by using the Eulerian-Eulerian two-phase model. The governing equations of mass, momentum and energy are discretized using the finite volume method and the modified SIMPLE algorithm is used for solving the flow equations of both phases. Also, the effects of the nanoparticle size and volume concentration and Reynolds number on the average Nusselt number are investigated. Considering the importance of the microfluidic systems and the field of nanofluid heat transfer, the investigation of the flow and heat transfer of nanofluid in a microfluidic device like a sudden expansion microchannel can be so important and practical. In expansion microchannels, the effective hydraulic diameter decreases because of the vortex creation that increases the resistance of the fluid flow. In fact, the sudden expansion microchannels can be used as a microfluidic rectifier, like: Tesla rectifier, simple nozzle/diffuser structures and cascaded nozzle/diffuser structures. These rectifiers are valveless and do the rectification with different flow resistance in the forward flow and the reverse one. On the other hand, two-phase models are proper substitution for a single phase model. Among two-phase models, the Eulerian-Eulerian model has good efficiency because of considering the relative velocity and temperature between phases and the nanoparticle volume concentration distribution.
Numerical investigation of nanofluid flow and heat transfer in a sudden expansion microchannel by using the Eulerian-Eulerian two-phase model and applying the modified SIMPLE algorithm is a new work in the field of nanofluid heat transfer in a microfluidic system.
2. Governing Equations
In this study, the geometry is the sudden expansion microchannel with different expansion ratios. Figure 1 is a schematic of the sudden expansion microchannel. The laminar flow of the mixture of water and copper nanoparticles, enters the microchannel with uniform velocity and temperature and the heat transfer
Figure 1. The sudden expansion microchannel.
occurs between the nanofluid and the isothermal walls of the microchannel. The upstream length and height of the microchannel, the downstream length and height of microchannel and the reattachment length after the expansion region ( ) are shown in Figure 1.
Considering the Eulerian-Eulerian two-phase model and the steady state, laminar and two-dimensional flow of nanofluid, the governing equations for the base liquid and nanoparticle phases can be written as follows.
2.1. Continuity Equations
Continuity equations for solid and liquid phases in the cartesian coordinate system are as follows   :
where , ρ and φ are velocity vector, density and volume concentration, respectively. Subscripts l and p show the base liquid and nanoparticle phases, respectively. Also, for volume concentration of the phases we have
2.2. Momentum Equations
Momentum equations for solid and liquid phases in cartesian coordinate system are as follows:
where , , , P and µ are the particle-particle interaction, drag, virtual mass (added mass) forces, pressure and viscosity, respectively. Here, the lift force between the phases can be neglected because of the small size of the nanoparticles and the gravitational force is neglected because of the small size of the microchannel.
The force terms in the momentum equations are defined as
The friction coefficient β, drag coefficient, the particle-particle interaction modulus G and the particle Reynolds namber are
Here, is the nanoparticle diameter.
Equation (9) is valid for . Also, it should be considered that Equations (8)-(14) are not obtained for nano-sized particles. But, they can be used for nanoparticles because of the lack of the proper correlations for nano-sized particles  .
2.3. Energy Equations
The base fluid and particle phases are considered as incompressible fluids and the viscous dissipation and radiation are neglected. So, the energy equations in cartesian coordinate system can be written as
Here, , , and are the temperature, effective thermal conductivity, heat capacity at constant pressure and volume interphase heat transfer coefficient, respectively. For spherical nanoparticles can be calculated as
Here, is the fluid-particle heat transfer coefficient that is obtained from experimental correlations.
where Pr is the base liquid Prandtl number.
The effective thermal conductivities for both phases are
For spherical nanoparticles
Equations (18)-(25) are not obtained for nano-sized particles. But these correlations are used because of the lack of the appropriate correlations for nano- sized particles  .
The Nusselt number is calculated from the temperature difference between the nanofluid mean temperature and the microchannel walls:
Here, , and are the wall convective heat transfer flux, convective heat transfer coefficient and the microchannel hydraulic diameter, respectively. The mean temperature of the nanofluid is calculated as follow
where the integration is applied on the microchannel cross section.
According to the Equations (15) and (16), the wall convective heat transfer flux can be obtained as
Considering the local heat transfer coefficient, the average heat transfer coefficient is
The non-dimensional parameters using in converting the governing equations are as follows:
where i = l, p indicates the liquid and particle phases.
2.5. Boundary Conditions
At the inlet, the mixture of water and copper nanoparticles enters the sudden expansion microchannel with the same uniform axial velocity. At the outlet, the velocity boundary condition is considered for the outflow of the both phases. In this study, the non-slip boundary condition at the walls is assumed for both phases.
For thermal boundary conditions, the microchannel walls are isothermal and the outflow temperature boundary condition is considered for both phases. Because the length of the microchannel is too long, the variations of the outflow velocity and temperature are assumed as follow
3. Numerical Method
In this study, the non-dimensional forms of the governing equations are discretized by the finite volume method. The first order upwind scheme is used for discretizing the convection-diffusion term. After discretization, the governing equations are converted to a set of algebric equations that is solved iteratively and the SIMPLE algorithm is used for the pressure-velocity coupling   . For using this algorithm, the pressure correction equation is obtained from the combination of the continuity equations of the base fluid and particle phases. In fact, combining the continuity equations for both phases causes a new pressure correction equation that modifies the SIMPLE algorithm and the under-relaxa- tion coefficients are used for velocity and pressure in order to accelerate the convergence of the algorithm. In this algorithm, the source term of the pressure correction equation is considered as the convergence criterion. In this study, the iteration is continued until the convergence criterion becomes smaller than 10−6 in all cells.
4. Grid Independence Study
To investigate that the results are independent from the number of grid points, the average Nusselt number is calculated for various numbers of grid points and Reynolds numbers. Table 1 and Figure 2 show the average Nusselt numbers for four different numbers of grid points and Re = 100. It can be seen that the average Nusselt numbers between the 3 and 4 grids have a negligible difference, so the grid 3 is selected for this study.
Figure 2. The average Nusselt number versus number of grid points.
Table 1. Grid-independency study results for, and.
Because of the lack of the experimental studies for nanofluid flow in a sudden expansion microchannels, the reattachment length for the pure water () at different Reynolds numbers is calculated and compared to Scott et al.  results in order to investigate the accuracy of the code. Table 2 and Figure 3 show the reasonable coincidence between the different amounts of the reattachment length of these two studies for ER = 3 and the maximum deviation is less than 4%.
6. Results and Discussion
6.1. The Viscosity of the Solid Phase
One of the important parameter in this two phase study is the solid viscosity ( ). In fact, the solid viscosity is the virtual viscosity that must be defined correctly in this model. There is not much experimental data for this viscosity in a solid-liquid two phase mixture, so the trial and error method is used for obtaining the proper value for the solid viscosity. This viscosity appears in the Reynolds number of the solid phase, thus is calculated for the different Reynolds number of the solid phase until the velocity profile of this numerical solu-
Figure 3. Comparison between reattachment length of this study and Scott’s one for different Reynolds number.
Table 2. Comparison between reattachment length of this study and Scott’s one.
tion has an appropriate coincidence with it’s analytical solution profile. Consequently, is obtained from the trial and error method and the Reynolds number of the solid phase. Table 3 shows the sensitivity of the Nusselt number to the solid viscosity for, , and.
It can be seen that the changes of the Nusselt number is small when changes. So, the amount of the solid viscosity is not so effective on the results and it is not necessary to find the exact amount of this viscosity. Also, Saffaraval et al.  reported the same result about the particle viscosity. The amount of the solid viscosity is considered to be 0.089 Pa.s in this study.
6.2. Force Terms in the Momentum Equations
The momentum equation has three interphase forces including virtual mass, particle-particle interaction and drag forces. In Table 4, the average Nusselt number is investigated for, and different conditions of the force terms. It can be seen that the particle-particle interaction and vitual mass forces are not effective on the average Nusselt number, however, the drag force changes the average Nusselt number slightly. Neglecting the drag term causes the average Nusselt number increase, especially in the higher nanoparticle
Table 3. Average Nusselt number on the solid viscosity for, , and.
Table 4. Effect of the force terms on the average Nusselt number for, and.
volume concentrations. According to Table 4, by neglecting the drag term at, the average Nusselt number increases 0.03 and it increases 0.24 at. Consequently, the particle-particle interaction and virtual mass forces can be neglected in the mathematical models of nanofluid.
6.3. Enhancement of the Nanofluid Heat Transfer in Comparison to the Pure Water
Figure 4 shows the increase of the average Nusselt number of the nanofluid in comparison to the pure water for different nanoparticle volume concentration. The results show the non-linearly increase in the nanofluid heat transfer enhancement when the nanoparticle volume concentration increases. In fact, the nanofluid heat transfer enhancement is caused by the presence of the copper nanoparticles that increases the thermal conductivity coefficient of the fluid. At it has been observed that the percentage enhancement in the average Nusselt number at Re = 100 is 7.6% more than the corresponding percentage at Re = 50.
The enhancement of the nanofluid heat transfer in comparison to the pure water can be seen from Figure 5 that the coefficient of performance (COP) increases with an increase in the nanoparticle volume concentration for every Reynolds number. Another important observation in Figure 5 is that the lower
Figure 4. Comparison of percentage enhancement in average Nusselt number with respect to pure water for different volume concentrations and ER = 3.
Figure 5. Coefficient of performance of the nanofluid versus Reynolds number for different nanoparticle volume concentration, and.
Reynolds number affects more on increasing the COP for a definite nanoparticle volume concentration. Here, the COP of the nanofluid is obtained as
where subscripts nf and pw show the nanofluid and pure water, respectively.
6.4. Expansion Ratio Effect on the Nanofluid Heat Transfer and Flow
Figure 6. Non-dimensional pressure drop versus Reynolds number for different expansion ratio and.
pressure drop increases with decrease of the expansion ratio of the microchannel but it is small for all Reynolds numbers.
Figure 7 depicts the effect of Reynolds number and the expansion ratio of the microchannel on the average Nusselt number for. It has been observed that the average Nusselt number increases with increase in Reynolds number and decrease in the expansion ratio of the microchannel. Reynolds number increases the convection effects but the increase of the expansion ratio decreases the pressure drop and makes the larger recirculation zone that decreases the convection effects and causes the lower heat transfer.
The increasing of the Reynolds number increases the expansion pressure. The pressure forces dominate the viscous force in the critical Reynolds number where the viscous force cannot hold the flow structure in a symmetric couple of vortexes anymore. In this case, the reattachment length of the upper wall decreases when the reattachment length of the lower one increases and the bifurcation of the reattachment length occurs in the critical Reynolds number. Figure 8 shows that the critical Reynolds number decreases with increase in expansion ratio of the microchannel and the bifurcation point moves to the left side. It can be seen from Figure 8 that the trinity occurs in the second critical Reynolds number for the reattachment length with more increase in the Reynolds number.
The trend of the first critical Reynolds number versus expansion ratio is a straight line. In Figure 9, a correlation for the first critical Reynolds number and expansion ratio is obtained by using the linear regression.
The first and second critical Reynolds numbers for different expansion ratios are shown in Table 5.
Figure 7. Average Nusselt number versus Reynolds number for different expansion ratio and.
Figure 8. Comparison of the bifurcation for different expansion ratios.
Figure 9. Linear regression for the first critical Reynolds number.
Table 5. The first and second critical Reynolds numbers for different expansion ratios.
A two-phase model is considered to simulate the nanofluid flow inside the microchannel and among two-phase models, Eulerian-Eulerian model is very efficient because of considering the relative velocity and temperature of the phases and the nanoparticle concentration distribution. The governing mass, momentum and energy equations for both phases are solved using the finite volume method. In solving the flow equations for both phases, the SIMPLE algorithm is modified for the coupling of the velocity and pressure and the continuity equations for both phases are combined in order to create the pressure correction equations. However, the Eulerian-Eulerian modeling results show higher heat transfer enhancement in comparison to pure water, so that for a 2% copper-water nanofluid, it has been observed a 35% increase of the heat transfer. The average nusselt number increases with increase in Reynolds number and nanoparticle volume concentration and decrease of the expansion ratio of the microchannel. For a constant volume concentration, the lower Reynolds number causes the larger average Nusselt number ratio. The critical Reynolds number for bifurcation decreases with increase in the expansion ratio of the microchannel. The importance and developments of microfluidic devices, like expansion microchannel, has made the investigation of the flow and the heat transfer of nanofluid in sudden expansion microchannel so important and practical. On the other hand, the two-phase models can be used instead of single phase model very well.