Table Of Content

Linear and Nonlinear Surface Waves in Electrohydrodynamics 5 M. Hunt, J-M. Vanden-Broeck, D. Papageorgiou, E. Parau 1 0 2 January 13, 2015 n a J Abstract 2 Theproblemofinterestinthisarticlearewavesonalayeroffinitedepth 1 governed by the Euler equations in the presence of gravity, surface tension, and vertical electric fields. Perturbation theory is used to identify ] canonical scalings and to derive a Kadomtsev-Petviashvili equation with n y an additional non-local term arising in interfacial electrohydrodynamics. d Whenthe Bond numberis equalto 1/3, dispersion disappears and shock - waves could potentially form. In the additional limit of vanishing elec- u tricfields,anewevolutionequationisobtainedwhichcontainsthird-and l f fifth-orderdispersion as well as a non-local electric field term. . s c i 1 Introduction s y h ClassicalwaterwavemodelsleadingtoequationssuchastheKadomtsev-Petviashvili p (KP) equation have long since been used to understand the nonlinear phenom- [ ena and interactions of waves and are useful in the theoretical basis for further 1 studies. Themaininteresthereistoextendtheresultsin[8]tothethreedimen- v sionalcase. Theimportanceofinterfacialelectrohydrodynamicsphenomenahas 3 been highlighted in many different cases. 8 7 2 2 Set-Up and Governing Equations 0 . 1 Consider a perfectly conducting, inviscid, irrotational and incompressible fluid 0 (region 1) bounded below by a wall electrode at y = h and bounded above 5 − 1 by a free surface y = η(t,x,y), here h is the mean depth of the surface. The : fluid motionis describedby a velocitypotential ϕ(t,x,y,z)satisfying Laplace’s v equation in region 1. Surface tension with coefficient σ and gravity, g, are i X included. The region y > η(t,x,y), denoted by region 2, is occupied by a r hydrodynamically passive dielectric having permittivity ǫp. It is assumed that a there are no free charges or currents in region 2 and therefore the electric field can be represented as a gradient of a potential function, E= V. A vertical −∇ electrical field is imposed by requiring that V E y as y , where E 0 0 ∼ − → ∞ is constant. The voltage potential satisfies the Laplace equation. On the free surface the Bernoulli equation holds: 2 2 2 ∂ϕ 1 ∂ϕ ∂ϕ ∂ϕ p + + + +gη+ =C (1) ∂t 2 ∂x ∂y ∂z ρ "(cid:18) (cid:19) (cid:18) (cid:19) (cid:18) (cid:19) # 1 The pressure p is obtained through the Young-Laplace equation: [nˆ T nˆ]1 =σ nˆ (2) · · 2 ∇· Where the stress tensor is given by: 1 T= pδ +ǫ E E δ E E (3) ij p i j ij k k − − 2 (cid:18) (cid:19) The unit normal is given by: (∂ η,∂ η, 1) nˆ = x y − (4) (1+(∂xη)2+(∂yη)2)21 The governing equations are then: ∂2ϕ ∂2ϕ ∂2ϕ + + =0 on h6z 6η(t,x,y) (5) ∂x2 ∂y2 ∂z2 − ∂2V ∂2V ∂2V + + =0 on z >η(t,x,y) (6) ∂x2 ∂y2 ∂z2 ∂η ∂ϕ∂η ∂ϕ∂η ∂ϕ + + = on z =η(t,x,y) (7) ∂t ∂x∂x ∂y ∂y ∂z ∂ϕ 1 ∂ϕ 2 ∂ϕ 2 ∂ϕ 2 P 1 σ + + + +gη+ nˆ Σ nˆ = nˆ+C (8) ∂t 2"(cid:18)∂x(cid:19) (cid:18)∂y(cid:19) (cid:18)∂z(cid:19) # ρ −ρ · · ρ∇· ∂V ∂η∂V + =0 on z =η(t,x,y) (9) ∂x ∂x ∂z ∂ϕ =0 on z = h (10) ∂z − V E z as z (11) 0 ∼− →∞ 3 Linear Theory Wemoveintoareferenceframewhichismovingwiththefluid. Thiswillremove the time derivatives from the governing equations. 3.1 Finite Depth The scaling we use for the finite depth case is: x=hxˆ, y =hyˆ, z =hzˆ, η =hηˆ ρh3 σh σ (12) t= tˆ, ϕ= ϕˆ, V =E hVˆ, P = pˆ 0 r σ s ρ h 2 Therearetwodimensionlessparameterswhichcomeoutofthenon-dimensionlisation process: ρh2g ǫE2ρh B = , E = 0 (13) b σ σ The parameters are called the Bond number and electric Bond number respec- tively. The variables are expanded as: p = εp +o(ε) (14) 1 η = εη +o(ε) (15) 1 ϕ = Ux+εϕ +o(ε) (16) 1 V = z+εV +o(ε) (17) 1 − The variables V ,η ,ϕ ,p are written as inverse Fourier transforms: 1 1 1 1 1 V = eiklx+ilyA(k,l,z)dkdl (18) 1 (2π)2 R2 Z 1 ϕ = eiklx+ilyB(k,l,z)dkdl (19) 1 (2π)2 R2 Z 1 η = eiklx+ilyC(k,l)dkdl (20) 1 (2π)2 R2 Z The expression for the free surface is then given by: 1 µei(kx+ly)pˆtanhµ η(x,y)= dkdl, (21) 4π2 R2 k2U2 (Bµ Ebµ2+µ3)tanhµ Z − − where µ=√k2+l2. The stability relation is given as: tanhµ U2 =(Bµ E µ2+µ3) (22) − b k2 Typicalfree surfaceprofiles predictedby (21) for E =1.5, B =2 and U =0.5 b in figure (3.1). 3.2 Infinite Depth The difference between the finite depth and the infinite depth is the boundary conditions, these new boundary conditions are given by: ∂ϕ U as z (23) ∂x → →∞ ∂ϕ 0 as z (24) ∂z → →−∞ Using the expansion: ϕ = Ux+εϕ +o(ε) (25) 1 η = εη +o(ε) (26) 1 V = E z+εV +o(ε) (27) 0 1 − p = εp +o(ε) (28) 1 3 Figure 1: Linear Waves Profiles - Finite Depth The same approach can be done for the infinite case, by writing the variables as inverse Fourier transforms, it is possible to write the free surface as: 1 µpˆρ−1 η = ei(kx+ly)dkdl (29) 1 4π2 R2 k2U2 gµ+ ǫE02µ2 σµ3 Z − ρ − ρ Where pîs the Fourier transform of the pressure. Equation (29) can be made clearer by the introduction of the following variables: α=U2g−1k, β =U2g−1l, xˆ=xgU−2, yˆ=ygU−2, η =gηÛ−2 (30) 1 To obtain: 1 νe−ν202ei(αxˆ+βyˆ) ηˆ (xˆ,yˆ)= dαdβ (31) 1 8π2 ZR2 α2−ν+µ1ν2−µ2ν3 where: ǫE2 σg µ = 0 µ = (32) 1 ρU2 2 ρU4 Figure (2) show the waves for infinite depth for values µ = 1 and µ = 2. 1 2 When β =0, the denominator in the integrandcan be arrangedto be a perfect square provided that µ and µ satisfy a certain relation 1 2 1+µ =2√µ . (33) 2 2 4 Weakly Nonlinear Theory Previous work in this area has been carried out by Katsis and Alylas ([11]), wherebyatwofluidscenariowithinterfacewasgivenbyz =η(t,x,y). Analysis 4 Figure 2: Profile of 2D Infinite Depth Wave showedthattheresultingequationwasthesameasthatobtainedinthissection, despite therebeingadifferentmethodanddifferentcoefficients. Forthe weakly nonlinear theory, use the following scaling ([1],[8],[7]): λ gλa x=λxˆ, y =µyˆ, t= tˆ, η =aηˆ ϕ= ϕˆ (34) c c 0 0 V =λE Vˆ z(1) =hzˆ z(2) =λZˆ, (35) 0 Define three parameters ([8],[1]): a h2 λ2 α= , β = , γ = (36) h λ2 µ2 The governing equations become: 1 ∂2ϕˆ ∂2ϕˆ ∂2ϕˆ +γ + = 0 16zˆ6αηˆ (37) β ∂z2 ∂yˆ2 ∂xˆ2 − ∂2Vˆ ∂2Vˆ ∂2Vˆ +γ + = 0 αηˆ6Zˆ 6 (38) ∂zˆ2 ∂yˆ2 ∂xˆ2 ∞ 2 2 2 E ∂Vˆ ∂Vˆ ∂Vˆ b Σ = γ (39) 33 2  ∂zˆ! − ∂x! − ∂y!    ∂Vˆ ∂ηˆ∂Vˆ +α β =0 on z =αηˆ. (40) ∂xˆ ∂xˆ ∂yˆ p 5 1∂ϕˆ ∂ηˆ ∂ϕˆ∂ηˆ ∂ϕˆ∂ηˆ = +α +αγ on z =αηˆ (41) β ∂zˆ ∂tˆ ∂xˆ ∂xˆ ∂yˆ∂yˆ E ∂ϕ E βBˆ nˆ b = +η+p+ bnˆ σ nˆ+ ∇· − 2α ∂t α · · 1 ∂ϕ 2 ∂ϕ 2 ∂ϕ 2 + α +αγ + (42) 2 ∂xˆ ∂y ∂z " (cid:18) (cid:19) (cid:18) (cid:19) (cid:18) (cid:19) # Bβ ∂2ηˆ(1+α2βγ(∂ ηˆ)2)+γ∂2ηˆ(1+α2β(∂ ηˆ)2) 2α2βγ∂ ηˆ∂ ηˆ∂ ∂ ηˆ xˆ yˆ yˆ xˆ − xˆ yˆ xˆ yˆ h (1+α2β(∂xˆηˆ)2+α2βγ(∂yˆηˆ)2)32 i (43) Where: σ ǫ E2 B = , E = p 0 (44) ρgh2 b ρgh The non-dimensional parameter B is the Bond number. Note that E can be b thought of as the ratio of two things, ǫ E2/ρ and gh. The quantity gh has the p 0 units of speed2 and therefore so does ǫ E2/ρ. The speed ǫ E2/ρ is therefore p 0 p 0 characteristic to the system of interest, the quantity √E can be thought of as b p the electric Foude number. The constant was calculated by using the solution ϕˆ = ηˆ= 0 and Vˆ = zˆto show that C = E /2α. To begin with restrict the b − − attention to the classical shallow water scaling by taking: α=β =γ =ε 1, T =εtˆ, X =xˆ tˆ (45) ≪ − From here on in, drop the hats. The asymptotic expansions for the variables are: ϕ = ϕ +εϕ +ε2ϕ +o(ε2) (46) 0 1 2 η = η +εη +o(ε) (47) 0 1 3 3 V = z+ε2V1+o ε3 (48) − p = εp (cid:16) (cid:17) (49) 1 Theexpansionfor(48)comesfromexamining(40). Equation(37)canbesolved as: (z+1)2∂2ϕ (z+1)4∂4ϕ (z+1)2∂2ϕ ϕ=ϕ ε 0 +ε2 0 0 +o(ε2) (50) 0− 2! ∂X2 4! ∂X4 − 2! ∂y2 (cid:20) (cid:21) Where the ε−1 equation shows that ϕ = ϕ (T,X). The electric term in the 0 0 Bernoulli equation (42) is: Eb 1 3∂V1 ε2 (51) ε −2 − ∂Z (cid:20) (cid:21) The first part of this equation cancels the Bernoulli constant and then this 3 leaves a term of order ε2, so to include this into the order ε equation, the 6 electric Bond number is scaled according to E = Eˆ √ε. This makes the O(ε) b b Bernoulli equation become: 1 1 η =B∂2 ∂ ϕ + ∂3ϕ p +Eˆ ∂ V (∂ ϕ )2 (52) 1 X − T 0 2 X 0− 1 b z − 2 x 0 Movingontothefreesurfaceequation(41),theO(1)equationgivesη =∂ ϕ , 0 X 0 which was known previously and the O(ε) equation is: 1 1 ∂4ϕ ∂2ϕ η ∂2ϕ =∂ η ∂ η +∂ ϕ ∂ η (53) 6 X 0− y 0− 2 0 X 0 T 0− X 1 X 0 X 0 Inserting the expression for η coming from (52) shows that: 1 1 1 1 3 ∂4 ϕ ∂2ϕ =2∂ ∂ ϕ + B ∂4ϕ +∂2p + ∂2(η2)+Eˆ ∂2 ∂ V =0 6 X 0− y 0 T X 0 2 3 − X 0 X 1 2 X 0 b X Z (cid:18) (cid:19) (54) In terms of η the equation is: 0 ∂ ∂η 1 1 ∂3η ∂p 3 ∂η Eˆ ∂2V 1∂2η 0 0 1 0 b 0 + B + + η + + =0 ∂X "∂T 2(cid:18)3 − (cid:19) ∂X3 ∂X 2 0∂X 2 ∂X∂Z# 2 ∂y2 (55) The O(1) equation for (38) needs to be solved, it is given by: ∂2V ∂2V 1 1 + =0 (56) ∂X2 ∂Z2 This equation requires a boundary condition in order to write down a solution. All boundary conditions are zero except the one at Z = 0. The boundary condition can be found by expanding (40) to yield: ∂V ∂η 1 =0 (57) ∂X − ∂X Integrating this equation and assuming that the free surface dies off at infinity shows that V = η on Z = 0. Equation (56) is solved by use of a Green’s 1 function. The Green’s function for a 2D Laplace equation in the upper half plane is: 1 (X X′)2+(Z Z′)2 ′ ′ g(X,Z X ,Z )= log − − (58) | 2π (X X′)2+(Z+Z′)2 (cid:20) − (cid:21) and use of Green’s second identity: g 2V V 2gdτ = g(nˆ )V V (nˆ )gdΣ (59) 1 1 1 1 ∇ − ∇ ·∇ − ·∇ ZD Z∂D Inserting the boundary condition into (59) shows that the solution is given by: 1 Zη 0 ′ V = dX (60) 1 −π R (X X′)2+Z2 Z − Then the required result is: ∂V ∂η 1 =H 0 (61) ∂Z (cid:12) ∂X (cid:12)Z=0 (cid:18) (cid:19) (cid:12) (cid:12) (cid:12) 7 So inserting (61) into (55) shows that: ∂ ∂η 1 1 ∂3η ∂p 3 ∂η Eˆ ∂2η 1∂2η 0 + B 0 + 1 + η 0 + bH 0 + 0 =0 ∂X "∂T 2(cid:18)3 − (cid:19) ∂X3 ∂X 2 0∂X 2 (cid:18)∂X2(cid:19)# 2 ∂y2 (62) In terms of dimensional variables the equation is: ∂ 1 ∂η 3 ∂η h2 1 ∂3η 1 ∂p E ∂2η + η + B + + bH + ∂x c ∂x 2h ∂x 2 3 − ∂x3 ρg∂x 2 ∂x2 (cid:20) 0 (cid:18) (cid:19) (cid:18) (cid:19)(cid:21) 1∂2η + =0 2∂y2 (63) Whenp=E =0,(63)reducestothe standardKPequation. Inthe2Dcase,it b has been shownthat if B <1/3,then “generalised”solitary wavesare possible. The equation which we have derived is exactly the same as equation (2.7) in [13] but with Eˆ replaced with the ratio of densities of the fluids. So there b appears to be a link with interfacial flows, at least on the level of equations of motion. There is different behaviour around B = 1/3 and the next section details the derivation of a new equation in this case. Note that a complete analytical derivation of all terms in the equation has been given whereas they simply add the nonlinear term in [13]. The dispersion relation can be easily derived for the problem, it is given by: ǫE2 σ ω2 =µ g 0µ+ µ2 tanhhµ (64) − ρ ρ (cid:18) (cid:19) where µ = √k2+l2 and k and l are the wavenumbers. Equation (64) can be expanded in powers of k and l to obtain the linear portion of (63). 5 Analysis around B = 1/3 In this section, the scaling in (45) is replaced by: α=ε2, β =ε, γ =ε2, T =ε2tˆ, X =xˆ tˆ (65) − The expansions are now given by: ϕ = ϕ +εϕ +ε2ϕ +ε3+o(ε3) (66) 0 1 2 η = η +εη +ε2η +o(ε2) (67) 0 1 2 5 5 V = Z+ε2V1+o ε2 (68) − p = ε2p +o(ε2) (cid:16) (cid:17) (69) 1 1 B = +εB +o(ε) (70) 1 3 The governing equation for ϕ is solved in exactly the same fashion as in the previous case, the solution is given by: (z+1)2∂2ϕ (z+1)4∂4ϕ ϕ=ϕ (T,X) ε 0 +ε2 0 (71) 0 − 2! ∂X2 4! ∂X4 8 Only the derivative of the third term needs to be calculated, this is: ∂ϕ 1 ∂6ϕ ∂2ϕ 3 0 0 = (72) ∂z −5! ∂X6 − ∂y2 The O(1) Bernoulli equation (42) once again shows that ∂ ϕ = η , the O(ε) X 0 0 part of (42) is: 1∂2η ∂ϕ 0 1 = +η (73) 3∂X2 −∂X 1 The term ϕ has been calculated in (71), the O(ε2) part of the free surface 1 equation (41) is: ∂ η ∂ η +∂ ϕ ∂ η = η ∂2ϕ +∂ ϕ (74) T 0− X 2 X 0 X 0 − 0 X 0 z 3 To find η , the O(ε2) Bernoulli equation (42) is used. The same rescaling as 2 before is used with the electric Bond number, Eb as Eb = Eˆbε23 in order to include the electric term in. The equation becomes: 1 1 B ∂2 η + ∂2η =∂ ϕ ∂ ϕ +η +p + (∂ ϕ )2 Eˆ ∂ V (75) 1 X 0 3 X 1 T 0− X 2 2 1 2 X 0 − b Z 1 Inserting this into (74) shows that: 1 2∂ η B ∂3η +3η ∂ η + ∂5η + T 0− 1 X 0 0 X 0 45 X 0 (76) +∂2ϕ +∂ p +Eˆ ∂ ∂ V =0 y 0 X 1 b X z 1 Inthe previoussection,the term∂ ∂ V wasexaminedbefore andthe resultis X y 1 just quoted, ∂ ∂ V =H(∂2 η ), making the final equation: X y 1 X 0 ∂ ∂η 1 ∂5η 3 ∂η B ∂3η 0 0 0 1 0 + + η + ∂X ∂T 90∂X5 2 0∂X − 2 ∂X3 (cid:20) (77) Eˆ ∂2η 1∂p 1∂2η + bH 0 + 1 + 0 =0 2 ∂X2 2∂X 2 ∂y2 # (cid:18) (cid:19) In dimensional variables the equation becomes: ∂ 1 ∂η ∂η h4 ∂5η 3 ∂η h2 1 ∂3η + + + η B + ∂x c ∂t ∂x 90∂x5 2h ∂x − 2 − 3 ∂x3 (cid:20) 0 (cid:18) (cid:19) (78) E h ∂2η 1 ∂p 1∂2η + b H + + =0 2 ∂x2 2ρg∂x 2∂y2 (cid:18) (cid:19) (cid:21) 6 Fully Nonlinear Results in Infinite Depth Only the infinite depth case in considered for simplicity, but an extension to include the finite depth is possible. The conditions at infinity are V E z, for z . 0 ∼− →∞ ∂φ 1, for z . ∂x → →−∞ 9 TheinterestedisinsteadywavestravellingwithaconstantspeedU,sochoosea framemovingwith the waveandnon-dimensionalisethe equationsbyusing the σ unit velocity U and the unit length L = . Also V is non-dimensionalised ρU2 using E L. The non-dimensional parameters are 0 ǫE2 σg Eˆ = 0 α= . b ρU2 ρU4 These are the exact same parameters obtained in the linear case with µ = Eˆ 1 b and µ = α. This is very useful as it can be used to determine the rage of 2 validity of the linear solution. Set V =0 at the free surface z =η(x,y) and introduce Vˆ =V +V z, which 0 satisfies Vˆ 0 as z . The Bernoulli’s equation becomes in this steady → → ∞ frame 2 1 ∂φ 2 ∂φ 2 ∂φ 2 1 Eˆ ∂Vˆ 1 Eˆ b b + + +αη K + + =0, 2"(cid:18)∂x(cid:19) (cid:18)∂y(cid:19) (cid:18)∂z(cid:19) #−2 − − 2 ∂nˆ 1+ηx2+ηy2 2  q (79) wherenˆ isthedownwardsunitnormal(seeequation(4)),andK isthecurvature given by η η x y K = + .  1+η2+η2  1+η2+η2 x y x y x y q  q  By applying the second Greens identity for φ(x,y,ζ(x,y) x in the region z <η(x,y)andforVˆ inthe regionz >η(x,y), weobtainthe b−oundaryintegral equations ∂ 1 1 ∂ 2π(φ(Q) x )= (φ(P) x ) dS (φ(P) x )dS , − Q − P ∂nˆ·∇ R P− R ∂nˆ − P P Z ZSF (cid:18) PQ(cid:19) Z ZSF PQ (80) ∂ 1 1 ∂Vˆ 2πη(x ,y )= η(x ,y ) dS + dS , Q Q − P P ∂nˆ ·∇ R P R ∂n P Z ZSF (cid:18) PQ(cid:19) Z ZSF PQ (81) whereP andQarepointsonthefree-surfaceS withcoordinates(x ,y ,ζ(x ,y )) F P P P P and (x ,y ,ζ(x ,y )) and Q Q Q Q R = PQ = (x x )2+(y y )2+(ζ(x ,y ) ζ(x ,y ))2. PQ P Q P Q P P Q Q | | − − − q The equations (79)-(81) are desingularized and projected on the Oxy plane, andthenintegratednumerically(seeLandweber&Macagno1969,Forbes1989, Pa˘r˘au & Vanden-Broeck 2002 for details). For most of the computations 40 grid points are chosen in x direction and 40 points in y direction, with a grid interval of dx = 0.8 and dy = 0.8. The fully-localised solitary waves near the minimum of the dispersion relation have been computed which corresponds to the region 2 1+Eˆ b α . ≥ 2 ! 10