# Series Solutions for Marangoni Convection on a Vertical Surface

Authors: Gabriella Bognár; Krisztián Hriczó
Journal: Mathematical Problems in Engineering
 (2012)

---

## Abstract

The problem of steady, laminar, thermal Marangoni convection flow of Newtonian fluid over a flat surface is investigated. The boundary layer equations for the momentum and energy equations are transformed with the similarity solutions to ODEs to obtain the analytical approximate solutions. The analysis assumes that the temperature variation is a power law function of the location. The approximate solutions to the similarity equations are obtained by exponential series. The effects of the power law exponent and Prandtl number on the velocity and temperature profiles are presented.


---

## Body

## 1. Introduction

For many metallurgical and materials processing applications in space, it has been evidenced that if a free liquid surface (or liquid-liquid interface) is present, a surface tension gradient of this free surface will cause liquid movement at this free surface which alternatively drives fluid movements in the corresponding phase. This phenomenon is referred to as capillary motion. When a free liquid surface is present, the surface tension variation resulting from the temperature or concentration gradient along the surface can also induce motion in the fluid called solutal capillary, and thermocapillary motion, respectively.The study of liquid movement resulting from thermocapillarity (or so-called Marangoni) convection is very important for a liquid system either in microgravity or in normal gravity [1]. Under normal gravity, liquid movement is mainly driven by buoyancy force because of the temperature-dependent density while the liquid is exposed to a temperature gradient field. As the size of the liquid system decreases especially having the size decrease in the direction of gravity, the buoyancy effect begins to diminish and the Marangoni effect will then dominate the system as the main driving force for liquid interface movement. In the absence of gravity, Marangoni convection always plays a main role in the determination of the fluid movement because of varying liquid surface tension in a temperature gradient field regardless of the system’s size. It has significance in the processing of materials, especially in small-scale and low-gravity hydrodynamics [2] to investigate the interaction between the heat and mass transfer.Marangoni convection appears in many industrial processes and space technologies, for example, in the study of semiconductor crystallization processes, in crystal growth melts where the flow produces undesirable effects (see [3–6]) and it occurs around vapor bubbles during nucleation [5]. In several papers authors investigate Marangoni driven boundary layer flow in nanofluids. These fluids can tremendously enhance the heat transfer characteristics of the base fluid and have many industrial applications such as lubricants, heat exchangers, and coolants. Nanofluids are studied when different types of nanoparticles such as oxide ceramics, nitride ceramics, and aluminium, are used [7–9]. Marangoni flow has also significance in welding, semiconductor processing, and other fields of space science. Its mathematical model is studied in [7, 10, 11]. Marangoni boundary layers may occur along the liquid-liquid or liquid-gas interferences.These phenomena have been investigated for various substances in geometries with flat surfaces by similarity analysis (see [3, 7, 12–16]). Arafune and Hirata [3] presented a similarity analysis for just the velocity profile for Marangoni flow when the surface tension variation is linearly related to the surface position. Christopher and Wang [4] studied Prandtl number effects for Marangoni convection over flat surface and presented approximate analytical solutions for the temperature profile. They showed that the calculated temperature distribution in vapor bubble attached to a surface and in the liquid surrounding the bubble was primarily due to the heat transfer through the vapor rather than in liquid region and the temperature variation along the surface was not linear but could be described by a power-law function [12]. Using the similarity transformation, the governing system of nonlinear partial differential equations are transformed into a pair of similarity nonlinear ordinary differential equations, one for the stream function and one for the temperature. The velocity and temperature distributions can be given by numerically using the Runge-Kutta method [5, 7, 12, 17], analytical approximate solutions can be determined for these problems by using Adomian decomposition method and Padé technique [15, 16, 18, 19] or by power series method [20].This paper investigates a similarity analysis for Marangoni convection inducing flow over a flat surface due to an imposed temperature gradient. The analysis assumes that the temperature variation is a power law function of the location and the surface tension is assumed to depend on the temperature linearly.We first present the derivation of the equations and show how the boundary layer approximation leads to the two points boundary value problem and the similarity solutions. The new model, written in terms of stream function and temperature, consists in two strongly coupled ordinary differential equations. Its analytical approximate solutions are represented in terms of exponential series. The influence of various physical parameters on the flow and heat transfer characteristics are discussed.


