LOCAL HYDRODYNAMICS OF FLOW IN A PULSATING HEAT PIPE: A REVIEW

Global thermal performance modeling of Pulsating Heat Pipes (PHPs) requires local, spatio-temporally coupled, flow and heat transfer information during the characteristic, self-sustained thermally driven oscillating Taylor bubble flow, under different operating conditions. Local hydrodynamic characteristics such as velocities, lengths, shapes and profiles of bubbles and slugs, their dynamic contact angles, thickness of the liquid film that surrounds the bubbles, enhanced mixing/ flow circulation within the liquid slugs and net pressure drop along the flow, etc., are needed to predict local heat transfer and thus, the global thermal performance. In this paper, we systematically review the experimental, theoretical/analytical, and modeling methodologies to predict these hydrodynamic properties in unidirectional two-phase Taylor bubble flows, in the context of Pulsating Heat Pipes. Indeed, there is little literature available for oscillating Taylor bubbles flows. In view of the state-of-the-art, we therefore recommend some directions and perspectives for furthering research on understanding and modeling PHPs.


INTRODUCTION
With the advances in mini/micro fabrication techniques, many applications have emerged which utilize single-phase or two-phase flows in such mini/micro geometries to achieve specific goals. Finding and delineating transport mechanisms in such systems is a major challenge and research impetus has greatly increased in this direction so as to achieve systems with better control and improved efficiency. One such device which utilizes self-excited thermally driven twophase flow oscillations for enhanced passive heat transfer is a Pulsating Heat Pipe (PHP), also referred to as an Oscillating Heat Pipe (OHP). This device has attracted considerable research in recent past due to its simplicity, ease of manufacture, favorable operating conditions under specified boundary conditions, vis-à-vis other conventional heat pipes. The transport mechanism of this 'simple' device also poses an excellent opportunity to understand its complex internal two-phase thermohydrodynamics [Akachi, 1990;Akachi et al., 1996;Khandekar, 2004;Vasiliev, 2005;Zhang and Faghri, 2008].
The basic structure of a typical pulsating heat pipe consists of a meandering capillary tube having no internal wick structure, as shown in Figure 1-a. The inside diameter of the tube is sufficiently small (Bo crit ~ 1.84; see Section 2) such that surface tension dominates and capillary slug flow is maintained Zhang and Faghri, 2008]. It can be designed in at least three ways, namely: The entire essence of a PHPs thermo-mechanical physics lies in the closed (constant volume), two-phase, bubble-liquid slug system, spontaneously formed inside the tube bundle, at the time of filling the device, due to the dominance of surface tension. This tube bundle receives heat at one end and is cooled at the other end. Temperature gradients give rise to temporal and spatial pressure disturbances due to resulting phase-change phenomena, i.e. generation and growth of bubbles in the evaporator and simultaneous collapse of bubbles in the condenser. The bubbles act as pumping elements, transporting the entrapped liquid slugs in a complex oscillating-translating-vibratory fashion, resulting in self-sustained thermally driven flow oscillations and ensuing highly efficient passive heat transfer thereof. In addition to the latent heat, considerable amount of sensible heat transfer also occurs in a PHP. While sweeping the evaporator section, a liquid slug accumulates heat, which is eventually transferred to the condenser.
The fundamental transport processes that occur inside the PHP can be visualized, based on Figure 1-c, which suggests the various forces, and heat and mass transfer processes acting on a unit cell of liquid-slug vapor-bubble, as formed inside the PHP. These primary processes are: • The flow pattern in a PHP may be broadly categorized as capillary slug flow, which is characterized by: (a) the flow pattern is 'generally' axisymmetric, at least in vertical flows (in horizontal flows, there will be some asymmetry depending on the Bond number, as explained later in Section 4), (b) the velocity of vapor bubbles relative to liquid slug is somewhat faster [Fairbrother and Stubbs, 1935;Harmathy, 1960;Suo and Griffith, 1964;Wallis, 1969;Fukano and Kariyasaki, 1993;Thulasidas et al., 1995;Angeli and Gavriilidis, 2008].

Frontiers in Heat Pipes
Available at www.ThermalFluidsCentral.org variations of a PHP (c) Taylor bubble train in a vertical 2 mm tube and the forces and transfer processes acting on a unit-cell [Khandekar, 2010].
• Due to the capillary dimensions of the PHP tube, a train of liquid slugs and vapor bubbles having menisci on their edges is formed due to surface tension. Usually a liquid thin film exists around the bubbles (see Section 4.1). The contact angle of menisci, liquid thin film stability and its thickness depend on the fluid-solid combination and the operating parameters (see Section 3.2).
• Liquid slugs and vapor bubbles move against the gravity vector, in its direction or at an angle to it, depending on the global PHP orientation and location of slugs/ bubbles in the up-header or down-header tubes.
• The liquid slugs and vapor bubbles are subjected to pressure force 1 F  and 2 F  (see Figure 1), from the adjoining slugs/bubbles. These are not only caused due to phase-change mass transfer but also due to capillary forces. This is further explained in Section 4.1.1.
• The liquid slugs and vapor bubbles experience internal viscous dissipation as well as wall shear stress. Their relative magnitude decides the predominant force to be considered (see Section 4.1).
• The liquid slugs and vapor bubbles may receive heat, reject heat, or move without any external heat transfer, depending on their location in the evaporator, condenser or the adiabatic section, respectively. Most thermal transport occurs through the thin film and its dynamics plays a crucial role in the overall thermal transport (see Section 5.2).
• In the evaporator, liquid slugs receive heat, which is simultaneously followed by evaporation mass transfer to the adjoining vapor bubbles or breaking up of the liquid slug itself, with creation of new bubbles in between, as a result of nucleate boiling in the slug flow regime; P sat and T sat thus increase locally. Vapor bubbles may frequently come in direct contact with the internal tube surface of the evaporator. In this case, saturated vapor bubbles receive heat via the liquid thin film surrounding them, which is simultaneously followed up by evaporation mass transfer from the film as well as the adjoining liquid slugs. Heat transfer under such conditions is strongly dependent on local film geometry as highlighted later in Section 5.3. The processes in the evaporator are repeated in a reverse direction in the condenser.
• In the adiabatic section, while passing from the evaporator to the condenser, the train of vapor bubbles and liquid slugs is subjected to a series of complex heat and mass transfer processes. Essentially nonequilibrium thermodynamic conditions exist whereby the high pressure, high temperature saturated liquid slugs/vapor bubbles are brought down to low pressure, low temperature saturated conditions existing in the condenser. If ideal adiabatic conditions are maintained, with no axial conduction of heat through the tube wall/ fluid itself, then an inherently irreversible isenthalpic process can accomplish this task. Internal enthalpy balancing in the form of latent heat takes place by evaporation mass transfer from the liquid slugs to the vapor bubbles whereby near saturation conditions tend to be imposed on the system during the bulk transit in the adiabatic section. In real systems, this transit is certainly much more complex with non-equilibrium metastable conditions existing throughout. It must be noted that there occurs no 'classical steady state' in PHP operation as far as the internal hydrodynamics is concerned. Instead, pressure waves and pulsations are generated in each of the individual tubes, which interact with each other, generating secondary and ternary reflections with perturbations. Many operating variables control these self-exited thermally driven oscillations [Khandekar, 2004].
There are several important and unresolved issues with the present understanding of PHPs [Zhang and Faghri, 2008;Das et al., 2010-a]: (i) There is no comprehensive mathematical model to predict the PHP thermal performance for a given boundary condition.
(ii) The understanding of heat transfer and pressure drop under selfexcited thermally driven oscillating two-phase flow inside capillary tubes is not convincing. The complete transport phenomena in the unitcell need to be resolved to predict global heat transfer parameters.
(iii) Multiple unit-cells also interact with each other mutually; merger and coalescence of liquid slugs, breakage of Taylor bubbles under the impact of inertia and surface tension, nucleation inside liquid slugs, confined bubble formation, condensation on liquid films, instabilities, surface capillary waves, etc. are additional complexities for formulation of predictive tools for PHPs.
In this background, we review the present understanding and highlight the important parameters, which need consideration in our global understanding of Pulsating Heat Pipes. Primarily, the flow pattern in these systems is characterized by thermally driven pulsating/ oscillating Taylor bubbles. Therefore, we focus our attention on local hydrodynamic transport behavior of Taylor bubbles. A detailed literature review is presented to highlight the present state of understanding on uni-directional as well as oscillating Taylor bubble flows, in the context of further understanding of PHPs.

