Numerical study on phase change of water flowing across two heated rotating circular cylinders in tandem arrangement

Abstract

The problems of fluid flow and heat transfer phenomena over an array of cylinders are quite prominent in fluid dynamics and industry applications. The current work focuses on fluid flow and heat transfer analysis over two heated rotating cylinders arranged in tandem. The flow of water over heated cylinders faces a phenomenon of phase change from liquid (water) to vapor phase (steam). The mechanism of this phase change is studied through a numerical simulation supplemented with verification of the code and validation. The problem is simulated when flows from two cylinders in a tandem arrangement become interacting and non-interacting. The Eulerian model is used during simulation to comprehend the multiphase phenomena. The volume fractions of both the phases such as water and vapor and heat transfer coefficients of both the cylinders have been computed and presented as findings of the problem. The mass and heat transfer mechanism is unidirectional from one phase to the other phase. The vapor fraction of each phase is to be observed and compared when three different rotations are given to the two cylinders immersed in a turbulent flow of water.

Introduction

The problems of fluid flow and heat transfer phenomena over an array of cylinders are quite prominent in fluid dynamics and industry applications. ¹ These problems give rise to some of the important aspects in fluid dynamics theory such as fluid flow interaction, interferences in flow, vortex dynamics and a variety of engineering applications such as compact heat exchangers, cooling of electronic equipment, nuclear reactor fuel rods, cooling towers, chimney stacks, offshore structures, hot-wire anemometry, and flow control. The mentioned structures are subjected to air or water flows and therefore, experience flow-induced forces which can lead to their failure over a long time.

Some of these situations of natural convection (and hence phase change) exist in the food industry where hot surfaces, such as retorts which may be vertical or horizontal cylinders, are exposed with or without insulation to colder ambient air. It occurs when food is placed inside a chiller or freezer store in which circulation is not assisted by fans.

The essence of the situation is to develop an understanding of convective heat transfer from the hot surface to the flowing fluid and hence its phase change. Basically, with respect to the free stream direction, the configuration of two cylinders can be classified as tandem, side-by-side and staggered arrangements. Quite a few studies on these problems have been carried out analytically, experimentally, and numerically, especially under the configuration of two tandem cylinders for simplicity. Some of the outstanding research activities in the above field have focused on the effect of spacing between the cylinders on the flow characteristics and heat transfer around them. It is observed, therefore, that the nature of the flow relies upon the arrangement of cylinders.^1–3

Numerical simulations of flow over a pair of circular cylinders have been carried out by applying different methods that are mainly based on finite-element formulation. Mittal et al. ² have investigated the problem numerically using a stabilized finite element method and reported their study for the Reynolds numbers of 100 and 1000 in tandem and staggered arrangements for different spacings and concluded that at Re = 1000 and L/D = 2.5, unlike Re = 100, in which the flow converged to initial steady state after some transience, the shear layers caused instability due to the increased velocity of flow. Increasing the gap to L/D = 5.5, the flow at Re = 100 showed unsteady behavior. It was observed that the Strouhal numbers that are associated with the vortex shedding of the twin cylinders could take on the same value.

Numerical investigation of the characteristics of two-dimensional heat transfer in a steady laminar flow around two rotating circular cylinders in a side-by-side arrangement and validation for the large gap-spacing between cylinder surfaces is found in literature. ⁴

In this respect, wake interaction between two circular cylinders in tandem and side-by-side arrangements was studied experimentally by some researchers such as Zhang and Melbourne, ⁵ Liu et al., ⁶ and Ryu et al. ⁷ At low Reynolds numbers, Liu et al. ⁶ employed also the unstructured spectral element method to investigate the flow pattern of two side-by-side cylinders for different spacings.

The flow pattern for tandem arrangement has recently been studied numerically by Mahir and Altac, ⁸ Singha and Sinhamahapatra, ⁹ Ding et al., ¹⁰ and Kitagawa and Ohta ¹¹ for both laminar and turbulent regimes. A numerical study on three-dimensionality effects in the wake of two fixed tandem cylinders at Re = 220 has been performed by Deng et al. ¹² They used the virtual boundary method to apply the no-slip condition. Like 2D case, they found the critical spacing range for which instability occurred at 3.5 ≤ L/D ≤ 4, implying that for L/D ≤ 3.5, the flow wake maintained a 2D state while for L/D ≥ 4, three-dimensionality effects appeared in the wake.

As a matter of fact, a few numerical experiments have been done on flow over a pair of cylinders at high Reynolds numbers. One of the most recent significant studies involves a 3D simulation of flow over two tandem cylinders at the sub-critical Re = 2.2 × 10⁴ by Kitagawa and Ohta. ¹¹ They changed the gap from 2D to 4D and analyzed the interference effect and vortex interaction of two cylinders. Their results showed good agreement with the experimental data at the same Reynolds number. Of noticeable experimental works at subcritical Reynolds numbers, one can mention the studies of Jungkrona et al., ¹³ at Re = 2 × 10⁴ and Moriya et al. ¹⁴ at Re = 6.5 × 10⁴. They thoroughly investigated the flow characteristics of two tandem cylinders.