## 2. Boundary Layer Equations

Consider the steady laminar boundary layer flow of a viscous Newtonian fluid over a flat surface in the presence of surface tension due to temperature gradient at the wall. Assuming that the surface is impermeable, the surface tension varies linearly with temperature and the interface temperature is a power-law function of the distance along the surface. The governing equations for two-dimensional Navier-Stokes and energy equations describing thermocapillary flows in a liquid layer of infinite extent are considered. The layer is bounded by a horizontal rigid plate from one side and opened from the other one. The rigid boundary is considered as thermally insulated. The physical properties of the liquid are assumed to be constant except the surface tension. This investigation are based on the balance laws of mass, momentum and energy. In the domainx>0,y>0, these can be written in the form [2]:(2.1)∂u∂x+∂v∂y=0,(2.2)u∂u∂x+v∂u∂y=K∂2u∂y2,(2.3)u∂T∂x+v∂T∂y=αf∂2T∂y2,
where the x and y axes are taken parallel and perpendicular to the surface, u and v are the parallel and normal velocity components to the surface, respectively, and αf denotes the thermal diffusivity, K=μ/ρ is the kinematic viscosity and density ratio of the ambient fluid.Marangoni effect is incorporated as a boundary condition relating the temperature field to the velocity. The boundary conditions at the surface (aty=0) are(2.4)μ∂u∂y|y=0=-σT∂T∂x|y=0,v(x,0)=0,T(x,0)=T(0,0)+Axm+1,u(x,∞)=0,∂T∂y|y=∞=0,
as y→∞, where σT=dσ/dT, A denotes the temperature gradient coefficient, m is a parameter relating to the power law exponent. The case m=0 refers to a linear profile, m=1 to the quadratic one. The minimum value of m is −1 which corresponds to no temperature variation on the surface and no Marangoni induced flow.Introducing the stream functionψ by(2.5)u=∂ψ∂y,v=-∂ψ∂x,
equation (2.2) is reduced to(2.6)∂ψ∂y∂2ψ∂y∂x-∂ψ∂x∂2ψ∂y2=K∂3ψ∂y3.Applying similarity functions(2.7)ψ=C1xaf(η),Θ=T-T(0,0)Axm+1,
and similarity variable η=C2xby with C1=(m+1)μσTA/ρ23,C2=(m+1)ρσTA/μ23,a=(m+2)/3 and b=(m-1)/3 one can obtain from the system of partial differential equations (2.5), (2.6) one single ordinary differential equation of the third order(2.8)f′′′-2m+13f′2+m+23ff′′=0,
and boundary conditions (2.4) become(2.9)f(0)=0,f′′(0)=-1,f′(∞)=0.
For (2.3) by the similarity temperature function Θ with the corresponding boundary conditions we get(2.10)(m+1)f′Θ-m+23fΘ′=1PrΘ′′,(2.11)Θ(0)=1,Θ′(∞)=0,
where Pr=μ/(ραf) is the Prandtl number. For the dimensionless stream function f(η) and the temperature field Θ(η), the system (2.8), (2.10) is derived and the primes denote the differentiation with respect to η.Now, the velocity components can be expressed by similarity functionf as follows(2.12)u(x,y)=∂ψ∂y=κ2ρμ3x(2m+1)/3f′(η),v(x,y)=-∂ψ∂x=-κ3x(m-1)/3[(m+2)μρ3f(η)+(m-1)κρμ3x(m-1)/3yf′(η)],
where κ=(m+1)σTA/ρ3, σT= constant.It should be noted thatu and v are proportional to x(2m+1)/3 and x(m-1)/3, respectively. It means, that for m=-1/2 the velocity component u is a constant on the upper surface of the boundary layer. If m=1 then η=(m+1)ρσTA/μ23y. In the case of m>1,v is proportional to x(m-1)/3 and is strictly monotone increasing to infinity as x tends to infinity, which is not accepted in physics. Therefore, we restrict our investigations for the interval -1<m≤1.We note that the special casem=1 do admit explicit solution. In [21, 22] the solution to (2.8), (2.9) is given by(2.13)f(η)=1-e-η,
and easy computation shows that(2.14)Θ(η)=Φ(Pr)-Ψ(Pr)e-η+Ω(Pr)e-2η,
with Φ(Pr)=[1/[(Pr-1)/Pr+Pr/(Pr-2)-2]][(Pr-1)/Pr], Ψ(Pr)=2/[(Pr-1)/Pr+Pr/(Pr-2)-2], Ω(Pr)=[1/[(Pr-1)/Pr+Pr/(Pr-2)-2]][Pr/(Pr-2)] is the solution to (2.10), (2.11).Due to the inherent complexity of such flows, to give exact analytical solutions of Marangoni flows are almost impossible. Exact analytical solutions were given by Magyari and Chamkha for thermosolutal Marangoni convection when the wall temperature and concentration variations are quadratic functions of the location [6].Our goal is to present approximate exponential series solution to the nonlinear boundary value problem (2.8), (2.9), moreover to (2.10), (2.11). Several values of the power law exponent and Prandtl number are considered. The influences of the effects of these parameters are illustrated.