CLASSIFICATION OF CHANNELS
It is to be appreciated that PHPs operate in the regime when surface tension forces are dominant in the tubes. For normal heat transfer fluids, this suggests the use of small, mini or micro ducts or channels for the construction of PHP. Due to significant differences of transport phenomena in mini-micro systems as compared to normal sized channels, it is important to define what is meant by the terms mini, meso and micro sizes of the channels. Merely categorizing on the basis of tube geometry is not enough, especially in the context of two-phase flows. The criterion ought to be based on the combination of size and thermo-physical properties of the fluid under consideration.
Channels are generally regarded small when surface tension forces dominate the gravitational forces. A measure of the relative value of the two forces is given by the Eötvös or Bond number: Researchers have suggested different values of the Bond number, ranging from 0.88 to 2π 2 , below which channels can be considered small [e.g. see White and Beardmore, 1962;Brauner, 1992;Kandlikar, 2002]. Bretherton [1961] stated that, for Bo < 1.84, the channel can be considered as small.
According to  and later by Serizawa et al. [2002], the flow channels having D hyd less than or equal to Laplace constant, λ, may be considered as microchannels, where: Kew and Cornwell [1997] give another interpretation of the same quantity in terms of the confinement of bubbles in a channel. They reported that channels having two-phase flow can be considered as mini systems when the Confinement number, Co, as defined below, is greater than 0.5.
Under such conditions the two-phase flows typically have a confined bubble flow regime, which is unique to mini systems, as also is the case with PHPs. As the hydraulic diameter goes on decreasing, the effect of Bond number on various hydrodynamic parameters goes on increasing. In fact, all classification numbers, Bo, Eö, Co, or a constant λ, consider the effect of fluid density, surface tension and channel size on two-phase flow.
If we estimate the relevant threshold critical diameters i.e. transition diameters for considering a channel as 'micro', for common working fluids respectively, for the entire range of operating pressures, we notice that there are considerable differences in the threshold diameters, as defined above by various criteria; so far, no unified methodology is available. In this paper we will confine our discussion to mini or smaller systems where the Bond number is definitely less than 2. Such tubes/channels will also be called as 'capillary' systems.
The foregoing discussion has important implications in defining a PHP. The definition of a threshold diameter needs to be decided so that distinct liquid plugs and vapor slugs can be formed without separation, stratification or agglomeration, at least under adiabatic conditions. For diabatic cases (such as bubble pumps), it is well known that slugbubble systems can form for higher diameters than the critical diameters defined herein [Delano, 1998]. Thus, rather than a certain fixed diameter which classifies the boundary between classical macroscale systems and PHPs, there is a finite transition zone wherein PHPs can effectively function as truly pulsating/oscillating heat transfer devices ].

CONTACT LINES: STATICS AND DYNAMICS
To predict the shape, length and pressure drop under Taylor bubble flows, it is important to understand the basic physics of three-phase contact lines. No sooner a second phase (gas or vapor) is introduced in a parent liquid phase, there is a possibility of formation of three-phase contact lines between the parent fluid, the dispersed phase and the solid container/duct in which the system is stored or is flowing. These threephase contact lines are in addition to the gas-liquid two-phase interfaces. Bubbles are formed when a continuous gas-liquid interface forms a closed volume entrapping the gas or vapor inside it. Gas or vapor can also be trapped between two (or more) distinct three-phase contact lines. Possibility of formation of three-phase contact lines in a PHP is high, especially in the evaporator section where bubble formation and growth of confined bubbles takes place. Elsewhere in the PHP, the Taylor bubbles usually travel with a liquid film around them.
It is therefore prudent to review the available literature on the statics and dynamics of contact lines although exhaustiveness in this section cannot be claimed due to the sheer vastness of the available material. We only highlight the important issues, which ought to be considered while developing models of PHPs.

Three-phase contact lines: Statics
Three-phase triple contact lines are said to be formed when materials in different phases, e.g. solid, liquid and gas (or vapor) intersect. Common examples are a liquid droplet spreading on a solid surface or a liquid meniscus in a capillary tube. The angle between the tangents at the liquid / gas interface and the liquid / solid interface is known as the contact angle θ of the liquid on the solid surface. Typically, when a liquid spreads on a solid surface in the presence of the third phase (gas or vapor), two distinct equilibrium states can be achieved: a partial wetting state or a complete wetting state, as shown in Figure 2-a. This wetting phenomenon of a solid, by a given liquid, and the ultimate equilibrium state is intrinsically connected to complex interplay of physical chemistry, statistical physics, short and long-range forces and to fluid dynamics. The condition for static equilibrium of a triple contact line involving an ideal solid (perfectly smooth and chemically homogeneous), liquid and a gas/vapor is stated as the Young's equation [Carey, 2007]: Consequently, the contact angle, which the liquid mass makes with a solid surface, depends not only on the thermophysical properties of the liquid, but also on the physico-chemical structure of the solid substrate. If all the respective interfacial tensions are known, the equilibrium static contact angle can thus be estimated. The SV and SL tensions are not only difficult to estimate but predictive theories are also quite weak. Moreover, validity of Eq. (4) is restricted to ideal surfaces and planer interfaces, which are rarely found in engineering practice. This makes the Young's equation rather impractical and frequently, careful experiments are needed to estimate the static 'equilibrium' contact angle. It must be noted here that true 'thermodynamic' equilibrium (mechanical, chemical and thermal) may also be extremely difficult to achieve, especially for non-volatile fluids; at least, mechanical equilibrium (force balance) is mandatory in the definition of θ eq . It is in this light that Young's equation is best derived by considering a reversible change in contact position, using global energetic arguments, i.e. the nature of the contact line region, over which intermolecular forces are acting, is not considered. Accordingly, θ eq is understood to be measured macroscopically, on a scale above that of long range intermolecular forces [Bonn et al., 2009]. Depending on this equilibrium contact angle, solid substrates, for a given fluid, may behave as lyophilic (0° ≤ θ eq ≤90°) lyophobic (90° ≤ θ eq ≤150°) or super-lyophobic (150° ≤ θ eq ≤ 180°).
The equilibrium wetting conditions of a liquid on a solid surface can also be stated in terms of the equilibrium Spreading Coefficient, S eq , which represents the surface free energy per unit area σ SV relative to its value for complete wetting, i.e.: For complete wetting, the equilibrium Spreading coefficient S eq is zero while for partial wetting conditions it is a negative quantity. As noted above, systems may not achieve complete equilibrium, at least in the time frame of the experiments, and to deal with such conditions, the Initial Spreading Coefficient is defined as follows, The initial spreading coefficient is a more relevant quantity of interest, especially under dynamical conditions of the contact line. Note that σ S0 is the surface tension of the 'dry' solid surface, contrary to σ SV which is the surface tension of the solid surface in equilibrium with the gas or vapor. It is known that this latter quantity is always lower than the former because of adsorption of gas/vapor molecules which lowers the surface tension of the solid. Thus, S i is always larger than S eq . For situations where equilibrium is not reached, and S i is negative (which means that S eq is also more negative), the contact line will display a finite static contact angle, θ i , which will not be the same as the equilibrium contact angle, θ eq . When S i is positive, although the droplet will spread and tend to flatten, there is no surety whether S eq will necessarily be zero or negative. Thus, the drop may initially spread and if enough time is given for achieving true equilibrium, it may either become a case of complete wetting (S eq = 0) or it may retract back to a case of partial wetting (S eq < 0). In general, though, if S i is large and positive, the equilibrium state is usually characterized by S eq = 0, however there can be exceptions when S eq < 0 [de Gennes, 1985;Bonn et al., 2009]. It is difficult to get values of S eq and even more difficult to estimate S i ; careful experiments are the only meaningful alternatives, at least for engineering practice.
Application of the above theory to real surfaces which are neither smooth nor chemically homogeneous requires some corrections. Under such situations, the static contact angle turns out to be non-unique and it is observed that it depends on the way the surface was "prepared" [Huh and Mason, 1977;Léger and Joanny, 1992]. Experiments involving the estimation of wettability are extremely sensitive to physical and chemical heterogeneities of the substrate. One of the first attempts at understanding the explicit influence of physical roughness (but chemically homogeneous surface) on wetting is due to Wenzel [1936], who proposed: * cos cos eq r θ θ = (7) Here θ * is the apparent contact angle, r is the roughness (r = 1 is a smooth surface and r > 1 is a rough surface). It embodies two types of behaviors. If 90 eq θ <° (lyophilic substrate), we will have * eq θ θ < since r > 1. Likewise, if 90 eq θ >°(lyophobic substrate), we will have * eq θ θ > . It is known that these angles can be sufficiently tuned by varying solid roughness [Shibuichi et al., 1996]. The same line of reasoning can be applied to a surface that is physically homogeneous but having n chemically heterogeneous species on its surface, which gives the Cassie-Baxter's relation: Here, θ * is the apparent contact angle on the chemically inhomogeneous surface, θ eq-i are the respective contact angles on the n chemically homogeneous surfaces respectively, i f  is the fractional area occupied by the i th specie. Therefore, the apparent angle θ * (which is restricted to the interval (θ eq-1, θ eq-n ) is given by an average involving the angles characteristic of each constituent, but the average is applied to the cosines of these angles. Thus, we see that the contact angle of a liquid on a substrate depends on the roughness and chemical homogeneity of the surface; the three-phase contact line of the drop is deformed due to the physical and/or chemical topographic heterogeneities (this is true whenever the length scale of the threephase contact line is larger than the typical length scale of the physical/ chemical in-homogeneity).
Real time engineering surfaces are both physically and chemically heterogeneous and tubes used of making PHPs will be no exceptions. Therefore, we need to restore to systematic experimental work to get the real picture of actual contact angles. In doing such experiments and subsequent post-processing of the data, there are two more related issues which need careful consideration (i) 'Apparent' nature of the contact angle and its hysteresis, and, (ii) pinning of contact lines.
The scale at which we look at the three-phase contact line manifests a dramatic change in our perception of the contact angle, which the liquid-gas interface makes with the liquid-solid interface. This is explained in Figure 2 Bonn et al., 2009]. From a macroscopic 'distal' viewpoint, wherein the length scale is of the order of capillary radius (or the droplet radius, in case of spreading droplets), the angle which the liquid-gas interface makes with the liquid-solid interface is termed as the 'apparent' contact angle ('apparent dynamic contact angle', if the contact line is moving). As we move to smaller scales from the 'outer' macroscopic regions, and approach towards the contact line, we pass through an 'intermediate' or proximal region where there is interplay of viscous and surface tension forces. This region is characterized by a rapid change in the interface slope [Huh and Scriven, 1971;de Gennes et al., 1990]. A point of inflection is evident in this region where the interface profile turns from concave to convex. Moving further closer to the contact line, we eventually reach the 'inner region'. This region is of the order of nanometers and characterized by interaction of molecular long-ranged forces. The interface here is constantly under thermal and diffusive fluctuations and it is not possible to uniquely determine the contact angle [Thompson and Robbins, 1989;Thompson and Troian, 1997]. Thus, for all practical purposes, we shall use and refer to the 'apparent' contact angle (and drop repetitive use of the word 'apparent'), unless stated otherwise. In general, for experiments involving submersion of a solid plate in free surface of a liquid with very slow moving constant velocity, two distinct values of the apparent contact angles are seen [Young and Davis, 1987]. These angles are known as apparent advancing θ a and receding θ r angles respectively, depending on the direction of motion, as shown in Figure 3-a. The difference between the advancing and receding contact angles is known as contact angle hysteresis. It is generally acknowledged to be a consequence of three factors: (i) surface in-homogeneity, (ii) surface roughness, and (iii) impurities on the surface. This can be more clearly understood by considering Figure 3-b which shows the sketch of the wetting behavior on a substrate with a continuously varying physical topography (chemical texture is homogeneous) and continuously varying chemical texture (physical topography is uniform). Figure 3-b (i) shows that if the drop size is smaller than the topography, then the drop shape is not affected by the topography. But if the drop size is bigger than the physical topography then the global shape of the drop will be affected due the deformation of the three phase contact line due to the topography. Similarly, for a substrate with a gradient in wettability, as shown in Figure 3-b-(ii), the drop shape will be deformed due to the variation of the contact angles around the three-phase contact line of the drop. If, on the other hand, the substrate has sharp topography or wettability pattern then situation is very different, as shown in Figure 3-c. At such sharp discontinuities, Young's equation (Eq. (4)) becomes ill-defined. As a result, the three-phase contact line gets locally immobilized which is known as the pinning of contact line. The contact angle at the discontinuous boundary can have any value in between the smaller angle θon the more lyophilic part and the larger value θ + on the more lyophobic part. The contact line then becomes immobilized and its position gets fixed to the line of discontinuity as long as the contact angle falls into the free range between θand θ + [Krishnacharya, 2007].