The characteristics of the flow in both laminar and turbulent regimes ¹⁵ have been numerically simulated for two-dimensional viscous flow around two circular cylinders in a tandem arrangement.

The flow visualization parameters, the Strouhal numbers, and drag and lift coefficients are comprehensively presented and compared for different cases in order to reveal the effect of the Reynolds number and gap spacing on the behavior of the flow. The above literature review has encouraged the authors to contribute toward understanding the phase change phenomenon while water is flowing over two circular cylinders of equal diameters and in tandem arrangement by changing the inter-cylinder spacing between the two heated cylinders and Reynolds number of the flow.

The fluid flow (either laminar or turbulent) over these cylinders with certain heat flux is more relevant as far as the challenges of phase change are concerned. In the present case, it is deemed appropriate to attempt to simulate numerically the phase change phenomenon. This attempt would throw light upon the various challenges and opportunities of phase change from water to steam in a simple case like two-cylinder arrangement in tandem. This is kept as one of the objectives for a faster understanding of the physics of phase change. The second objective of the existing problem is to understand the effect of rotation of one cylinder on the overall phase change process happening due to the flow of water over two heated circular cylinders arranged in tandem. This situation is quite evident in many applications. One such application is the cooling off of nuclear fuel rods by immersion in water and making a provision to circulate water over the cylindrical rods. Another application of rotating cylinders is observed in the landing gears of an airplane and the flow of air over them. Since there is scarcity of earlier studies where rotation of a cylinder is considered in the phase change phenomenon, the attempt made in solving the above-mentioned problem would definitely bring some findings for appreciating the physical scenario in convective heat transfer in the presence of rotation. Due to the lack of sufficient literatures, it is, therefore, considered to analyze the existing literature on fixed cylinders in tandem arrangement as outlined below. After validating the output of a workhorse problem like flow over a circular cylinder when the Reynolds number is 200, the simulation is carried out for the study of phase change of water flowing over two cylinders arranged in tandem and, during the study, rotation is provided to one cylinder(the one that touches water first, called upstream cylinder) and then the overall phase change (in terms of volume fraction) is observed.

The fluid flow over an isolated cylinder has been substantially studied and the results of which have been expressed through relevant literatures in the past. Observing the flow patterns from their origin to growth within a flow over an isolated cylinder, researchers have become quite inquisitive in investigating the situation where there is more than one cylinder. Before embarking on solving a problem where a two-cylinder arrangement is considered, it is felt that the approach and essence of this research would be incomplete without solving a flow over a single cylinder. To validate a few outcomes of the flow over a single cylinder, leading literatures were reviewed and results of laminar flow (Re = 200) over a single cylinder is presented here in this paper.

The above literature review has encouraged the authors to contribute toward understanding the phase change phenomenon while water is flowing over two rotating circular cylinders of equal diameters and in tandem arrangement by changing the inter-cylinder spacing between the two heated cylinders and Reynolds number of the flow.

The fluid flow (either laminar or turbulent) over these cylinders with certain heat flux is more relevant as far as the challenges of phase change are concerned. In the present case, it is deemed appropriate to attempt to simulate numerically the phase change phenomenon. This attempt would throw light upon the various challenges and opportunities of phase change from water to steam in a simple case like two-cylinder arrangement in tandem. This was kept as one of the objectives for a faster understanding of the physics of phase change.

Physical problem and mathematical modeling

Set up

A sketch of the numerical set-up (shown in x–y plane) for two rotating circular cylinders in tandem arrangement with boundary conditions is shown in Figure 1. The necessary dimensions of the fluid domain are expressed in terms of the diameter (D) of the cylinder. The diameters of both the cylinders are the same. The components of the fluid velocity (U) are u, v in the x- and y-directions, respectively. Here T is the temperature which has a free stream value at the inlet of the flow (left to right in the flow domain) and it is equal to the wall temperature specified at the boundary of the solid surfaces of both cylinders. Water enters at the left inlet and comes in contact with rotating cylinder 1 and rotating cylinder 2 with a distance of “L” between them.

Figure 1.

Sketch definition of numerical set-up (shown in x–y plane) for two circular cylinders in tandem arrangement with boundary conditions.

Mathematical modeling

The Eulerian model has been adopted in the current problem. ¹⁶ In this model, a set of n momentum and continuity equations for each phase have been solved. Coupling is achieved through the pressure and interphase exchange coefficients. The manner in which this coupling is handled depends upon the type of phases involved; granular (fluid–solid) flows are handled differently than non-granular (fluid–fluid) flows. For the current problem, coupling has been handled by considering non-granular flows.