## 3. Exponential Series Solution

First, our aim is to determine the approximate local solution off(η) to (2.8), (2.9). We replace the condition at infinity by one at η=0. Therefore, (2.8), (2.9) is converted into an initial value problem of (2.8) with initial conditions(3.1)f(0)=0,f′(0)=ζ,f′′(0)=-1.In view of the third of the boundary conditions (2.9), let us take the solution of the initial value problem (2.8), (3.1) in the form(3.2)f(η)=α(A0+∑i=1∞Aidie-αηi),
where α>0, A0=3/(m+2), Ai(i=1,2,…) are coefficients and α>0 and d are constants. Conditions in (2.9) yield the following equations:(3.3)α(A0+∑i=1∞Aidi)=0,α3∑i=1∞i2Aidi=-1.
It may be remarked that the classic Briot-Bouquet theorem [23] guarantees the existence of formal solutions (3.2) to (2.8), (3.1), the value of A0, and also the convergence of formal solutions.Let us introduce the new variableZ such as(3.4)Z=de-αη.It is evident that the third boundary condition in (2.9) is automatically satisfied. From differential equation (2.8) with (3.2) we get(3.5)-∑i=1∞i3AiZi+m+23(A0+∑i=1∞AiZi)∑i=1∞i2AiZi-2m+13(∑i=1∞iAiZi)2=0.Equating the coefficients of like powers ofZ one can obtain the expressions for coefficients A2,A3,… with m and A1:(3.6)A2=-112A12(m-1),A3=1216A13(m-1)(m-2),A4=-115552A14(m-1)(4m2-15m+17),A5=14665600A15(m-1)(62m3-371m2+757m-610),A6=-146656000A16(m-1)(32m4-257m3+810m2-1171m+730),A7=1740710656000A17(m-1),×(25742m5-263609m4+1108202m3-2419211m2+2737856m-1383380),⋮
From system (3.3) with the choice of A1=1 the parameter values of d and α can be numerically determined. By these parameters the complete series solution (3.2) is reached. The calcutaled values of d, α and f′(0) are shown in Table 1, and the variation of f′(0) with m on Figure 1.Table 1mdαζ=f′(0)−0.7−3.6470382351.1515955552.124598444−0.6−2.9657609801.1274158341.983315576−0.5−2.6377576811.063879191.732325541−0.4−2.3761728621.0333540731.593916052−0.3−2.1623107101.0144144561.494034266−0.2−1.9840743281.0018200701.415321059−0.1−1.8331837710.99339785011.3506758060−1.7037580500.98793949661.2961852350.1−1.5914983540.98467330131.2493678420.2−1.4931868630.98307107321.2085321220.3−1.4063657450.98275608581.1724721170.4−1.3291245510.98345178041.1402996280.5−1.2599554230.98495053901.1113434380.6−1.1976520640.98709361081.0850853410.7−1.1412377580.98975771031.0611178970.8−1.0899131100.99284580341.0391156680.9−1.0430174650.99628062091.0188150711−111Figure 1
Variation ofζ with m.The series forms forf(η) and f′(η) are given below for some special values of the exponent m(m=-0.5;m=0;m=1):(3.7)m=-0.5:f(η)=2.127758380-2.806255505(e-1.06387919η)+0.9252777520(e-1.06387919η)2-0.3389803468(e-1.06387919η)3+0.1266709685(e-1.06387919η)4-0.04756400678(e-1.06387919η)5+0.01788183745(e-1.06387919η)6-0.006724299058(e-1.06387919η)7+0.002528494356(e-1.06387919η)8-0.0009506669013(e-1.06387919η)9+0.0003573924970(e-1.06387919η)10,f′(η)=2.985516834(e-1.06387919η)-1.968767491(e-1.06387919η)2+1.081902410(e-1.06387919η)-0.5390504296(e-1.06387919η)4+0.2530117851(e-1.06387919η)5-0.1141446885(e-1.06387919η)6+0.05007689285(e-1.06387919η)7-0.02152010022(e-1.06387919η)8+0.009102552598(e-1.06387919η)9-0.003802224402(e-1.06387919η)10,m=0:f(η)=1.481909245-1.683209870(e-0.9879394966η)+0.2389818639(e-0.9879394966η)2-0.04524080827(e-0.9879394966η)3+0.009099650359(e-0.9879394966η)4-0.001854352461(e-0.9879394966η)5+0.0003780882937(e-0.9879394966η)6-0.00007689156390(e-0.9879394966η)7+0.00001559224546(e-0.9879394966η)8-0.000003154096351(e-0.9879394966η)9+0.0000006368276490(e-0.9879394966η)10,f′(η)=1.662909512(e-0.9879394966η)-0.4721992446(e-0.9879394966η)2+0.1340855441(e-0.9879394966η)3-0.03595961598(e-0.9879394966η)4+0.009159940185(e-0.9879394966η)5-0.002241170151(e-0.9879394966η)6+0.0005317494905(e-0.9879394966η)7-0.0001232335611(e-0.9879394966η)8+0.00002804450725(e-0.9879394966η)9-0.000006291471870(e-0.9879394966η)10,m=1:f(η)=1-e-η,f′(η)=e-η.
It can be seen that for the case m=1 the obtained solution coincides with the exact solution (2.13). The effect of the exponent m on the velocity profiles f′(η) is illustrated in Figure 2. The values of f′(0)=ζ decrease as m is changing from negative values to positive ones.Figure 2
Variation off′ with η.Applying the series solution forf the second-order linear differential equation (2.10) for Θ can be solved similarly, which presents the temperature distribution. Here we define Θ(η) as the series(3.8)Θ(η)=B0+∑i=1∞Bidie-αηi,
with coefficients Bi(i=0,1,2,…) and hence the individual coefficients will be determined from differential equation (2.10) with (3.2) as follows(3.9)B1=A1B0PrPr-1(m+1),B2=112A12B0Pr(Pr-1)(Pr-2)(3m2Pr+m2+6mPr+3Pr-1),B3=-1216A13B0Pr(Pr-1)(Pr-2)(Pr-3)F(Pr,m),F(Pr,m)=((m3-m)(3Pr2-19Pr-2)+(m2-1)(4Pr2-20Pr+4)),⋮Remark that these coefficients as expressions ofB0 can be calculated only for noninteger values of the low Prandtl numbers. In (2.11) the second boundary condition is automatically satisfied, and from the first condition coefficient B0 is to be determined, that is, from the equation(3.10)B0+B1d+B2d2+B3d3+⋯=1,
together with (3.9) (see Table 2).Table 2B0Pr∖m−0.5010.270.6993601030.6177162890.5563436130.70.2618673400.1874641880.1444444442.5−6.8688903240.7323522340.1666666675.5−0.010101191−0.0203388152.10000004970−129.591844090.92275412391/6298−521.5866457−253.81722067326/25ForΘ(η) with Prandtl number Pr =298 and three values of m(m=-0.5;m=0;m=1) the first ten terms are given below(3.11)m=-0.5:Θ(η)=-521.5866457+690.2257842(e-1.06387919η)-228.3499003(e-1.06387919η)2+83.59931856(e-1.06387919η)3-31.32427417(e-1.06387919η)4+11.72818950(e-1.06387919η)5-4.437848713(e-1.06387919η)6+1.649735001(e-1.06387919η)7-.6357069408(e-1.06387919η)8+.2268403841(e-1.06387919η)9-0.09549145798(e-1.06387919η)10,m=0:Θ(η)=-253.8172206+433.8991704(e-0.9879394966η)-185.8554249(e-0.9879394966η)2+23.32308704(e-0.9879394966η)3-11.64354979(e-0.9879394966η)4-.6456997006(e-0.9879394966η)5-1.774772066(e-0.9879394966η)6-.8708402162(e-0.9879394966η)7-.6992264953(e-0.9879394966η)8-.5109496856(e-0.9879394966η)9-.4045737228(e-0.9879394966η)10,m=1:Θ(η)=732625-4410475(e-η)+2220175(e-η)2.It may be noted that the Prandtl numberPr=298 corresponds to the power transformer oil. We point out that for the case m=1 the solution Θ(η) coincides with the exact solution (2.14).The effects of the power law exponentm and the Prandtl number are exhibited in Figures 3–14, where Pr =0.27 corresponds to the mercury and Pr=0.7 corresponds to the air. Figures 3–5 illustrate the influence of the Prandtl number on the temperature Θ for m=1. It can be observed in Figure 3 that for low Prandtl numbers 0.27≤Pr≤1.00001 the maximum value of Θ decreases as Pr increases and for a high Prandtl numbers 2.5≤Pr≤7.00001 and 70≤Pr≤298 the maximum value of Θ increases as Pr increases. In all three cases the boundary layer thickness increases as Pr increases. Figures 6–8 depict the effect power exponent m for fixed values of Pr. It can be observed in Figures 6 and 7 that the boundary layer thickness increases as m increases and the maximum value of Θ decreases as m increases for Pr=0.27,2.2, while for high Prandtl number (Pr=298) the reverse effect of m on the maximum of Θ can be seen. Figures 9–11 illustrate the effect of Pr on Θ′ for m=1 and Figures 12–14 represent the effect of m for Θ′ for Pr=0.27,2.2,298.Figure 3
Variation ofΘ with Pr (0.27≤Pr≤1.00001) for m=1.Figure 4
Variation ofΘ with Pr (2.5≤Pr≤7.00001) for m=1.Figure 5
Variation ofΘ with Pr (70≤Pr≤298) for m=1.Figure 6
Variation ofΘ for Pr=0.27.Figure 7
Variation ofΘ for Pr=2.2.Figure 8
Variation ofΘ for Pr=298.Figure 9
The effect of the Prandtl number onΘ′ for m=1 (0.27≤Pr≤1.00001).Figure 10
The effect of the Prandtl number onΘ′ for m=1 (2.5≤Pr≤7.00001).Figure 11
The effect of the Prandtl number onΘ′ for m=1 (70≤Pr≤298).Figure 12
The effect ofm on Θ′ for Pr=0.27 (m=-0.5;0;1).Figure 13
The effect ofm on Θ′ for Pr=2.2 (m=-0.5;0;1).Figure 14
The effect ofm on Θ′ for Pr=298 (m=-0.5;0;1).


## 4. Conclusion

In this paper the incompressible flow and heat transfer over a flat impermeable plate has been investigated. The resulting governing equations have been transformed into a system of nonlinear ordinary differential equations by applying suitable similarity transformation. For these equations approximate exponential series solutions are determined and the effects of the power exponent and the Prandtl number are illustrated in Figures1–14, and it is observed that the values of f′ decrease as power exponent m increases (see Figures 1 and 2), moreover the boundary layer thickness increases as m or Pr increases. The temperature profiles are exhibited in Figures 3–8, it is observed that for low Prandtl number the temperature Θ decreases as Pr increases and for high Prandtl numbers the influence of Pr is opposite.



---

*Source: 314989-2012-04-10.xml*