Three-phase contact lines: Dynamics
Most classical work on dynamics of spreading is confined to two major arrangements (i) forced flow in a capillary tube and (ii) spontaneous spreading of a droplet on a smooth horizontal surface [Hocking, 1976;de Gennes, 1985;Bonn et al. 2009]. It is also important to note that these classical theories of contact line dynamics are usually restricted to cases where there is (a) no hysteresis of contact angles (b) gravitational effects are negligible, Bo → 0, and, (c) practically no inertia effects, i.e. Weber number tending to zero. Experiments with these constrains require patient efforts; a general review on the significance of the experiments have been given by Dussan [1979] and de Gennes [1985] and more recently, an exhaustive overview has been provided by Bonn et al. [2009]. Another limitation of these experiments, which have tried to isolate the interplay of purely viscous and surface tension forces, is that most experiments have been performed on dry surfaces for fluids with S i > 0. Further more, for the experiments with forced flow in capillary tubes, the effect of capillary diameter for a given liquid, i.e., its Bond Number, on the dynamics has not been extensively studied although there are clear indications that such dependence exists [Ngan and Dussan, 1982]. In case of an operating PHP none of these assumptions may be valid and most hydrodynamics is under pre-wetted conditions. It is well known that Bond number also affects the thermal performance of the PHP system. Nevertheless, the available theory of spreading provides fundamental insights into the phenomena under consideration. This should help improve closed form equations for the global model of the PHP.
The contact angle formed between a flowing liquid meniscus or a front and a solid surface is achieved by a balance between the capillary forces and the viscous forces. These two forces can be conveniently scaled by the Capillary number defined as: The viscous drag force on the liquid slug inside a capillary tube and the capillary force on the slug, scale respectively as: Thus, we deduce the scaling ratio of the drag and the capillary force as: It is clear that increasing liquid slug length L s increases the net drag force vis-à-vis the capillary forces. Multiple bubbles in the PHP system increase the capillary contribution to the total pressure drop, as will be seen later in Section 4.1.1.
Estimation of capillary forces require explicit knowledge of the dynamic apparent contact angles, both advancing and receding, when a unit-cell (Taylor bubble along with the adjoining liquid slug) train is passing through a capillary tube. One of the earliest mention of a 'universal' law between the velocity of spreading and the apparent advancing contact angle is highlighted by the work of Rose and Heins [1962] who could correlate cos(θ a ) with the velocity U of the contact line. The range of experiments was very narrow. Later, for capillary tubes which are already pre-wetted with a thin layer of liquid, the following equation was found to be more suitable: The log (tan(θ a )) vs. log (Ca) has a slope of 1/3 and it is also obvious that Eq. (12) cannot be valid for values of θ a approaching 90º [Schwartz and Tejada, 1972].
The geometry chosen by Hoffman [1975] was also a capillary tube of diameter ~ 2.0 mm with variation of five orders of magnitude of the velocity of the meniscus. Two sets of experiments were conducted: (a) with silicon oils having complete wetting, i.e. S i > 0 and (b) with liquids having non-zero static contact angle. For the former range of experiments and Ca → 0, Hoffman [1975] found a rather universal relation between θ a and Ca, given by: where q1 = 3.0 ±0.5. In the experiments at low Ca, θ a first increases as (Ca) 1/3 but ultimately θ a → π as Ca →1. For the second series of experiments also, the data could be represented by the form: While there are many assumptions and underlying limitations in the above 'laws', it is remarkable to note that these relations hold for liquids having different spreading coefficient S i (S i > 0), i.e., the magnitude of spreading coefficient apparently does not have any influence. This has been well documented, both experimentally and theoretically [de Gennes, 1985]. A systematic theoretical and experimental study in a dry capillary was conducted by Tanner [1979] and a simple cubic relation between the apparent dynamic advancing contact angle and Ca, θ a ~ K 1 (Ca) 1/3 was established. Optical observations of spreading silicone oil drops by Tanner showed that the droplet profile had an inflection near the edge, consistent with the existence of a pressure gradient arising from surface tension, which caused the spreading motion. The edge profile was predictable and the analysis showed that the edge velocity was approximately proportional to the cube of the slope at the inflection. The relationship given by Eq. (13) is therefore referred to as Tanner's Law or Hoffman-Tanner Law. Figure 4 shows the typical variation of the Tanner Law curve by taking K 1 = 1.92 [Kavehpour et al., 2003]. Other improved forms, by taking a somewhat stronger power law with the exponent of Ca = 0.39 and constant = 3.4 and more detailed treatment by incorporating a weak logarithmic dependency of the constant K 1 which gives θ a ~ K 3 (Ca·ln(K 4 Ca 2/3 ) 1/3 is also depicted along with the experimental data of Kavehpour et al. [2003]. All these variants of Hoffman-Tanner law are valid for low Ca ranging from 10 -6 -10 -2 , under the assumptions discussed earlier in this section.
For dynamic systems with non-zero static contact angle, the Hoffman-Tanner Law can be slightly modified and extended as: If Ca is small, linearization of Eq. (15) leads to [Berthier, 2008]: These equations inform us that the advancing contact angle is larger than the static contact angle and a receding contact angle is smaller than the static contact angle. This observation has been frequently reported qualitatively in the context of PHPs also.