The description of multiphase flow as interpenetrating continua incorporates the concept of phasic volume fractions, denoted here by ${\alpha }_q$ . Volume fractions represent the space occupied by each phase, and the laws of conservation of mass and momentum are satisfied by each phase individually.

The volume of phase $q\text{'}$ is defined by $V_q$

${\mathbf{V}}_{\mathbf{q}}=\int {\alpha }_{\mathbf{q}}\text{d}V$

(1)

where

$\sum _{\mathbf{q}=1}^{\mathbf{n}}{\alpha }_{\mathbf{q}}=1$

(2)

The effective density of phase q is ${\rho \wedge }_{\mathbf{q}={\alpha }_{\mathbf{q}}}{\rho }_{\mathbf{q}}$

where ${\rho }_{\mathbf{q}}$ is the physical density of phase $\mathbf{q}$ . The volume fraction equation may be solved either through implicit or explicit time discretization.

The Eulerian multiphase model allows for the modeling of multiple separate, yet interacting phases. A set of conservation equations for momentum, continuity and (optionally) energy is individually solved for each phase.

The basic set of governing equations used to solve the multiphase flow problem is given below.

The general conservation equations for conservation of mass, conservation of momentum and energy are presented below.

Conservation of mass. The continuity equation for phase q is

$\frac{\partial }{\partial \mathbf{t}}({\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}})+\nabla ·({\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}}{v\to }_{\mathbf{q}})=\sum _{p=1}^n({\mathbf{m}·}_{\mathbf{pq}}-{\mathbf{m}·}_{\mathbf{qp}})+{\mathbf{S}}_{\mathbf{q}}$

(3)

where $v\to$ _q is the velocity of phase q and $\mathbf{m}·$ _pq characterizes the mass transfer from the p th to q th phase, and ${\mathbf{m}·}_{\mathbf{qp}}$ characterizes the mass transfer from phase q to phase p . The source term S_q on the right-hand side of the above equation is zero, but it can be specified as a constant or user-defined mass source for each phase.

Conservation of momentum. The momentum equation for phase q is

$\frac{\partial }{\partial \mathbf{t}}({\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}}{v\to }_{\mathbf{q}})+\nabla .({\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}}{v\to }_{\mathbf{q}}{v\to }_{\mathbf{q}})=-{\alpha }_{\mathbf{q}}\nabla \mathbf{p}+\nabla .{¯¯\tau }_{\mathbf{q}}+{\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}}g\to +\sum _{\mathbf{p}=1}^{\mathbf{n}}({R\to }_{\mathbf{pq}}+{\mathbf{m}·}_{\mathbf{pq}}{v\to }_{\mathbf{pq}}-{\mathbf{m}·}_{\mathbf{qp}}{v\to }_{\mathbf{qp}})+({F\to }_{\mathbf{q}}+{F\to }_{\mathbf{lift},\mathbf{q}}+{F\to }_{\mathbf{vm},\mathbf{q}})$

(4)

where q th phase stress–strain tensor is

${¯¯\tau }_{\mathbf{q}}={\alpha }_{\mathbf{q}}{\mu }_{\mathbf{q}}(\nabla {v\to }_{\mathbf{q}}+\nabla {v\to }_{\mathbf{q}}^{\mathbf{T}})+{\alpha }_{\mathbf{q}}({\lambda }_{\mathbf{q}}-\frac{2}{3}{\mu }_{\mathbf{q}})\nabla .{v\to }_{\mathbf{q}}¯¯\mathbf{I}$

(5)

Here ${\mu }_{\mathbf{q}}$ and ${\lambda }_{\mathbf{q}}$ are the shear and bulk viscosity of phase q , ${F\to }_{\mathbf{q}}$ is an external body force, ${F\to }_{\mathbf{lift},\mathbf{q}}$ is a lift force, ${F\to }_{\mathbf{vm},\mathbf{q}}$ is a virtual mass force, ${\mathbf{R}\to }_{\mathbf{pq}}$ is an interaction force between phases, and p is the pressure shared by all phases. Then ${v\to }_{\mathbf{pq}}$ is the interphase velocity, defined as follows; if ${\mathbf{m}·}_{\mathbf{pq}}$ (i.e. phase p mass is being transferred to phase q), ${v\to }_{\mathbf{pq}}$  =  ${v\to }_{\mathbf{p}}$ ; if ${\mathbf{m}·}_{\mathbf{pq}}<0$ (i.e. phase q mass is being transferred to phase p), ${v\to }_{\mathbf{pq}}$  =  ${v\to }_{\mathbf{q}}$ . Likewise, if ${\mathbf{m}·}_{\mathbf{qp}}>0$ then ${v\to }_{\mathbf{qp}}$  =  ${v\to }_{\mathbf{q}}$ , if ${\mathbf{m}·}_{\mathbf{qp}}<0$ then ${v\to }_{\mathbf{qp}}$  =  ${v\to }_{\mathbf{p}}$ .

