Abstract
The effects of wavy geometry on natural convection in an enclosure with two wavy vertical walls and filled with fluid saturated porous media are investigated numerically by using finite element method. The wavy enclosure is transformed to a unit square in the computational domain and the finite element formulations are solved in terms of -coordinate based on iterative method. In order to investigate the effects of interested parameters, the values of wave amplitude ( = 0.05 and 0.1) and number of undulations (n = 1 and 2) are chosen with constants Ra = 105, Da = 10−3, and Pr = 0.71. It is found that the increase in number of undulations has small effect on natural convection inside the enclosure whereas the increase in wave amplitude reduces the strength of convection because higher wave volume plays a barricade role.
1. Introduction
The study of natural convection in an enclosure filled with fluid saturated porous media has received attention because of many applications in science and engineering, for example, geophysics, geothermal reservoirs, insulation of building, heat exchanger design, building structure, and so on. These applications motivate many researchers to perform the numerical simulation for investigating the flow pattern, temperature distribution, and heat flow. Literature reviews on natural convection inside the rectangular and nonrectangular enclosures having different temperature boundary conditions and filled with porous medium are available in [1–4].
The geometrical complexity is increased when the boundary of the enclosure becomes wavy. Wang and Chen [5] analyzed the rates of heat transfer for the flow in a complex wavy channel transformed into a parallel-plate channel by a simple coordinate transformation and solved the obtained equations by the spline alternating direction implicit (SADI) method. Liaw and Char [6] also used a simple coordinate transformation to transform a complex wavy microtube into a smooth circular tube in the computational domain and solved the transformed governing equations by the finite difference method to investigate the electroosmotic flow and heat transfer. The study of natural convection in the enclosure having two wavy walls (top and bottom) and containing internal heat sources at different wave ratios is found in [7]. The enclosure geometry is transformed into computational domain using algebraic coordinate transformation and the transformed governing equations are solved using finite volume method. Dalal and Das [8] were interested in a cavity having three flat walls and wavy right vertical wall consisting of one, two, and three undulations and the heat function equation is utilized in [9] to visualize the heat transfer by fluid flow. The governing equations are solved using finite volume method in the nonorthogonal body-fitted coordinate system. Non-Darcian effects on natural convection in the enclosure with left wavy wall were studied by Khanafer et al. [10]. The governing equations on the physical domain are discretized by using finite element method based on the Galerkin weight residual method with biquadratic interpolation functions and the equations resulting from the discretization are solved using iterative method called segregate solution algorithm. Hasan et al. [11] employed the six-node triangular elements to develop the finite element equations for a heated enclosure with different sinusoidal corrugation frequencies of the top surface and Newton's method is used for solving the nonlinear residual equations. Hasan et al. [12] used the finite element method with the six-node triangular elements for solving unsteady natural convection problem in the heated enclosure with sinusoidally corrugated side walls.
In this study, the effects of wavy geometry on natural convection in the enclosure with wavy side walls are investigated. The numerical solution procedure is based on the finite element method and the finite element formulations are solved in terms of -coordinate in the computational domain. The physical model, boundary conditions of wavy enclosure, mathematical formulation, and solution procedure are described in the next section and the results of a variation in wavy geometry, number of undulations, and wave amplitude are then presented.
2. Mathematical Formulation
2.1. Physical Model and Boundary Conditions
Figure 1 shows the physical model and boundary conditions of the enclosure having two wavy vertical walls. The expressions of the left and the right walls are given by respectively. The fluid saturated porous media considered as incompressible and Newtonian is contained inside the enclosure and its properties are constant except the density variation. The enclosure is nonuniformly heated () on the bottom wall and cooled () on the top wall while the vertical walls are adiabatic. The temperature value assigned to the bottom wall is where is the length of the enclosure, and are distances (m) in horizontal and vertical directions, respectively, is wave amplitude, is number of undulations, and and are temperatures (K) of hot and cold walls, respectively. It is assumed that there is no slip on a boundary. Thus, the velocity components in and directions are and for all boundaries.
(a)
(b)
(c)
2.2. Governing Equations
The governing equations for steady two-dimensional natural convection flow in a porous enclosure can be expressed by the differential equations of the conservation of mass, the conservation of momentum, and the conservation of energy [13]: where is the density (kg ), is the pressure (Pa), is the kinematic viscosity (), is the permeability of the porous medium, is the acceleration due to gravity (m), and are volume expansion coefficient () and thermal diffusivity (), respectively, and is the temperature (K). By using the following change of variables: the governing equations in terms of dimensionless are where and are dimensionless distance in horizontal and vertical directions, respectively, and are dimensionless velocity components in and directions, respectively, is dimensionless pressure, is dimensionless temperature, and , , and are Prandtl, Darcy, and Rayleigh numbers, respectively. The transformed boundary conditions are changed to and for all boundaries and the temperatures are The pressure is eliminated by using the penalty finite element method with penalty parameter such that , where [14]. Substitution of into (7) yields In order to study the flow characteristic, the fluid motion is displayed by using the stream function which is defined as and [15] and yields an equation with the boundary conditions for all boundaries. To visualize the heat transfer by fluid flow, heat function for a two-dimensional convection problem in dimensionless form can be defined as [16] such that Heat transfer rate can be measured in terms of the local Nusselt number () defined by [13] where represents the normal direction on a plane. The local Nusselt number computed along heated bottom wall is defined by and the average Nusselt number at the bottom wall is
2.3. Finite Element Equations
In order to obtain the finite element equations of (8) to (12), the wavy enclosure in -plane called physical domain is transformed to a unit square in -plane called computational domain and the finite element equations are formulated in terms of -coordinate. Let be a coordinate transformation in vector form mapping any point (,) in the physical domain to the points (,) in the computational domain such that and . By using an algebraic method, transfinite interpolation [17], we have the transformation , where and are related to and by the following equations: which transforms the wavy enclosure to a unit square in the computational domain. The computational domain is divided into equal square elements as shown in Figure 1(b) in which nine-node biquadratic elements are used in this study. Using the Galerkin weight residual method with the approximations of , , , , and as where represents column matrix and are shape functions defined for an individual element with node numbering shown in Figure 1(c), yields the element equations in matrix form as follows: energy equation: momentum equations: stream function equation: heat function equation: Due to the fact that the solutions (18) are approximated on a computational domain, the change of variables in double integrals is used and , can be computed by the following relations: where is Jacobian: and (). The shape functions for biquadratic element with node numbering shown in Figure 1(c) are defined as where These finite element equations are assembled to form the global equations and boundary conditions are applied prior to solving for new solutions at nodal point.
2.4. Solution Procedure
The finite element equations formulated in terms of -coordinate in previous section are solved by iterative method. The values of velocity components at the nodal points are first initiated and the new temperature values are obtained by solving the energy equation with the initiated velocity. The new temperature values are substituted into the nonlinear momentum equations and Newton-Raphson method is applied to obtain the new nodal velocity components. The nodal velocity components are then updated and again substituted into the energy equation for new temperature values in the next iteration. The iterative process is continued until the relative error is satisfied: represents the unknown variables , , and at the th iteration. After the process is terminated with the convergence criteria (27), the values of , , and are substituted into the stream function and heat function equations to obtain the values of and , respectively.
3. Results and Discussion
The finite element method and solution procedure proposed in previous section are performed on the computational domain based on biquadratic elements to obtain numerical solutions presented by streamlines, isotherms, and heatlines. Figures 2–4 show the effects of the variations in wave amplitude ( and 0.1) and undulations ( and 2) on the flow field, temperature distribution, and heat transfer within the enclosure having two wavy vertical walls for , , and . Streamlines presented in Figure 2 illustrate the effects on fluid motion when and are changed. It can be seen that two main circulations are formed inside a cavity and symmetric with respect to the vertical center line. As seen from values of streamlines (), the left half is positive and the right half is negative indicating that flow circulations rotate in anticlockwise and clockwise directions, respectively. Because the enclosure is heated from below and the buoyancy force generated is strong enough to initiate fluid convection, fluid moves upwards and impinges to the top of the enclosure and then flows down along the wavy side walls. Due to surface waviness on side walls, the shape of fluid circulation near that region is corrugated. It is observed that the shape of fluid flow for with is squeezed due to higher wave volume in middle portion of the enclosure and flow intensity is found to be lesser compared to that of . The strength of flow circulation for can be decreased by increasing surface waviness to and the magnitude of circulation becomes smaller compared to the other cases.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
The effects of wave amplitude and surface waviness on temperature distribution are presented in Figure 3. As seen from temperature contour, a variation of temperature is high in lower region and less in upper region. For , temperature distribution with shows similar trend with , the temperature contours which disperse from middle portion of bottom wall push the isotherms with towards the upper portion of the enclosure due to the influence of convection. The temperature lines also appear in a small area of bottom corners. As is increased to 0.1, the temperature distribution is smooth curve because flow intensity is obstructed by wave volume and the area of temperature distribution in lower portion of the enclosure is found to be smaller because of the change in size of flow circulation.
Figure 4 shows heatlines pattern for different wave amplitudes and undulations. It can be seen that heat transfer occurs from the bottom towards the side walls and heatlines pattern is symmetric with respect to the middle line resulting from the symmetry of flow pattern and temperature distribution. The results of increasing number of undulations for show similar trend with in which heatlines along wavy side walls are little distorted, whereas the increase of wave amplitude from 0.05 to 0.1 causes smaller region of heat transfer.
The local Nusselt number in Figure 5 is plotted to show the heat transfer rate along bottom wall for different wave amplitude and number of undulations. The distribution of local Nusselt number is symmetric with respect to the line and the value of local Nusselt number is negative near the left and right edges of the bottom wall for all cases. Figure 5 illustrates that heat transfer rate increases from both edges towards the central region. For with and , heat transfer rate is maximum at and and then slightly decreases to the center. It can be seen that the local Nusselt numbers for all cases are almost the same at . Table 1 shows the average Nusselt number () for different wave amplitude and number of undulations calculated for heated bottom wall. The values of indicate that heat transfer is decreased when is increased from to for both and . It is observed that the increase in surface undulations has small effect on .
4. Conclusion
The effects of wavy geometry on natural convection in an enclosure having two wavy vertical walls and filled with fluid saturated porous media have been investigated numerically by using finite element method. The wavy enclosure was transformed to a unit square in computational domain and the finite element formulations in terms of -coordinate were derived and solved by iterative method in which Newton-Raphson method was used for nonlinear equations system. The numerical solutions are presented by streamlines, isotherms, and heatlines. In order to investigate the effects of interested parameters, the values of wave amplitude ( and 0.1) and number of undulations ( and 2) were chosen with constant , , and . From the study results, it is found that(i)as wave amplitude is constant, the flow intensity can be decreased by increasing number of undulations;(ii)the increase in wave amplitude reduces the strength of convection inside the enclosure because higher wave volume plays a barricade role;(iii)temperature distribution in lower portion of the enclosure and region of heat transfer becomes smaller when wave amplitude is increased;(iv)the local Nusselt number is found to be larger for and compared to the other cases due to lower surface waviness; there is a little change in average Nusselt number for increasing number of undulations;(v)the increase in number of undulations has small influence on natural convection compared to the increase in wave amplitude.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgment
This paper is supported by the Centre of Excellence in Mathematics, the Commission on Higher Education, Thailand.