TWO-PHASE FLOW IN CAPILLARIES
With this background of basic wettability issues and formation of contact lines as discussed in the previous section, let us now focus our attention on two-phase flows in capillary tubes.
Transport mechanisms of heat, momentum and species under twophase flow conditions in mini/micro systems are greatly affected by the distribution of phases in the channel (or the flow patterns that emerge). One of the important multi-phase flow pattern frequently occurring in two-phase flow systems is the motion of a gas/vapor bubble entrapped between liquid slugs as was shown earlier with reference to PHPs (in Figure 1-c), usually termed as slug flow. Slug flow belongs to a class of intermittent flows that has very distinctive features. There is an inherent random/statistical component for such a flow because of intermittent phase reversal which results in fluctuations of pressure and velocity. The complications arise due to the unique flow structure which is neither periodic in space nor in time [Ajaev and Homsy, 2006]. The two-phases may be composed of a either gas-liquid two-component system or a single component, vapor-liquid system.
Capillary slug flow or Taylor bubble flow is one of the sub-classes of conventional slug flow; as noted, this special pattern appears when surface tension dominates over gravitational body force, as discussed in Section 2. Due to this, capillary slug flows exist in horizontal as well as in vertical orientations. Also, flow is essentially laminar and predominantly viscous; the liquid plugs are free of smaller bubbles and breakage and coalescence of air slugs is virtually absent. This is often seen in the adiabatic section of PHPs. The understanding of species transport under such a flow configuration is a challenging problem. Other than PHPs, in recent times, the attention to two-phase Taylor/slug flows has increased due to the development of mini/micro scale systems in diverse branches ranging from bio-medical, bio-chemical to thermal management of electronics, water management of fuel cells, micro-two-phase heat exchangers and reactors, nuclear rod bundles, DNA separation and analysis, digital microfluidics, microscale mass transfer processes, lab-on-chips, micro-fluidic devices, loop heat pipes, etc.; in all these cases of emerging technologies, Taylor bubble flow is the dominant flow pattern. [for example, see Ghiaasiaan and Abdel-Khalik, 2001;Devesenathipathy, 2003;Spernjak et al., 2007;Steijn et al., 2008;Moharana et al. 2010].