Equation (4) must be closed with appropriate expressions for the interphase force, ${R\to }_{\mathbf{pq}}$ . This force depends on the friction, pressure, cohesion, and other effects, and is subject to the conditions that ${R\to }_{\mathbf{pq}}=-{R\to }_{\mathbf{qp}}$ and ${R\to }_{\mathbf{qq}}=0$ .

Here a simple interaction term is used and stated below

$\sum _{p=1}^n{R\to }_{\mathbf{pq}}=\sum _{p=1}^n{K}_{\mathbf{pq}}({v\to }_{\mathbf{p}}-{v\to }_{\mathbf{q}})$

(6)

where $K_{\mathbf{pq}}(=K_{\mathbf{qp}})$ is the interphase momentum exchange coefficient.

For fluid–fluid flows, each secondary phase is assumed to form droplets or bubbles. This has an impact on how each of the fluids is assigned to a particular phase. For example, in flows where there are unequal amounts of two fluids, the predominant fluid should be modeled as the primary fluid, since the sparser fluid is more likely to form droplets or bubbles. The exchange coefficient for these types of bubbly, liquid–liquid or gas–liquid mixtures can be written in the following general form:

$K_{\mathbf{pq}}=\frac{{\alpha }_{\mathbf{q}}{\alpha }_{\mathbf{p}}{\rho }_{\mathbf{p}}\mathbf{f}}{{\tau }_{\mathbf{p}}}$

(7)

where f is the drag function, is defined differently for the different exchange-coefficient model and ${\tau }_{\mathbf{p}}$ , the “particulate relaxation” time, is defined as

${\tau }_{\mathbf{p}}=\frac{{\rho }_{\mathbf{p}}{\mathbf{d}}_{\mathbf{p}}^2}{18{\mu }_{\mathbf{q}}}$

(8)

where ${\mathbf{d}}_{\mathbf{p}}$ is the diameter of the bubbles or droplet of phase p .

The drag function f is not computed separately in the problem. However, it is given here as a part of the explanation of the basic equations involved in the entire computing work of the project.

Nearly all definitions of f include a drag coefficient (C_D) that is based on the relative Reynolds number (Re). It is this drag function that differs among the exchange-coefficient models.

For all these situations, $K_{\mathbf{pq}}$ should tend to zero whenever the primary phase is not present within the domain. To enforce this, the drag function f is always multiplied by the volume fraction of the primary phase q , as is reflected in equation (7).

The Schiller and Naumann model ¹⁷ is the default method, and it is acceptable for general use for all fluid–fluid pairs of phases and hence f is expressed as

$\mathbf{f}=\frac{{\text{C}}_{\text{D}}\text{Re}}{24}$

(9)

where