Taylor bubble flows
Taylor slug flow occurs over a wide range of parameters, for all pipe inclinations and for a wide range of gas and liquid flow rates. Taylor slug flow conditions are typically characterized by a sequence of long bubbles which are trapped in between liquid plugs. The diametrical size of these bubbles is nearly comparable to the pipe diameter while their axial length scale can be several times the pipe diameter. A thin liquid film usually always separates the bubbles from the channel wall. The thickness of this film varies from millimeter scale to micro scale depending upon the dimensions, geometry, flow velocity, orientation of the channel and the thermophysical properties of the fluid and channel used. The intermittent liquid plugs may also have small diameter bubbles (much smaller than the pipe diameter) entrapped inside them, as seen in the third image on Figure 5; this tendency reduces as the Bond number decreases. Taylor bubbles are characterized by strong geometric constraints. The definitive pioneering work, for configurations involving circular conventional channels, which has spawned a large literature, is that of Bretherton [1961], who observed that a Taylor bubble does not rise spontaneously in a water filled vertical capillary under the effect of gravity for Bo < 1.835; this value can also be used to define small channels, as noted in Section 2. Contemporary researchers of Bretherton, who also addressed this problem were Harmathy [1960] and White and Beardmore [1962].
The presence of the film that separates the bubble from the wall means that the bubble velocity is not equal to the liquid one [Fabre and Lińe, 1992]. The presence of bubbles in front and at the back of the slugs, modifies the flow field in the liquid slug compared with singlephase flow and toroidal vortices extending the length of the slug can form [Thulasidas et al., 1995[Thulasidas et al., , 1997. The recirculation patterns within the liquid slugs improve heat and mass transfer from liquid to wall and interfacial mass transfer from gas/vapor to liquid. Thus, Taylor flow offers many advantages for carrying out heat mass and momentum transfer operations compared with other patterns and to single-phase laminar flow. The film surrounding the bubbles is the only means of communication between two successive slugs, and in the majority of cases, its thickness is only a very small percentage of the tube diameter [Moriyama and Inoue, 1992]. For practical estimation of species transfer coefficients, properties such as thickness of film that surrounds the bubbles, bubble shape and velocity, bubble and slug length, flow patterns in the liquid slug, and pressure drop, are some of the primary parameters of interest [Devesenathipathy, 2003;Nogueira et al., 2003Nogueira et al., , 2006Angeli and Gavriilidis, 2008;Steijn et al., 2007Steijn et al., , 2008.
The reduction in the length scale poses challenging tasks to model the Taylor bubble flows in mini/micro systems. The flow and heat transfer are influenced by the dynamics of a very small group of bubbles or typically, isolated slug bubbles, rather than the averaged behavior of a large population of bubbles. This, in principle, improves the prospects of developing mechanistic models. Such development needs to be facilitated by localized experimental observations of slugbubble systems with synchronized measurements of the resulting fluctuations in local conditions such as temperature, pressure and wall heat flux [King et al., 2007]. Mass transfer characteristics are also affected by the local hydrodynamic properties of the flow [Berčič and Pintar, 1997]. Knowledge of local hydrodynamic characteristics of the unit cell during Taylor slug/bubble flow is vital for complete understanding of the behavior and improving the performance of micro-thermofluidic and micro-chemical systems that operate in this regime. increasing superficial velocity [Khandekar, 2004].
The main forces acting on slug-bubble system are due to surface tension, viscosity, inertia, gravity, applied tangential shear stress and the force associated with the disjoining pressure at the molecular level. This gives rise to primary non-dimensional numbers of interest, i.e. Reynolds number, Capillary number, Bond number and Weber number plus geometric aspect ratios, if any. Marangoni number (M) and dimensionless Hamaker constant (ε) are also important depending on the domain of interest. When dealing with additional time scales in the system, for example arising due to externally imposed oscillations, Strouhal or Womersley number need to be considered for appropriate scaling. Additional non-dimensional numbers appearing for nonisothermal cases, without phase-change, are Peclet number (Re·Pr) and Nusselt number. If the axial conduction in the substrate is important under given boundary conditions, the Biot Number of the wall material also needs to be considered. The issues become more complicated with phase-change phenomena of evaporation and condensation occurring simultaneously in the Taylor bubble system.
It has also been found that the surface energy (surface wetting capability) of the wall may have great impact on the interfacial film dynamics when a confined bubble is flowing along with liquid in a mini/microchannel [Aussilous and Quéré, 2000]. This indicates that the friction of Taylor slug flows may strongly depend on the wall property, effective dynamic contact angles and hysteresis in a mini/microchannel [Taha and Cui, 2004;Berthier, 2008]. Modeling multi-phase fluid flow at the micro-scale with the assumption of simple disjoining pressure and liquid surface tension may not be accurate [Akbar and Ghiaasiaan, 2006]. The former will depend on the interactions between the fluid, gas/vapor, and the surface properties of the wall. Furthermore, due to the complexities of micro scale thermal capillary effects and the aspect ratio of the mini/micro channel, the relationship between the interfacial liquid film thickness, velocity profiles, friction, plug/slug menisci shape, dynamic contact angle hysteresis, is also not well understood [Chamarthy, 2008].
While in many systems capillary flows are uni-directional in nature or are employed in this way, there are many instances where the flow is either oscillating or pulsating in nature, as in the case of PHPs. These oscillations/pulsations may be externally controlled, thermally driven or alternatively an effect of the dynamic instabilities which are inherent part of two-phase boiling/condensation systems. The review by Angeli and Gavriilidis [2008] clearly suggests that bulk of the work has been carried out for unidirectional fully developed 'steady-state' Taylor bubble flows in capillaries. Two-phase flow in non-circular channels is also not commonly addressed in literature. There is scarcely any work done on oscillatory Taylor bubble flows.
As was noted earlier, exhaustive mathematical modeling of the phenomena occurring inside a PHP is yet an unaccomplished task. The models which are available for pressure drop and heat transfer coefficients are indeed very much simplified; they not only ignore the complex physics of the 'unit-cell', the singular characteristics of oscillating hydrodynamics of Taylor bubbles also remain completely unaddressed.  and Dobson [2004Dobson [ , 2005 studied the governing mechanism of the PHP using highly simplistic models. They studied a U-shaped miniature tube (i.e. single bend PHP) with a single liquid plug or a vapor bubble. The evaporation/condensation rate was assumed to be proportional to the difference of the temperatures of the vapor and the walls in contact with it. In the vapor bubble evolution equation, this leads to terms analogous to those of sensible heat transfer between the vapor and the tube walls. The approach by  has been extended by Shaffi et al. [2002] to model both looped and unlooped PHPs with multiple vapor bubbles, liquid plugs and tube bends. This model has been used later by another team [Sakulchangsatjatai et al., 2004], also for multi-bubble PHP modeling. It is well known from general considerations of thermal resistance that during the meniscus evaporation, an important contribution to the heat and mass transfer comes from thin liquid films that may cover the interior of the capillary (see Section 5.3). The local two-phase equilibrium exists at the interface of microscopically thin films so that the interface is at saturation or very close to it. This effect was completely neglected in the above mentioned modeling approaches. Dobson [2004Dobson [ , 2005 has included a film in his single bubble model. However, the film mass exchange in his model was not related to the liquid-vapor equilibrium and the mass exchange was proportional to the difference of temperatures of the vapor and the wall, just like in the other works that did not treat the films at all. The single-bubble model of Faghri [2002, 2003] has taken a step forward by rigorously showing that most part of heat and mass exchange occurs via the films in the PHP. The shape of the curved meniscus including the film has been calculated.
Globally, the existing models describe oscillations of only small amplitudes. During these oscillations, the Taylor bubble trains (or their menisci) are located almost all the time in either the evaporator or the condenser section. This contradicts most experimental results of PHPs where strong amplitude meniscus oscillations are observed, i.e., the pulsating flow is superimposed with strong unidirectional flows. During such large amplitude oscillations, the Taylor bubble menisci may sweep both the entire length scale of the pulsating heat pipe. It is imperative that for comprehensive understanding and subsequent modeling ventures of PHPs, there is a need to systematically study the thermohydrodynamics and transport mechanism of oscillating/pulsating Taylor bubble flows in mini/micro channels.

Parameters of interest:
The main parameters which impact the transport behavior in the context of unidirectional flow of Taylor bubbles flowing in circular geometries are:

(a) Liquid film thickness
It is known from quite sometime now that, when a gas bubble displaces a wetting fluid which is filled in a capillary, a liquid film is deposited on the capillary between the gas bubble and the inside wall [Fairbrother and Stubbs, 1935]. Early on it was known that the film thickness follows the scaling (δ/R) ~ Ca 1/2 for Ca ranging from 10 -5 to 10 -1 . Seminal analytical work, with some simplified assumptions by Bretherton [1961] established the scaling as (δ/R) ~ Ca 2/3 , the theory being valid for creeping flows at low Ca ranging from 10 -3 to 10 -2 and bubble Weber number <<1. Although, ironically the theory is based on low Ca assumption, it under-predicts the film thickness at the lower range of Ca. This discrepancy has been explained by various mechanistic models and arguments [Goldsmith and Mason, 1963;Schwartz et al. 1986], the most convincing being the effect of variable surface tension in the film and the bubble cap region, which in turn increases the film thickness at low Ca [Ratulowski and Chang, 1989]. Many experiments have confirmed the scalings, some results being close to Bretherton's theory while some more closely following that proposed by Fairbrother and Stubbs [1935]. Various numerical studies, especially at high Ca, have been done not only to predict the bubble film thickness but also the Taylor bubble shape [Edvinsson and Irandoust, 1996;Giavedoni and Saita, 1997;Heil 2001]. The effect of inertia (average flow Re) is also felt on the bubble shape and film thickness; it is non-monotonic but not profound. With increasing Re the film thickness first decreases slightly and then somewhat increases. In addition, ripples and capillary waves can be seen on the interface. Moreover, the Froude number also affects the bubble shape as the significance of body forces increase; this is especially true while comparing upward and downward Taylor bubble flows. Increasing Ca also leads the shape of the back of the Taylor bubble to change from convex to concave [Taha and Cui, 2006]. Figure 6 shows the variation of film thickness with Capillary number by using four commonly available correlations, as listed therein, along with their validity range. The effect of high Ca region is well captured by the correlations of Aussilous and Quéré [2000], which at low Ca, overlaps with the Bretherton's theory. In the context of a PHP, looking into the experiments available so far with common working fluids like water, ethanol, acetone, and common refrigerants etc., typical Ca will be ranging from 10 -4 to 10 -1 [Xu et al., 2005]. Systematic PHP visualization experiments are needed to quantify the film thickness regimes likely to be encountered in practice under oscillatory conditions. This is especially important in the adiabatic section which has substantial contribution in the overall pressure drop. In the evaporator and condenser sections of the PHP the situation is more complicated due to phase-change processes.

(b) Bubble velocity and shapes
It is well known that the bubble velocity is somewhat faster than the average liquid slug velocity due to the presence of the liquid film around it. Thus, there is a relative motion between the bubble and the adjoining liquid, which is characterized by: As this drift or relative velocity depends on the film thickness surrounding the bubble, it is obvious that the scaling for film thickness vs. Ca should, in general, be also applicable for m. Indeed, as per Fairbrother and Stubbs [1935], m scales as ~ Ca 1/2 while as per the Bretherton theory it scales as ~ Ca 2/3 . Again, as expected, when inertia forces become dominant and the liquid film thickness gets limited by the thickness of the boundary layer that develops in the liquid as it moves, the parameter m also approaches a constant value ~ 0.56-0.60 at large Capillary numbers [Taylor, 1961;Cox, 1964;Giavedoni and Saita, 1997]. Some comparative studies of vertical upward, vertical downward and horizontal or inclined tubes also indicate the effect of Bond number or Eötvös Number on the relative bubble velocity.
The possible streamline patterns in the liquid slug ahead of the gas bubble were first proposed in an early study by Taylor [1961]. For cases where Ca is low, m < 0.5, two possible flow patterns in the liquid slug were proposed (a) single stagnation points on the bubble tip and, (b) formation of a stagnation ring around the bubble cap with vortex pairs in the liquid slug; there is also a possibility of formation of one stagnation point on the bubble tip and a second stagnation point inside the liquid slug for a narrow range of 0.60<Ca<0.69. As Ca increases, the vortices in the liquid slug tend to disappear leading to bypass flows. Thus, at high Ca, for m > 0.5, only one stagnation point on the bubble tip was suggested, with complete bypass flow in the liquid slug [Cox, 1964]. More recent experimental and numerical studies by Thulasidas et al. [1997], Giavedoni and Saita [1999] and, Cui [2004, 2006] have confirmed these findings. Typical result of Taylor [1961] and Taha and Sui [2004] are given in Figure 7. Some PIV and visualization studies have also been undertaken which confirm the existence of toroidal vortices extending the length of the liquid slug [Thulasidas et al., 1995[Thulasidas et al., , 1997Kashid et al., 2005]. One of the parameters of interest in enhancement of heat and mass transfer is the typical velocity of the particles in the toroidal vortices. The inverse ratio of the time required by the liquid slug to travel its own length to the time it takes for a typical particle in the vortex to come from one end of the liquid plug to the other end, is defined as the dimensionless recirculation time. The ratio increases to infinity for very large Ca as there is complete bypass flow. For Ca < 0.2 this is found to be of the order of 2 [Thuslasidas et al. 1997]. Increasing inertia forces changes the shape of the bubble nose [Heil, 2001].

(c) Pressure drop across the unit cell
The nature of the Taylor bubble flow clearly suggests that flow pattern independent models, like the common Lockhart-Martinelli formulations [Lockhart and Martinelli, 1949] which are very successful on larger diameter channels, will not be successful. In fact, this observation was made by Lockhart and Martinelli themselves that their procedure is not applicable for slug flows having alternate liquid plugs and gas bubbles. Homogeneous flow models, wherein there is no interfacial slip, may be partially successful only under a very narrow range of slug flow operating conditions Kawahara et al. 2002, Chung andKawaji, 2004;Liu et al., 2005]. These common models fail when applied to Taylor bubble trains because of two primary reasons (i) dominance of surface tension in slug flows, which is neglected in these models (ii) additional and unique flow field attributes in terms of velocity patterns in the liquid slug and vapor bubbles, which is also ignored in conventional models. Chen et al. [2002] have shown that even after the inclusion of surface tension effects in terms of Bond number and Weber number, all the available experimental data sets on pressure drop in slug flows cannot be comprehensively predicted. This necessarily means that improved models should not only include the effect of surface tension but should also address to resolve the effect of flow fields in the unit cell on the overall pressure drop.
To develop a specific model applicable to Taylor bubble flows, the usual practice is to conveniently divide the net pressure drop in a unit cell into three components: (a) as proposed by Taylor [1961] and (b) as numerically computed by Taha and Cui [2004].
(i) Pressure drop in the liquid plug, (ii) Pressure drop around the ends of the bubble and, (iii) Pressure drop along the body of the bubble. The last contribution tends to zero or minimal when (i) gravity forces are small i.e., Bond number is small, and (ii) shear in the gas phase is neglected; viscosity and density of the gas phase is much smaller than the liquid phase. With these assumptions, the gas in the bubble is practically at constant pressure and for a stable constant liquid film thickness surrounding the bubble, its curvature is constant, necessitating that pressure drop cannot exist in the bubble body. The fact that, the terminal velocity of cylindrical bubbles rising vertically upwards in stagnant liquid is generally unaffected by the length of the bubble, also supports these assumptions. If we invoke this assumption of zero pressure drop in the bubble body, then the overall pressure drop curve may look something like what is depicted in Figure 8-a. An approximately triangular or saw-tooth alternating component of pressure drop gets superimposed on the average pressure gradient. It is interesting to note that this phenomenon, on its own, can excite oscillations and pressure perturbations in the flow [Wallis, 1969].  [Khandekar et al., 2002] (e) Explanation of the additional pressure drop due to capillarity.
In real-time engineering systems, contact angle hysteresis may be present wherein the advancing and receding angles of a typical bubble may be different, as noted earlier in Section 3.1. This gives rise to additional pressure drop due to capillary forces. This resistance is due to the capillary pressure difference between the plug front (in moving direction) and the plug end (against movement). If the bubble cap radii are assumed to be spherical in geometry, then the following simplified estimate of the capillary pressure gradient across a liquid-vapor interface in a tube of circular cross section, due to the Laplace equation, and because of the two different contact angles (advancing, θ a and receding, θ r ), can be formed: This scenario is depicted in Figure 8-b, c. If θ a is larger than π/2, there is a positive pressure drop associated with the advancing interface. If θ r is smaller than π/2, the receding front contributes positively to the pressure drop (Figure 8-c), and negatively in the opposite case (Figure 8-b). An example of this hysteresis occurring in a PHP operated by water and ethanol, respectively, is shown in Figure 8d while Figure 8-e explains the genesis of this additional pressure drop. The trailing slug curvature slips over a just wetted surface by the bulk liquid while the front meniscus has different conditions, as seen. This additional resistance is additive and gets amplified if several plugs are simultaneously located in a channel. These cumulative pressure gradients damp oscillations and therefore fluid/solid combinations having large dynamic contact angle hysteresis should be avoided for PHP application. Too many plugs in the capillary may even rapidly block the flow if the driving force is not sufficient. For N plugs, we have, where the index s stands for static contact angle; , s a θ and , s r θ are the two static contact angles, both being equal only if there is no static hysteresis, i.e., perfectly smooth and chemically homogeneous surface. As per the classical Bretherton's theory, for low Ca, the pressure drop across a gas bubble moving through a liquid filled circular channel is given by: As for the case of film thickness and relative bubble velocity, this formulation by Bretherton deviates from experimental results. Inertial effects also lead to a drift away from the Bretherton theory as explicitly demonstrated by Westborg and Hassager [1989], Giavedoni and Saita [1997], Heil [2001] and Fujioka and Grotberg [2005]. As the Capillary number increases, the film surrounding the bubble grows and subsequently occupies a significant portion of the tube around the bubble. At higher flow velocities there may be considerable interfacial resistance; flow fields develop in the gas bubble too, which ensue a non-zero pressure drop.
In a recent exhaustive work, Kreutzer [2005] has modeled the Taylor bubble flows extensively for a range of flow conditions. The general line of argument followed is that the Hagen-Poiseuille flow in the liquid slugs gets disturbed by the intermittent bubbles. Thus, the pressure drop calculation done as per the standard Hagen-Poiseuille formulation needs amendments to capture the additional pressure drop created by the presence of the bubble caps. For a given Ca, the flow disturbances in the liquid plug created by the bubble caps is confined to a small distance in its vicinity, and is therefore independent of the length of the liquid slug. The disturbance length scale is also a function of the thermophysical properties of the fluid/wall combination, the Bond number, and the average Reynolds number of the flow. Thus, the correction term in the fanning friction factor (f fanning ) was chosen to be inversely proportional to the length of the liquid slug. The modified fanning friction factor for liquid slug in the unit-cell suggested by Kreutzer [2005] is given by: In this equation, A and B are constants and the modified Reynolds number and Capillary numbers are defined on the basis of sum of the gas and liquid superficial velocities, j l and j g . While numerical and experimental values of constant B was reported to be equal to 0.33, the experimental value of constant A was higher (= 0.17) than the numerically predicted value (=0.07). In this formulation, Re * and Ca * are given by: To calculate the aggregate pressure drop, the statistical measure of the typical liquid slug lengths in a section of the tube under study needs to be known. All the data used for arriving at Eq. (23) were for Re b > 150. For lower bubble Reynolds numbers another modification of Eq. (23) is provided recently by Warnier et al. [2010].
All mathematical models for pulsating heat pipes use extremely simplistic pressure drop correlations, typically the Hagen-Poiseuille model. Further, the statistical measure of bubble/liquid slug lengths in the adiabatic section of the pulsating heat pipe is not yet fully resolved vis-à-vis the operating boundary conditions. Unless these issues are studied, proper pressure drop estimation in the adiabatic section of the device cannot be achieved. The fact that the Taylor bubble flow is also oscillating, further complicates the estimation. Very little information is available on oscillating Taylor bubble flows, which we review in the next section.

OSCILLATING/PULSATING MENISCI AND BUBBLES INSIDE CAPILLARIES
The subject of oscillating/pulsating motion of liquid slugs inside mini or micro-channel/pipe/capillary/duct has received attention over the last two decades, both due to the large number of practical applications in which it appears and the interesting scientific challenges it poses. There are many engineering systems where oscillating/ pulsating two-phase capillary flows occur. The oscillations may be externally controlled, thermally driven or alternatively an effect of the dynamic instabilities which are inherent part of two-phase boiling/condensation systems [e.g., see Spernjak et al., 2007, Tadrist, 2007. In unidirectional flow system, the only relevant time scale is the momentum diffusion time. Here, in contrast, the flow is characterized by a second imposed time scale that is due to the oscillatory pressure gradient; the form of the resulting flow depends on the ratio of the diffusion time to the imposed time scale, i.e., Strouhal number (St). The magnitude of St determines the importance of acceleration effect in the fluid relative to the viscous effects (diffusion of momentum). In this role, it is evident that St can be considered as a Reynolds number, based on a characteristic 'velocity' ~ f·D. In addition, during oscillatory flows, the velocity of the meniscus/bubble continuously changes, resulting in dynamic variation of Re, Ca, We etc.
Drawing an analogy from single-phase oscillating systems, heat, mass and momentum transfer coefficients are expected to be higher in case of oscillating two-phase systems. While a number of investigators have looked into the relevant issues for fully developed unidirectional steady state Taylor bubble flows [Ho and Leal, 1975;Martinez and Udell, 1990;Borhan and Pallinti, 1998], information on oscillating Taylor bubble flows is scarce [Graham and Higdon, 2000]. Young and Davis [1987] were among the first to study oscillatory contact line motion and they simplified Dussan's model [1979] for unidirectional creeping flow. Their analysis considered slow quasi-static oscillatory contact line motions of a solid plate being submerged in free surface of a liquid. Analytical solution was given when contact angle hysteresis was neglected. In studying a similar capillary gravity wave problem, Miles [1990] argued that the contact angle can be proportional to contact line velocity, but there will be a phase-lag in oscillatory motions. From a mathematical point of view, the description of the correct boundary condition at the contact line becomes an issue. The presence of capillarity adds an extra term to the dynamic-free surface boundary condition which is proportional to the free surface curvature and thus increases the order of the dynamic-free surface boundary condition [Perlin et al., 2004].

Oscillatory slug flows under adiabatic conditions
Oscillatory two-phase flows under 'adiabatic' conditions can be obtained by controlling the flow by external pressure variations.
Depending on the wettability of the fluid on the solid wall, in conjunction with the applied external surface and body forces, there can be three types of oscillatory contact line motion of the liquid-gas/vapor meniscus confined in a tube, as detailed in Figure 9, viz., (a) No contact line sticking -the average velocity of contact line scales with the average bulk liquid velocity, (b) Partial sticking -there is a phase-lag in the average contact line velocity and the bulk fluid velocity, and (c) Full sticking, the contact line is nearly motionless while the bulk liquid moves/ oscillates. Further, in cases (a) and (b), there can also be contact angle hysteresis in advancing and receding strokes (2D axisymmetric representation has been shown here). The corresponding schematic representation of the oscillatory movement of the liquid slug is shown in Figure 9-d, e, f. Recent exploratory studies by Xishi [2005, 2006] focus on finding the effect of oscillating frequency on the interfacial film thickness around Taylor bubbles. As we have seen in Section 4.1.1, many relations are available for film thickness of a Taylor bubble, under unidirectional quasi-steady flow conditions, suggesting that the film thickness increases with increase in the velocity. However, during oscillatory conditions, the bulk flow velocity varies from 0 to U max in a prescribed manner. Inertia forces may also be dominant depending on the Weber number. Xishi [2005, 2006] used a novel laser based optical fringe scattering technique to minimize the 'Mirage' effect due to the cylindrical nature of the tube cross section in which the flow takes place. The schematic of the experimental set-up is shown in Figure 10-a. The preliminary result for the change in the mean film thickness with the oscillating frequency is plotted in Figure 10-b. The film thickness increased rapidly with increase in the frequency from 0 Hz to 4.0 Hz; thereafter an increase from 4 to 8 Hz only increased the film thickness by 0.1 micrometer. Further increase of frequency introduced some periodicity in film thickness. Looking at the experimental results, they suggest changes in the steady flow model of film thickness to suit the oscillatory conditions. Shekhawat et al. [2009] have studied the interfacial contact line behavior of the single oscillating meniscus formed between a long liquid slug and air, inside a square capillary tube (2.0 mm x 2.0 mm). An eccentric cam follower system, as shown in Figure 11-a, was fabricated to provide simple harmonic oscillations of fluid meniscus (water and silicon oil). Dynamic apparent contact angle measurements were carried out for water at two oscillating frequencies, 0.25 Hz and 0.50 Hz using high speed videography. Contact angle of water at static condition varied from 19.8 o to 21.8 o , with average value of 21 o . Figure 11-c shows the images of the oscillating meniscus at respective locations marked on the displacement-time graph, i.e. Figure 11-b. Points a, b, c, and d show the advancing meniscus and e, f, g, and h show the receding meniscus respectively. Points d and e are at the top dead end of the stroke, and points h and a are at bottom dead end of the stroke respectively. A clear change in the curvature of the meniscus as it oscillates was noted; to quantify this change the contact angle was measured by magnified videography. The relevant parameters of the experiment are noted in Figure 11-d.  Xishi [2005, 2006] used for measurements of liquid film thickness in an oscillating vapor plug and, (b) mean film thickness as a function of oscillating frequency.  [Shekhawat et al., 2009]. Figure 12-a shows the variation of contact angle along the stroke length at an oscillating frequency of 0.5 Hz and 0.25 Hz, static contact angle line is set at 21°. It was seen that the advancing contact angle was more than the receding contact angle. It can be seen that the difference in the advancing and receding contact angle at 0.5 Hz is more than that at 0.25 Hz. This effect was attributed to the inertial and viscous forces as indicated by higher Re and Ca. Figure 12-b represents the velocitytime graph super-imposed on the dynamic contact angle. Contact angle variation were quite in-line with the temporal velocity fluctuations.
Visualization of silicon oil meniscus in Figure 13 revealed a major difference between the behaviors of the two fluids. It was observed that the contact point of the silicon oil meniscus stuck to the extreme end of the stroke, and due to this a continuous film was formed when the bulk fluid was pulled away from the top dead center position. Also, there was a drastic difference in the curvature of meniscus during advancing and receding stroke (see Figure 13-a and 11-c for comparison). Figure 13-a shows the images of the meniscus oscillating corresponding to the marked positions, at an oscillating frequency of 0.75 Hz. The thickness of the film formed was also recorded at eight different oscillating frequencies, as noted in Figure 13-b, along with the corresponding values of dimensionless numbers. Figure 13-c shows the variation of average film thickness (see inset: at station A, which is at the mid point of the stroke length) with the change in frequency. It was observed that the film thickness increased with increase in oscillating frequency, the rate of increase, decreasing with increasing frequency. The study concluded that an increase in the oscillating frequency increased the difference in the advancing angle and receding angle of the meniscus. There was considerable difference in the hydrodynamics of oscillating menisci of silicon oil and water respectively.
During forced mechanical simple harmonic oscillations of a single meniscus of silicon oil-air interface by Tripathi et al. [2010] inside a square capillary tube of D hyd = 3.0 mm, it was also seen that the meniscus contact line got totally immobile at the extreme end of the stroke during oscillations. Figure 14-a, b shows the different stages of the oscillating motion of single meniscus at 0.5 Hz, inside the capillary kept in horizontal and vertical positions respectively. The contact line got fixed relative to the glass tube and did not move at all, in spite of the sinusoidal motion given to the fluid alongside. Due to this, a continuous liquid film was formed along the stroke length. The dominance of adhesive forces in the case of silicon oil was clearly visible. The static contact angle of silicon oil on clean glass was found to be ~14.2°. Re D in this case changed from 0 to 3.3.
For oscillations in horizontal tube orientation the effect of Bond number was clearly seen; the meniscus shape was not symmetrical, the length as well as the thickness of the bottom liquid film was always greater than at the top side. Also, the hysteresis of dynamic contact angle in the forward and the return stroke were clearly visible. In vertical orientation, there was considerable change of the radius of curvature in the upward and the downward strokes. The contact line completely stuck to the glass tube inner wall. This orientation provided symmetrical body forces on the meniscus geometry.    [Tripathi et al., 2010].
Experiments with water were also carried out at oscillating frequency of 0.5 Hz. In this case, contrary to the silicon oil case, the contact line moved along with the bulk fluid/meniscus, albeit with some phase shift as clearly seen in Figure 14-c, d. Slip type motion was observed; the contact line moved freely along the tube with a time varying dynamic contact angle. Figure 14-c shows the oscillations of single water meniscus inside a horizontal capillary. When it was in the receding stroke (t = 0 to t = 1 s), the effect of gravity, even at a relatively small Bond number of 1.2, was quite visible; the meniscus shape was not quite symmetrical. The meniscus showed a relatively small radius of curvature in the receding stroke. During its advancing stroke (t = 1 to t = 2 s) the meniscus shape became flattened. Figure 14e shows the close view of silicon-oil and water menisci, for both advancing and receding interfaces, highlighting clear differences in the behavior of the two fluids.
Recent studies by Lips and Bonjour [2007] and Lips et al. [2010] on adiabatic oscillations of menisci also brought to the fore the importance of dynamic contact angle in the Taylor bubble flows and the dissymmetry between the advancing and receding contact angles. Figure 15-a depicts the schematic of their test bench. It consisted of two tanks connected by a 2.4 mm ID and 4 mm OD glass capillary tube of length 300 mm. A speed controlled motor was used to impose a quasisinusoidal volume variation on the upstream tank (variable pressure reservoir) through an eccentric wheel, and a flexible membrane mounted with a return spring. The downstream tank (constant pressure reservoir) was filled with liquid pentane in equilibrium with its vapor, whereas the upstream tank contained only pentane vapor. All experiments were done at saturation conditions. The amplitude of the pressure oscillation (~ 0.1 bar) and the frequency (~ 4 Hz) were representative of their typical values experimentally obtained on a branch of full-scale PHPs [Bensalem, 2008;Khandekar et al., 2009]. Figure 15-b illustrates a set of successive slug photographs at 20 ms interval. Enlarged views (i), (ii) and (iii) alongside show a strong dissymmetry of left and right interfaces: the curvature radius of advancing menisci (right interface in (i) and left interface in (iii)) was smaller than that of respective receding menisci representing of the pressure difference between both sides of the liquid slug, as explained earlier in Section 4.1. In case (ii), the interface velocity was equal to zero and the liquid slug was then symmetric. Figure 15-c shows one representative time evolution of the curvature radii of the two interfaces which was obtained by image processing at a liquid slug velocity oscillating between -0.2 and +0.3 m/s (negative velocity relates to a receding meniscus, while a positive velocity relates to an advancing one); eight periods of oscillation were recorded and superposed, the average being represented by the dotted line.
A large number of experiments were done for various conditions: the liquid slug length ranged from 2.4 to 15.2 mm, the amplitude of the velocity oscillations could reach up to 0.65 m/s, and the frequency of the imposed sinusoidal evolution of the pressure ranged from 1.8 to 6.1 Hz, which is believed to be the correct order of magnitude in PHPs. The experimental results were presented as the superposition of curves showing the inverse of the dimensionless meniscus curvature radius r* vs. meniscus velocity U (Figure 15-d). No clear effect of the slug length or frequency was detected; all the curves were roughly correlated by a linear fit as, where, r is the meniscus curvature radius and D is the tube diameter.
Although this correlation is only valid for their specific experiment and is not general enough, it can be used to assess the importance of the slug dissymmetry in pressure drop calculations. For that purpose, as an example a liquid slug of length L s = 1 cm whose maximum velocity is U = 0.5 m/s was demonstrated. Using Eq. (25), the advancing dimension-less radius of curvature comes out to be 1/r* = 0.465 (i.e. for U = 0.5 m/s), and that of the receding meniscus becomes 1.295 (i.e. for U = -0.5 m/s). Therefore, the capillary pressure difference reaches 10.9 Pa for pentane at 20°C (σ = 15.2 mN/m) in the 2.4 mm inner diameter tube used in the experiments. Besides, the frictional pressure drop for the liquid slug is calculated as: With ρ l = 625.7 kg/m 3 and μ l = 2.37.10 -4 Pa·s at 20°C, the frictional pressure drop for the considered 1 cm long slug reaches 13.5 Pa, which is of the same order of magnitude as the capillary pressure difference. This clearly demonstrates that the dynamic contact angle difference between receding and advancing menisci is vital for PHP pressure drop calculations (see inset of Figure 15-d). The capillary pressure difference will be as important as the frictional component in the momentum equation.

Oscillatory slug flows under diabatic conditions
Very few quantitative studies are available describing controlled oscillatory flows under diabatic conditions. Das et al. [2010Das et al. [ -a, 2010 have recently undertaken a parametric study of a two-phase oscillating flow in a capillary tube. Figure 16-a shows the schematic diagram of their experimental set up. This system represents the simplest version of a PHP; it consisted of a capillary tube of ID = 2 mm connected between two reservoirs. This tube was heated at one end (evaporator section = 15 cm; copper) and cooled at the other end (condenser section = 25 cm) by a transparent glass heat exchanger (for visualization) supplied by a constant temperature coolant, with an adiabatic length (variable = 1 cm, 2.5 cm and 5 cm respectively; glass) in between. An absolute pressure sensor was connected to one end of the evaporator. One end of the glass tube was connected to the evaporator and the other end was connected to a reservoir filled with saturated pentane. A heating coil was wound around this reservoir to fix its temperature, and thus experiments at different pressures could be performed. Visualization was performed and the oscillations of the liquid-vapor meniscus and the vapor pressure were observed. It was shown that such a system can develop instability that leads to thermally driven meniscus oscillations.
The oscillations in the system (i.e., its instability) appeared when the difference between the temperatures T e and T c exceeded a threshold value. At small T c (0 -10ºC), the meniscus did not move out of the condenser (towards the fluid reservoir) which, incidentally, was convenient for visualization. Figure 16-b shows meniscus oscillations and the vapor pressure P v , close to the instability threshold. Reservoir pressure P res was 90 kPa. The amplitude of pressure was about 15 kPa. The meniscus displacement curve was truncated from below because the adiabatic and evaporator sections were opaque and the meniscus displacement could not be measured inside them. The evolution was not exactly periodic; the average period was about 0.28 s. Note the intermittency in the height of the minima of the pressure curve: each second or third minimum is higher than the others. For T e = 65 ºC and T c = 0 ºC, the maximum meniscus displacement was about 39 cm as shown in Figure 16-c; the displacement variation is truncated not only from below but also from above because of the limited camera view. Unlike the previous case, the pressure variation is regular, amplitude being about 50 kPa and it was steeper at expansion than at compression. The period of oscillation was about 0.34 s. Kinks in the displacement plot correspond to the sudden relocation of the meniscus.
This happens because when the rate of condensation onto the moving liquid film in the condenser was high enough, the liquid film became unstable. Liquid film formed multiple bridges and coalesced to form secondary bubbles, relocating the meniscus of the main bubble. This could occur during the movement, both towards or away from the evaporator. The secondary bubble quickly disappeared by coalescing with the main bubble, which caused another kink. Figure 17-a show a series of images corresponding to Figure 16-b. Dark part of the capillary tube to the left is the vapor bubble and the bright part to the right is the liquid plug. Note that in both cases, at each oscillation, the meniscus had penetrated into the evaporator section. Figure 18 clearly shows the liquid film instability leading to formation of secondary bubbles by multiple bridges being formed by the liquid film.
The meniscus displacement in the transparent condenser and the pressure variation were analyzed both close to the instability threshold and far from it. Close to the instability threshold oscillations were less regular and had smaller amplitude than far from it. In an attempt to better understand the existing PHP models and their simplified approach, i.e., considering the vapor phase as an ideal gas with sensible heating/ superheating, was also scrutinized. While such a simplification allows some important parameters (in particular, the oscillation frequency) to be obtained analytically, Das et al. [2010-a] explicitly show that only small amplitude oscillations can be described by such a simplification. More coherent ('evaporation/condensation') model that allows the large amplitude oscillations (as observed experimentally in PHPs) has been proposed. The model differs from the existing PHP models (as discussed earlier Zhang and Faghri, 2003;Dobson, 2004, Dobson 2005) by virtue of the accounting of the two-phase equilibrium that occurs locally at the vapor-liquid interface. It also introduces and accounts for the time varying wetting films (see Figure 17-c) through which major part of the heat and mass transfer occurs [Höhmann and Stephan, 2002]. The model proposed by Das et al. [2010-a] succeeds in describing almost all observable qualitative features of the oscillations.

Implications on heat transfer
It is well known that during phase-change heat transfer across menisci, an important contribution comes from thin liquid films forming part of the meniscus that may cover the interior of the capillary [for example, refer Morris, 2003]. The total heat transfer resistance is composed of the diffusional resistance of the liquid film and evaporation at vapor-liquid interface. Höhmann and Stephan [2002] reported the temperature distribution at the micro-region between the adsorbed, non-evaporating film and the macroscopic region of a liquid meniscus. It was observed that the temperature in the micro-region drastically drops and a high heat transfer is obtained at that region; up to 50 % of evaporation takes place in this small region. The length of the micro-region depends on the wetting characteristic which is represented by the local apparent contact angle. Thus, the dynamic apparent contact angle is a very important process parameter. In the context of PHPs, the effect of the film evaporation has been clearly highlighted by models proposed by  and more recent comprehensive extension by Das et al. [2010-a], as noted earlier.
These models now need further refining by taking into account the exact shape/ dynamics of the interface for different working fluids of PHPs, as has been highlighted in this paper. Due to limitation of space, we have not explicitly reviewed the heat transfer implications of unidirectional and oscillating Taylor bubble flows; in the context of heat transfer in a PHP, the importance of this topic justifies the need of a separate dedicated review of the available knowledge. Coupling of heat transfer with the Taylor bubble hydrodynamics brings in other relevant phenomena such as Rayleigh-Taylor and Kelvin-Helmholtz instability, film dry-out, local non-equilibrium thermodynamics, meniscus deformation due to thermo-capillary waves, Marangoni convection, rapid compressibility effects, phase-change induced fluid motion, and heat transfer models due to phase change and resulting flow regimes [Faghri and Zhang, 2006;Carey, 2007].

SUMMARY AND CONCLUSIONS
Mathematical modeling of PHP is a contemporary problem, exhaustive and comprehensive solution strategy of which remains yet un-developed. Many fundamental phenomena need to be understood to develop a complete PHP model. Understanding the local hydrodynamics of oscillating Taylor bubble flows forms an important element in the building block of the mathematical description of the flow in a PHP. Eventually, the wettability of the liquid-solid combination has profound effects on the overall transport behavior. Needless to say, while the issue of hydrodynamics and resulting thermal behavior are intrinsically linked, in this article we have focused our attention only on reviewing the hydrodynamics of Taylor bubble flows. The fundamental understanding of the wettability and contact line motion are also reviewed. The related heat transfer implications necessarily require another detailed complementary review, which will be addressed in the near future.
Present review clearly highlights the fact that many hydrodynamic aspects of steady unidirectional Taylor bubble flow, such as the thickness of the liquid film that surrounds the bubbles, bubble velocity, bubble and slug lengths, mixing and flow circulation in the liquid slugs, as well as pressure drop during Taylor flow, are quite well understood. On the other hand, oscillating Taylor bubble flows require immediate attention. Unless the nuances of oscillatory confined bubbles in small capillaries is well discerned, the net pressure drop correlations for a PHP cannot be exhaustively constructed. The challenging issues are high inertia, effect of wettability, contribution of surface waves, film thickness dynamics and bubble breakage and merger due to flow instability.
There are not enough studies so far which explicitly highlight the qualitative and quantitative trends of the oscillating/pulsating flow. Preliminary available results for oscillating Taylor bubble flows, either generated through externally imposed mechanical pressure variation or by thermally inducing them, indicate the following trends: (i) Wettability of the working fluid-tube material drastically changes the shape of the meniscus in the advancing and retarding strokes. (ii) Most classical literature on the effect of Ca on the apparent contact angle is confined to unidirectional flow at low Ca values. Oscillatory flow velocities are not thoroughly addressed. (iii) Frequency of the flow oscillation affects the liquid film thickness which is formed around the Taylor bubble. (iv) In spite of the fact that the tube diameter satisfies the critical Bo crit criteria, gravity does affect the shape of the menisci. (v) Unless the thin film evaporation and phase-change thermodynamics is included in the mathematical model, the large scale oscillations of interfaces, as experimentally observed in PHPs, cannot be effectively captured. (vi) Inertia forces tend to create instabilities resulting in strong interfacial waves which can eventually result in bubble breakage and formation of secondary bubbles, etc. Fundamental hydrodynamic phenomena (and the resulting thermal trends) must be included in constructing global models for describing flow and heat transfer in a PHP. More studies to understand adiabatic as well as non-adiabatic phase-change induced oscillatory Taylor bubble flows are thus needed.