${\text{C}}_{\text{D}}=\{24(1+0.15\mathrm{Re}0.678)/\text{Re}\text{Re}\le 10000.44\text{Re}>1000$

(10)

and Re is the relative Reynolds number. The relative Reynolds number for the primary phase q and secondary phase p is obtained from

$\frac{{\rho }_{\mathbf{q}|{v\to }_{\mathbf{p}}-{v\to }_{\mathbf{q}}|{\mathbf{d}}_{\mathbf{p}}}}{{\mu }_{\mathbf{q}}}$

(11)

Lift forces. For multiphase flows, the effect of lift forces on the secondary phase particles (or droplets or bubbles) is included. These lift forces act on a particle mainly due to velocity gradients in the primary-phase flow field. The lift force will be more significant for larger particles, but for simplified analysis, it can be assumed that the particle diameter is much smaller than the inter particle spacing. Thus, the inclusion of lift forces is not appropriate for closely packed particles or for very small particles.∇

The lift force acting on a secondary phase p in a primary phase q is computed from Drew and Lahey ¹⁸

${\mathbf{F}}_{\mathbf{lift}}=-0.5{\rho }_{\mathbf{q}}{\alpha }_{\mathbf{p}}({v\to }_{\mathbf{p}}-{v\to }_{\mathbf{q}})\times (\nabla \times {v\to }_{\mathbf{q}})$

(12)

The lift force ${F\to }_{\mathbf{lift}}$ will be added to the right-hand side of the momentum equation for both phases ( ${F\to }_{\mathbf{lift},\mathbf{q}}=-{F\to }_{\mathbf{lift},\mathbf{p}}$ ). In most cases, the lift force is insignificant compared to the drag force, so there is no reason to include this extra term. If the lift force is significant (e.g. if the phases separate quickly), it may be appropriate to include this term. ${F\to }_{\mathbf{lift}}$ is not included. The lift force and lift coefficient can be specified for each pair of phases, if desired.

The computation of lift force and drag force is computed for the cylinders by the solver. However, the effect of lift force is assumed insignificant compared to the drag force and hence not included.

Virtual mass force. For multiphase flows, the model includes the “virtual mass effect” that occurs when a secondary phase p accelerates relative to the primary phase q. The inertia of the primary-phase mass encountered by the accelerating particles (or droplets or bubbles) exerts a “virtual mass force” on the particles. ¹⁸

${F\to }_{\mathbf{vm}}=0.5{\alpha }_{\mathbf{p}}{\rho }_{\mathbf{q}}(\frac{{\mathbf{d}}_{\mathbf{q}}{v\to }_{\mathbf{q}}}{\mathbf{dt}}-\frac{{\mathbf{d}}_{\mathbf{q}}{v\to }_{\mathbf{q}}}{\mathbf{dt}})$

(13)

The term $\frac{d_q}{\mathrm{dt}}$ denotes the phase material time derivative of the form

$\frac{d_q(Ø)}{\mathrm{dt}}=\frac{\partial (Ø)}{\partial \mathbf{t}}+({v\to }_{\mathbf{q}}·\nabla )Ø$

(14)

The virtual mass force ${F\to }_{\mathrm{vm}}$ will be added to the right-hand side of the momentum equation for both phases ( ${F\to }_{\mathrm{vm},q}$  =  $-{F\to }_{\mathrm{vm},p}$ ). The virtual mass effect is significant when the secondary phase density is much smaller than the primary phase density. ${F\to }_{\mathrm{vm}}$ is not included.

Virtual mass force is not considered in the current problem because the relative acceleration between the secondary phase (vapor) and the primary phase is neglected.

Conservation of energy in Eulerian multiphase applications. The current problem uses this form of the equation

$\frac{\partial }{\partial \mathbf{t}}({\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}}{\mathbf{h}}_{\mathbf{q}})+\nabla ·({\alpha }_{\mathbf{q}}{\rho }_{\mathbf{q}}{u\to }_{\mathbf{q}}{\mathbf{h}}_{\mathbf{q}})=-{\alpha }_{\mathbf{q}}\frac{\partial {\mathbf{p}}_{\mathbf{q}}}{\partial \mathbf{t}}+{¯¯\tau }_{\mathbf{q}}:{u\to }_{\mathbf{q}}-\nabla ·{q\to }_{\mathbf{q}}\mathrm{\hspace{1em}\hspace{1em}\hspace{1em}}+{\mathbf{S}}_{\mathbf{q}}+\sum _{\mathbf{p}=1}^{\mathbf{n}}({\mathbf{Q}}_{\mathbf{pq}}+{\mathbf{m}·}_{\mathbf{pq}}{\mathbf{h}}_{\mathbf{pq}}-{\mathbf{m}·}_{\mathbf{qp}}{\mathbf{h}}_{\mathbf{qp}})$

(15)

where ${\mathbf{h}}_{\mathbf{q}}$ is the specific enthalpy of the q th phase, ${q\to }_{\mathbf{q}}$ is the heat flux, ${\mathbf{S}}_{\mathbf{q}}$ is a source term that includes sources of enthalpy (e.g. due to chemical reaction or radiation), ${\mathbf{Q}}_{\mathbf{pq}}$ is the intensity of heat exchange between the p th and q th phases, and ${\mathbf{h}}_{\mathbf{pq}}$ is the interphase enthalpy (e.g. the enthalpy of the vapor at the temperature of the droplets, in the case of evaporation). The heat exchange between phases must comply with the local balance conditions ${\mathbf{Q}}_{\mathbf{pq}}=-{\mathbf{Q}}_{\mathbf{qp}}$ and ${\mathbf{Q}}_{\mathbf{pq}}=0$ .

Turbulence models. In comparison to single-phase flows, the number of terms to be modeled in the momentum equations in multiphase flows is large, and this makes the modeling of turbulence in multiphase simulations extremely complex. In the present problem, $\mathbf{k}-\epsilon$ models have been employed which are a mixture turbulence model. The description of the mixture turbulence model is presented below.

The $\mathbf{k}-\epsilon$ mixture turbulence model. The mixture turbulence model is the multiphase turbulence model which has been used in the present computation. It represents the first extension of the single-phase $\mathbf{k}-\epsilon$ model, and it is applicable when phases separate, for stratified (or nearly stratified) multiphase flows, and when the density ratio between phases is close to 1. In these cases, using mixture properties and mixture velocities is sufficient to capture important features of the turbulent flow.

The $\mathbf{k}$ and ε equations describing this model are as follows:

$\frac{\partial }{\partial \mathbf{t}}({\rho }_{\mathbf{m}}\mathbf{k})+\nabla ·({\rho }_{\mathbf{m}}{v\to }_{\mathbf{m}}k)=\nabla ·(\frac{{\mathbf{\mu }}_{\mathbf{t},\mathbf{m}}}{{\mathbf{\sigma }}_{\mathbf{k}}}\nabla \mathbf{k})+{\mathbf{G}}_{\mathbf{k},\mathbf{m}}-{\rho }_{\mathbf{m}}{\epsilon }_{\mathbf{m}}$

(16)

and

$\frac{\partial }{\partial \mathbf{t}}({\rho }_{\mathbf{m}}\epsilon )+\nabla ·({\rho }_{\mathbf{m}}{v\to }_{\mathbf{m}}\epsilon )=\nabla ·(\frac{{\mu }_{\mathbf{t},\mathbf{m}}}{{\mathbf{\sigma }}_{\varepsilon }}\nabla \epsilon )\mathrm{\hspace{1em}\hspace{1em}\hspace{1em}}+\frac{\epsilon }{\mathbf{k}}({\mathbf{C}}_{1\epsilon }{\mathbf{G}}_{\mathbf{k},\mathbf{m}}-{\mathbf{C}}_{2\epsilon }{\rho }_{\mathbf{m}}\epsilon )$

(17)

where the mixture density ${\rho }_{\mathbf{m}}$ and velocity ${v\to }_{\mathbf{m}}$ are computed from

${\rho }_{\mathbf{m}}=\sum _{\mathbf{i}=1}^{\mathbf{N}}{\alpha }_{\mathbf{i}}{\rho }_{\mathbf{i}}$

(18)

and

${v\to }_{\mathbf{m}}=\frac{\sum _{\mathbf{i}=1}^{\mathbf{N}}{\alpha }_{\mathbf{i}}{\rho }_{\mathbf{i}}v\to i}{\sum _{\mathbf{i}=1}^{\mathbf{N}}{\alpha }_{\mathbf{i}}{\rho }_{\mathbf{i}}}$

(19)

The turbulent viscosity, ${\mu }_{\mathbf{t},m}$ is computed from

${\mu }_{\mathbf{t},\mathbf{m}}={\rho }_{\mathbf{m}}{\mathbf{C}}_{\mu }\frac{{\mathbf{k}}^2}{\epsilon }$

(20)

and the production of turbulence kinetic energy, G_k,m is computed from

${\mathbf{G}}_{\mathbf{k},m}{\mathbf{=}\mathbf{\mu }}_{\mathbf{t}\mathbf{,}\mathbf{m}}({\mathbf{\nabla }\mathbf{v}\to }_{\mathbf{m}}\mathbf{+}({\mathbf{\nabla }\mathbf{v}\to }_{\mathbf{m}})\mathbf{T}){\mathbf{:}\mathbf{\nabla }\mathbf{v}\to }_{\mathbf{m}}$

(21)

The constants in these equations are given below for the single-phase k–ɛ model. The same constants are used while soling the equations for each phase.

${\mathbf{C}}_{\mathbf{\mu }}=0.09,{\mathbf{C}}_{1\mathbf{\varepsilon }}=1.44,{\mathbf{C}}_{2\mathbf{\varepsilon }}\mathbf{=}1\mathbf{.}92\mathbf{,}{\sigma }_{\mathbf{k}}\mathbf{=}1\mathbf{.}0\mathbf{,}{\sigma }_{\varepsilon }\mathbf{=}1\mathbf{.}3$

Boundary conditions

Boundary conditions for the above set up are as follows

Inlet to the domain: Velocity inlet, U_∞ = 1 m/s;

Outlet from the domain: Gauge pressure outlet, p = 0 Pa;

Wall of the domain: No slip wall boundary (top and bottom) s;

Cylinder wall surface: Heat flux, q″ = 10,000 W/m² (for both the cylinders);

Cylinder 1: rotating wall;

Cylinder 2: rotating wall.

Solution procedure and mesh independence study

The continuity, momentum and energy equations are solved as per the Eulerian model for multiphase using Ansys^(R) Fluent 12.1 solver. ¹⁶ The phase coupled SIMPLE method has been chosen in the solver to compute the flow variables. The turbulent quantities k and $\epsilon$ are solved as per the mixture turbulence model. The descriptions of both are given in the above section. The iteration of all the steps ends when the full convergence is achieved. Residual values include the momentum equations for each phase, k and $\epsilon$ equations for each phase and pressure correction residual for continuity equation.

During the heat transfer from the walls of the cylinders to water, which is liquid, it starts to go to a new phase, steam. The mechanism of heat and mass transfer from water to steam is set during the multiphase flow settings in the Ansys^(R) Fluent 12.1 solver. The mass transfer from Phase 1 (water) to Phase 2 (steam) is computed using the unidirectional mass transfer mechanism option in Fluent. Momentum, energy, and turbulence are also transported with the mass that is transferred.

The convergence of the numerical solutions is obtained from the above-mentioned problem using the residuals of the values of variables such as continuity of the flow, velocities of two phases (primary and secondary), and the energy of each phase, turbulent kinetic energy and its rate of dissipation. In this work, convergence occurs when the values of total residual in all the above-mentioned equations become smaller than 10⁻⁵. All these values have reached their acceptable steady solutions during the simulation. The solutions are also independent of the mesh resolution. For the present simulation, initial mesh elements are 306,550 and the convergence of residual error are below the above-mentioned value. The mesh elements are increased to 1.5 times due to finer meshing. The simulation is again carried out and convergence criteria satisfied. The values of all field variables obtained from the two simulations are compared and found to be the same.

Grid quality has been satisfactory after checking the skewness as 0.2 and over all grid quality 0.92 (the highest value is 1.0). The automatic gridding method of the Ansys workbench has optimized the grid density appropriately to fit into flow domain and overall grid quality.

Results and discussion

Validation of the code

For validation purpose, experimental setup at NASA laboratory ¹⁹ has been used. In the past, numerous studies have been conducted on the tandem configuration due to its practical application in systems used for cooling, venting, and structural support. In most cases, the studies were performed at subcritical Reynolds numbers (less than 1.5 × 10⁵) and focused on fundamental issues such as boundary layer development, numerical prediction of flow interference, flow-induced vibrations, and characterization of the wake structure for various cylinder spacings. The above-mentioned research carried out at NASS’s GRC has inspired the authors to go ahead with a very specific application that is the cooling of heated cylinders in a tandem arrangement where there is a change of phase of working fluid due to the heat released from the heated surfaces.

Therefore, to validate the approach of the multiphase flow study over a two-cylinder arrangement in tandem, numerical study was completed for the spacing-to-diameter (L/D) ratios such as 1.435 and 3.7 as same as that by NASA laboratory. ¹⁹ Figures 2 and 3 show the comparison of stream lines in the gap. It can be seen from the figures that the numerical results obtained in the present simulation is qualitatively same as reported in the literature. ¹⁹ Qualitative comparison means that the contour plots obtained from numerical simulation are kept side by side with the contours obtained from the experiments to understand the flow features around the two cylinders. Though the L/D and Reynolds number used for the comparison are the same, the working fluid and the phase change mechanisms of numerical study and experimental work are different. It was, hence, deemed appropriate to capture the flow features in terms of stream traces to ensure that the physical phenomena between the cylinders and in the wake zone of downstream cylinder match well. Quantitative comparison is not attempted because the values produced during the experiment and numerical simulations are ought to be different for the reason mentioned above.

Figure 2.

Mean streamline in the gap region (L/D = 1.435): (a) present model and (b) NASA experimentation. ¹⁹

Figure 3.

Mean streamline in the gap region (L/D = 3.7): (a) present model and (b) NASA experimentation. ¹⁹

Verification of interacting and non-interacting flow

Mean drag coefficient ${\text{C}}_{\text{D}}^{\text{M}}$ and mean lift coefficient ${\text{C}}_{\text{L}}^{\text{M}}$ are compared with some of the leading references. From Table 1, it is observed that the computed values of mean drag coefficient, mean lift coefficient and the Strouhal number agree well with the values published by leading researchers in the field.

Table 1.

Comparison of $C_{\text{D}}^{\text{M}}$ , $C_{\text{L}}^{\text{M}}$ , and $\text{St}$ .

	CD	CL	St
Present model	1.41	0.692	0.1902
Liu et al. 20	1.337	0.685	0.1955
Rajani et al. 21	1.3380	0.4276	0.1936
Wang et al. 22	–	0.71	0.1950
Zhang et al. 23	1.34	0.66	0.1970
Linnick and Fasel 24	1.34–1.37	0.71	NA
Farrant et al. 25	1.36–1.39	0.71	NA
He et al. 26	1.36	NA	0.1978
Henderson 27	1.34–1.37	NA	0.1971

To understand interacting and non-interacting flows when two cylinders are kept in a tandem arrangement; a parameter u⁺ (ratio of mean stream wise velocity and free stream velocity) is used. This is a non-dimensional parameter and it is computed along a line passing straight from the bottom to the top of the flow domain through a point lying on the mid-line (passing through a point which is lying at an equal distance from each of the cylinders) across the cylinders. Figure 4 depicts the situation when the flows from each cylinder are interacting and non-interacting. This is plotted against the y⁺ (a non-dimensional quantity and ratio of y and maximum value of the ordinate of the flow domain). The objective of this plot is to find out the influence of spacing on the mean stream wise velocity and recommend the spacing (in terms of D) where the flows are interacting and non-interacting.

Figure 4.

Comparison of the u⁺ values for three different spacings with respect to y⁺ to find interacting and non-interacting flows.

It is observed that u+ values for the spacing 1.5D are influenced by the closeness of the two cylinders (and hence flows are interacting) whereas the values of u+ remain steady (particularly linear inside the domain covering spacing and hence flows are mostly non-interacting) for the spacing of 6D. After confirming that the flows from the two cylinders, while fixed, are not interacting, the objective is extended to investigate the effect of rotation of one of the heated cylinders on phase change process. The cylinders having the same diameters are kept in a tandem arrangement with a gap of 6D, where D is the diameter of the cylinder.

Effect of rotation of cylinders on phase change behavior

The focus of the study is to observe the effect of rotation on the volume fraction of water and vapor. So, the arrangement is to give rotation to both cylinders. Both the cylinders are heated by certain heat source and the heat flux is maintained at 10,000 W/m². The numerical experiment was conducted by giving three different rotations to cylinders. Each rotation is expressed in terms of rate of rotation (indicated by α) which is Dω/2U_∞, where U_∞ is free-stream velocity, D the cylinder diameter, and ω the angular velocity of the cylinder.

Before beginning the numerical simulation, the volume fractions of water and vapor are noted down when both the cylinders are non-rotating. For the case, when both the heated cylinders are fixed, the contours of V_f for the vapor phase is depicted in Figure 5. Figure 6 depicts volume fraction (V_f) of vapor phase at different locations on the surfaces of cylinder 1and cylinder 2 when both the cylinders are fixed at a gap of 6D. V_f of vapor for cylinder 2 at location 0° begins with a value of 0.55, increases to almost 1.0 at location of 180° and then falls below 0.55, fluctuates in the lower half of the cylinder (180°–360°) before ending with a modest value of 0.5 while the V_f of water phase does exactly the reverse on the upper and the lower half of cylinder 2. It is observed that the major phase change happens in the upper half (0°–180°) of cylinder 2. Over the cylinder 1, between 90° and 240°, V_f of the vapor phase has been constantly increasing and then is fluctuating before staying at a value of 0.5. So, cylinder 1 releases more heat to water during a shearing phenomenon which takes place just after its tip and this continues toward cylinder 2 where upper half of it faces a major phase change.

Figure 5.

Contours of volume fraction (V_f) of water vapor over for two fixed cylinders.

Figure 6.

Volume fraction (V_f) of water vapor at different location (in angle measured clockwise for fixed cylinder).

Then, both the cylinders were maintained at three different revolutions per minute (r/min) such as 1000, 3000, and 6000. The corresponding rates of rotations are 104.71, 314.15, and 628.31, respectively. The readings are taken from 0° to 360° at an interval of 30° around each cylinder in the clockwise direction. The contour plots of volume fraction (V_f) of vapor phase for both the rotating cylinders are given in Figure 7. Figures 8 and 9 show the vapor phase at different angles for cylinder 1 and cylinder 2, respectively.

Figure 7.

Contours of V_f for water and vapor over two heated rotating cylinders (a) α = 104.71, (b) α = 314.15, and (c) α = 628.31 of equal diameters and in tandem arrangement.

Figure 8.

Volume fraction (V_f) of water vapor at different location for cylinder 1 (in angle measured clockwise as shown above) when cylinders are (a) α = 104.71, (b) α = 314.15, and (c) α = 628.31.

Figure 9.

Volume fraction (V_f) of water vapor at different location for cylinder 2 (in angle measured clockwise as shown above) when cylinders are (a) α = 104.71, (b) α = 314.15, and (c) α = 628.31.

From the contours plots shown in Figure 7(a) and the graphs shown in Figures 8 and 9, it is seen that for α = 104.71 for cylinder 1, V_f of vapor decreases from a value of 1.0 at 0° to a value close to 0.0 at 90° and then after increases to almost 1.0 at 220°. After that it fluctuates till 360°. As the rotation increases, a fewer fluctuations are visible and at 6000 r/min, the vapor fraction remains almost at 1.0 from 200° to 360° (Figure 7(b) and (c)).

The effect of rotation is also felt in phase change over cylinder 2 but at a uniform way which is displayed in Figure 9. The major change in phase change in terms of V_f for both the water phase and the vapor phase is attributed to an unidirectional mass and heat transfer mechanism between two discrete phases during rotation.

Due to rotation, boundary layer thickness changes on the solid body and due to which shearing action leads to more heat addition to the phase change processes.

Conclusion

In this paper, the effect of rotation on the phase change of water flowing across two heated circular cylinders in a tandem arrangement was studied. The Eulerian model was used during simulation to capture the data. First of all, the findings of the current numerical simulation were compared and validated with the experimental data published in the literature ¹⁹ by considering L/D = 1.435 and 3.7. The difference does exist between the two separate findings because the properties of water used in the current problem differ from that of the air taken as the working fluid in the experiment. The physics of flow in the gap region has demonstrated that the findings of numerical simulation are quite satisfactory because it compares well with the experiment.

The novelty in the problem is the study of phase change of water when the cylinders rotate. At a specific location on the surface of cylinder 2 it is observed that the V_f for vapor phase goes beyond 50% higher at α = 104.71 than that at α = 628.31. This finding concludes that rotation of cylinder has significant effect on the phase change of water flowing across cylinders which are heated and rotating.