One‐dimensional models for slender axisymmetric viscous liquid bridges

One-dimensional models for slender axisymmetric viscous liquid bridges F. Javier Garcı ´a Departamento de Fı ´sica Aplicada, E.U.I.T.A., Universidad de Sevilla, Ctra. Utrera km 1, 41013 Sevilla, Spain Antonio Castellanosa) Departamento de Electro ´nica y Electromagnetismo, Facultad de Fı ´sica, Universidad de Sevilla, Avda. Reina Mercedes s/n, 41012 Sevilla, Spain A set of one-dimensional models, previously derived for liquid jets, is generalized to viscous liquid bridges by applying suitably modified boundary conditions at the anchoring disks. A linear analysis for small-amplitude perturbations around the cylindrical static solution is performed. The oscillation frequencies and growth factors so obtained are compared to the already known linear three-dimensional results for a wide range in both the slenderness and viscosity. The relative error of each model is studied in terms of the typical axial length. Good agreement is found for slender enough bridges. The existence of boundary layers for weakly dissipative liquid bridges in the context of one-dimensional models is also discussed. I. INTRODUCTION The stability and dynamics of liquid bridges, i.e. a body of fluid anchored between two solid supports and held by capillary forces, has received much attention in the literature. The stability problem, whose solution is independent of the viscosity of the liquid, may be solved exactly using a static analysis. The solution for the linear dynamics of cylindrical inviscid liquid bridges is also well known.1However, its solution for viscous liquid bridges on the basis of general three-dimensional ~hereafter 3-D!hydrodynamic equations involves great difficulties, and only partial results are known. Tsamopoulos et al.2have formally solved the linear problem for arbitrarily viscous liquid bridges. However, a tedious computation is necessary, and convergence is only guaranteed for Ohnesorge numbers Cunder 0.5. A boundary-layer treatment of this linear problem, valid for very small values of C, has been carried out by Borkar and Tsamopoulos,3and improved by Higuera et al.4On the other limit, Nicola ´s5has studied this linear eigenproblem for paramount viscosity. The finite-amplitude oscillations of viscous axisymmetric liquid bridges have been studied by Chen and Tsamopoulos.6The nonlinear dynamic equations that describe this 3-D problem cannot be solved analytically, and hard computation is necessary in order to obtain numerical solutions. Because of the above-mentioned difficulties, several authors have used one-dimensional ~hereafter 1-D!models, initially derived for liquid jets, to study the dynamics of axisymmetric viscous liquid bridges. Most of these models7 retain integrally the capillary pressure term. Recently, Bechtel et al.8have performed a self-consistent asymptotic 1-D analysis for slender inviscid and viscous liquid jets. Their perturbative approach develops all terms in the Navier–Stokes equations, including the capillary pressure one, in terms of a nondimensional wave number. The resulting leading-order models do not account for the finite instability cutoff of the exact linear 3-D analysis of Rayleigh.9 Finite, though approximate values for this cutoff are recovered by higher-order models. This lack of accuracy is due to the truncation of the capillary pressure. In the usual 1-D models, the approach is not fully consistent from a formal mathematical point of view, in the sense that all orders in the small parameter are retained in the capillary term. However, on the one hand, it is desirable to conserve the term which leads the stability and dynamics of the physical system under consideration, and this cannot produce worse results than the Bechtel et al. models. As a matter of fact, the former models reproduce more accurately the well known 3-D linear results. On the other hand, the complexity of the models so derived hardly increases, and its main virtue remains unchanged: they are one-dimensional and do not depend on the radial variable. Meseguer10 extended the inviscid slice model of Lee11 and the Cosserat model12 to liquid bridges. The inviscid slice model neglects radial momentum and viscous effects. Due to its simplicity, it has been extensively used since ~see Perales and Meseguer13 for references!. The correct viscous contribution has recently been derived by Eggers and Dupont14 and, independently, by Garcı ´a and Castellanos.7In order to deal with viscous liquid bridges, some authors 10,13,15–18 have used the more elaborate Cosserat model, which enhances the Lee model in the inviscid case, since it accounts partially for radial momentum effects. Unfortunately, its range of application reduces to small viscosities, due to an inconsistency in the viscous terms.13,7 Garcı ´a and Castellanos7have recently derived a set of 1-D models for slender liquid jets, which include the Lee, Eggers, and Cosserat models. For slender bridges, the above-mentioned new set of 1-D models, suitably modified to take account of the boundary conditions at the anchors, is expected to be useful to find the oscillation frequencies and damping rates for arbitrary values of the viscosity. Here the results obtained with these models are compared to the linear 3-D solutions given by Sanz,1 Higuera et al.,4Tsamopoulos et al.,2and Nicola ´s,5for invis- a!Author for correspondence. Fax: 34-5-4239434, Electronic-mail: [email protected] cid, weakly viscous, moderately viscous, and highly viscous liquid bridges, respectively. In this paper, a wide range in both the slenderness (1<L<10) and the Ohnesorge number (0<C,`) is studied, where the slenderness is the nondimensional length of the liquid bridge and the Ohnesorge number is the nondimensional viscosity ~to be defined precisely below!. Notice that liquid bridges of slenderness larger than p are considered. The analysis of such unstable columns is interesting for two reasons. First, it is possible to have liquid bridges of slenderness above p , using either electric19,20 or magnetic21 forces to stabilize the interface, or by melting a solid metal rod by means of an intense electric current passing through it.22 Second, it is interesting to predict how many main drops, aside from the satellite drops, are going to be produced in the breaking of such slender columns, as well as to determine the minimum critical length above which these liquid bridges behave as jets. II. EQUATIONS AND BOUNDARY CONDITIONS Let us consider an axisymmetric liquid bridge ~see Fig. 1!held between two parallel rigid circular disks of radius R, separated by a distance L, and coaxially aligned. Here the influence of the surrounding gas as well as gravity are neglected. Also density r , viscosity m , and surface tension s are supposed to be uniform and constant, and only axisymmetric configurations are considered. The equations for the 1-D models for liquid jets are derived in detail in Ref. 7. With the latter assumptions, the same dimensionless equations apply to liquid bridges, but now additional boundary conditions at the disks must be added. Let Vand Wbe the radial and axial components of the velocity, which depend on r,z, and t, while the shape of the interface is given by r5F(z,t). In nondimensional coordinates, the disks are placed at z56L, where L5L/(2R)is the slenderness. Impenetrability of each of the rigid disks and anchoring to their respective edges yield W~r,6L,t!50, ~1! F~6L,t!51. ~2! As long as CÞ0, where C5 m /( rs R)1/2 is the Ohnesorge number, the no-slip condition at the disks must be also satisfied: V~r,6L,t!50. ~3! Finally, the volume of the liquid bridge must be conserved, which gives E 2L LdzF2~z,t!5 E 2L LdzF2~z,0!.~4! The volume of the column is important for its stability properties. However, if the dynamics of the system is under study, the conservation of volume is guaranteed by fulfillment of the continuity equation and the kinematic condition. Therefore, Eq. ~4!can be considered as an initial condition. A. Mean-velocity models The Lee, Cosserat, and averaged models may be called mean-velocity models, since all of them have the mean axial velocity on a slice W ¯ and the shape of the interface Fas dependent variables. Their differential order in zis four for any value of the Ohnesorge number.7Thus, four independent boundary conditions at the disks must be satisfied. First, the anchoring conditions ~2!remain unchanged: F~6L,t!51. ~5! Second, as Meseguer10 has shown, the radial integration of Eq. ~1!on a slice implies that the mean axial velocity must be zero at the disks: W ¯ ~6L,t!50. ~6! Finally, introducing the latter four conditions in the averaged kinematic equation, which relates Fto W ¯ , another two conditions for W ¯ are obtained: W ¯ z~6L,t!50, ~7! which are dependent with respect to ~5!. It is important to note that the no-slip condition ~3!leads to ~7!again, as can be shown by averaging the continuity equation on a slice at z56L. Therefore the boundary conditions at the disks are always the same, independent of the value of the Ohnesorge number. Finally, notice that the above boundary conditions have been derived without the aid of any approximation. Thus, the relative error of the mean-velocity models, defined as the order of magnitude of the neglected terms divided by the retained ones, is expected to be the same as for liquid jets. FIG. 1. Schematic description of an axisymmetric liquid bridge. B. Parabolic model The dependent variables of the parabolic model are W0,W2, and F. In terms of these, the velocity field is given by W~r,z,t!5W0~z,t!11 2r2W2~z,t!,~8! V~r,z,t!52 1 2rW0z21 8r3W2z.~9! When the viscosity is not zero, the differential order in zis eight, while it changes to six in the inviscid case.7This imposes the number of boundary conditions in each case. The anchoring conditions at the edges of the rigid disks ~5!apply again. If ~8!and ~9!are introduced in the impenetrability and no-slip conditions, ~1!and ~3!, the following eight boundary conditions for the viscous parabolic model are obtained: W0~6L,t!50, W2~6L,t!50, ~10! W0z~6L,t!50, W2z~6L,t!50. ~11! The four conditions ~11!do not apply if C50. Instead, substituting ~8!and ~9!into the kinematic equation gives another two conditions. In this particular case, the irrotationality of the velocity field allows decoupling W2from the problem. Therefore the following six conditions for W0must be fulfilled by the inviscid parabolic model: W0~6L,t!50, W0zz~6L,t!50, ~12! W0z~6L,t!21 8W0zzz~6L,t!50. ~13! Conditions ~13!come from imposing the anchoring of the shape of the interface to the edges of the disks. However they allow the liquid to slip on their surface, as is to be expected as long as viscosity effects are absent. Contrarily to the mean-velocity models, the parabolic one conserves the sensitivity of the general equations and boundary conditions to the viscous or inviscid character of the problem. III. LINEAR ANALYSIS Let us consider small perturbations around the static cylindrical solution. The linear 1-D models previously obtained for liquid jets7are valid for liquid bridges, as long as the above discussed boundary conditions are added. To solve the linear dynamics associated with the 1-D models, an exponential dependence in time is tried for the variables appearing in each model: f~z,t!5Re@f ˆ~z!eVt#,~14! w ¯ ~z,t!5Re@w ¯ ˆ ~z!eVt#,~15! w0~z,t!5Re@w ˆ0~z!eVt#,~16! w2~z,t!5Re@w ˆ2~z!eVt#,~17! where Vis complex in general, and its real and imaginary parts are the growth factor a and the oscillation frequency v , respectively; and the lower case letters f,w ¯ ,w0, and w2denote the respective perturbations of F,W ¯ ,W0, and W2.Although the evolution of a liquid bridge is an initialvalue problem with respect to time, we do not consider initial conditions here. Instead, we address our interest to a modal analysis, from which a countable infinite set of eigenvalues Vmas well as their corresponding eigenfunctions are obtained. The final solution of a particular problem, characterized by its initial conditions, would be an appropriate superposition of these eigenmodes. The above linear problem has well-defined parity with respect to z, since the equations and boundary conditions are invariant under the change of sign of z. Therefore, the eigenmodes can be classified in either antisymmetric or symmetric, according to the parity of the shape of the interface with respect to z. Note that w ¯ ,w0, and w2have opposite parity to the shape of the interface, since applying a derivative with respect to zto any quantity changes its parity. Similar to the dependence on time, an exponential dependence on zcould be tried. However, the well-defined parity of the problem suggests proposing a combination of hyperbolic sines and cosines, which is equivalent. A. Mean-velocity models The substitution of the proposed dependence on both t and zin the equations of the mean-velocity models leads to the following general solution for w ¯ ˆ : w ¯ ˆ 5( j51 2 @Ajcosh~ k jz!1Bjsinh~ k jz!#,~18! where 6 k j~with j51,2) are the four complex roots of the biquadratic equation k 41c1 k 21c050. ~19! The coefficients c0and c1, which are functions of V, take different values for each model: for Lee, c052V2,c15126CV;~20! for Cosserat, c052V2 11~C/4!V,c15126CV21 4V2 11~C/4!V;~21! for averaged, c052V2,c15126CV21 4V2.~22! Notice that the antisymmetric and symmetric modes with respect to the interface are characterized by Bj50 and Aj50, respectively. Substitution of ~18!in the linear counterparts of ~6!and ~7!gives two systems of homogeneous algebraic equations for the antisymmetric and symmetric modes. In order to have nontrivial solutions, the determinant formed by the coefficients of such algebraic equations must be zero. For antisymmetric modes, this yields k 2sinh~ k 2L!cosh~ k 1L!2 k 1sinh~ k 1L!cosh~ k 2L! 50, ~23! while the analogous condition for symmetric modes is k 2cosh~ k 2L!sinh~ k 1L!2 k 1cosh~ k 1L!sinh~ k 2L! 50. ~24! Given the parameters Land C, this compatibility condition is satisfied by a countable infinite set of eigenvalues Vm, where mis odd for antisymmetric modes and even for symmetric ones. The numerical procedure implemented to find such complex roots is described in the Appendix. B. Parabolic model A similar treatment can be applied to the parabolic model. The dependence on zproposed for w0induces an analogous one for w2when substituted in the equations, although with different integration constants. The resulting system of two homogeneous equations lead to the general solution w ˆ05( j51 J @Ajcosh~ k jz!1Bjsinh~ k jz!#, ~25! w ˆ25( j51 J @Aj h jcosh~ k jz!1Bj h jsinh~ k jz!#, where h jdepends on Cand V, and is given by h j52 k j 2V2C~81 k j 2! 2V1C~1626 k j 2!,~26! and 6 k j~with j51,2,...,J) are the 2Jroots of c4 k 81c3 k 61c2 k 41c1 k 21c050, ~27! whose coefficients are functions of Cand V: c05128CV2116V3, c1564C1~82384C2!V296CV224V3, c2548C1~7196C2!V118CV2,~28! c35215C2~1114C2!V, c45C. Again, the antisymmetric and symmetric modes are characterized by Bj50 and Aj50, respectively. Notice that Jis the degree of the polynomial equation ~27!. Therefore J53ifC50, and J54 otherwise. As shown before, the boundary conditions are also different, depending on the Ohnesorge number. Consequently, the viscous and inviscid cases must be treated separately. If the general solution ~25!is introduced in the linear counterparts of either ~10!and ~11!if CÞ0, or ~12!and ~13! if C50, a system of homogeneous algebraic equations for the integration constants Ajand Bjis obtained. However, a suitable combination of these equations leads to two systems of equations, one for antisymmetric modes and the other for symmetric modes. Again, the determinants formed with the coefficients of these homogeneous systems of algebraic equations must be zero. In the viscous case, the corresponding conditions to find the antisymmetric and symmetric eigenvalues are, respectively, ( i,j,k,l51 4 Eijkl cosh~ k iL! k jsinh~ k jL! h k 3cosh~ k kL! h l k lsinh~ k lL!50~CÞ0!,~29! ( i,j,k,l51 4 Eijkl sinh~ k iL! k jcosh~ k jL! h ksinh~ k kL! 3 h l k lcosh~ k lL!50~CÞ0!,~30! where Edenotes the completely antisymmetric tensor. In the inviscid case, the analogous conditions for antisymmetric and symmetric modes are ( j,k,l51 3 Ejkl cosh~ k jL! k kcosh~ k kL! k l~82 k l 2! 3sinh~ k lL!50~C50!,~31! ( j,k,l51 3 Ejkl sinh~ k jL! k ksinh~ k kL! k l~82 k l 2! 3cosh~ k lL!50~C50!.~32! IV. RESULTS AND DISCUSSION As will be shown below, the relative errors committed in the eigenvalues provided by the linear analysis of the 1-D models are reasonably small. In comparison with the 3-D computations, which are very heavy and present convergence problems,2the 1-D ones are very fast and efficient. This allows a careful parametric study of the eigenvalues and eigenfunctions, which provides valuable information about the influence of the slenderness and the viscosity on the evolution of the bridge. Figure 2 shows the set of three eigenvalues associated with the index m51, for a typical viscous liquid bridge. For each m, there is a value Lm *of the slenderness under which two eigenvalues are complex conjugate, with a <0. For L.Lm *, the latter ones become real. One of them, labeled FIG. 2. Frequency v ~–––!and growth factor a ~—!versus the slenderness, as given by the parabolic model for m51 and C50.2. Vdm , increases monotonically as Ldoes. The critical slenderness Lcm for which Vdm becomes positive, and therefore unstable, determines the stability of the mth mode. The other eigenvalue, Vsm , decreases monotonically. The modes corresponding to these will be called dominant and subdominant, respectively, since the former dominates the evolution of the liquid bridge, with the latter dying down relatively quickly. Thus, while nothing is said Vmdenotes Vdm . The two modes discussed above, called nearly inviscid when C!1,4are the only ones appearing if C50. Owing to the greater differential order of the viscous equations, new modes appear as long as CÞ0. These viscous modes have real negative eigenvalues Vvmthat are approximately proportional to Cfor small Ohnesorge numbers, and are subdominant for realistic values of C~see Fig. 2!. The range of validity of the referred models is limited by the condition of large slenderness. Concerning the other limit, it could be thought unrealistic considering values greater than p , for which the liquid bridge is unstable. However, as we have commented before, such slender liquid columns can be obtained by applying suitable electric or magnetic fields, or by melting a solid metal rod. In particular, using electric fields,20 dielectric liquid bridges with L.5 have been obtained, and greater values of the slenderness are attainable with slight changes in the geometry of the system. If the electric field which stabilizes the dielectric liquid bridge is suddenly removed at t50, the early evolution of the bridge is governed by the linear equations considered here. For t.0 the equations are purely hydrodynamic, because the depolarization of the liquid is instantaneous compared to any mechanical time of interest. Therefore, we will deal with slenderness within the range 1<L<10. Although less realistic, much more slender liquid bridges are considered to study the asymptotic connection between liquid bridges and jets. No restrictions have been imposed on the value of viscosity. The inviscid and very viscous cases deserve special attention, since the above-cited exact linear 3-D results provide a good test of 1-D cases in both limits. As already pointed out, the limit of small viscosity is also considered in detail, in order to compare the behavior of the 1-D models with respect to the boundary layers. Finally, we also discuss the moderately viscous case. A. Inviscid case In Fig. 3, the module of the most significant value of length l,defined as the distance between consecutive nodes of the shape of the interface. lis a half of the wavelength in a jet, but does not necessarily coincide with Lfor slender enough bridges. Let lmax be the value of lfor which the maximum growth factor is attained in a jet, which is about 4.51 in the inviscid case. In a bridge and for a given value of L, the most unstable mode should have a similar value of l. However, this value can never be greater than L, which shows that the slenderness can be considered a good measure of the characteristic axial length if L,lmax ~for m51!. For more slender bridges, it is expected that l.lmax , independent of the value of L, and new nodes appear in the shape of the interface. In fact, this is what can be observed in Fig. 4, where the shape of the interface is plotted against z/Lfor increasing values of L,C50, and m51. For a given mode m,lmis approximately given by lm. H 2L m11 S if L,m11 2lmax D , lmax S if L.m11 2lmax D . ~33! In fact the shapes of the interface of the odd and even modes are well approximated by FIG. 3. Module of the eigenvalue Vof the first four modes as a function of the slenderness, for C50. The short dashed line represents the growth factor of the inviscid infinite jet. The long dashed horizontal line is the tangent at the maximum of the latter curve. FIG. 4. Shape of the interface of the first ~antisymmetric!mode versus z/L, for C50 and L52, 3, 5, 10, 20; as given by the 3-D analysis of Sanz ~Ref. 1!. Vm ~dominant modes! is plotted against L, for the first four modes and C50, as calculated by Sanz.1 The short dashed line corresponds to the growth factor of the infinite jet, taking k5 p /L as the dimensionless wavenumber, while the long dashed line is the tangent line at the maximum of the former. Notice that, for L@Lcm and C50, the growth factor does not tend to zero, as happens to liquid jets. Instead it approaches a constant value which is common to all the modes, and coincides with the value of the maximum growth factor of an infinite jet. This is also true for the 1-D models for any finite value of the Ohnesorge number, and it can be explained in terms of a characteristic nondimensional axial f2n21~z!5 5 cos S n p z 2L D sin~ p z/lmax!~n51,3, . . . !, sin S n p z 2L D cos~ p z/lmax!~n52,4, . . . !; ~34! f2n~z!5 5 cos S n p z 2L D cos~ p z/lmax!~n51,3, . . . !, sin S n p z 2L D sin~ p z/lmax!~n52,4, . . . !; ~35! where the subscript denotes the index of the corresponding mode. The typical axial length is responsible for the relative error of the 1-D models. On the one hand, we have seen that lis directly related to Lif the latter is not larger than lmax . Thus, the mentioned error is expected to increase as Ldecreases or mincreases, due to the decreasing distance between nodes. On the other hand, ltends to lmax if the bridge is slender enough. This implies that the error must tend to a constant value as Lincreases. Indeed, the above reasoning is confirmed by Fig. 5, where the data of Fig. 3 have been used to calculate the relative errors in the eigenvalue as given by the 1-D models, for C50 and m51 and 2. The error of the Cosserat and averaged models, which coincide in the inviscid limit, is clearly of the same order than the one of the parabolic model. In contrast, the error of the Lee model is about twenty times greater than the previous ones. This confirms the theoretical expectations about the relative order of each of these models. Concerning the mean axial velocity and the shape of the interface, we have found that their relative error, measured as the deviation from the exact results divided by their amplitude, is very similar to that shown above for the eigenvalues. The same is not true for the velocity components themselves, as we show below. The different 1-D approximations of the velocity field can be compared to the linear 3-D solution given by Sanz1in the linear inviscid case. Radial profiles of the axial and radial velocities are shown in Figs. 6~a!and 6~b!, respectively, for C50 and L53. Two approximations are plotted. On the one hand the one-term averaged one, given by Meseguer,10 approximates Was its mean value on each slice, and Vas a linear function of r. Although the averaged model has been chosen, the results are the same within the plot accuracy for any mean-velocity model, the relative error being l22.On the other hand, the parabolic model provides a two-term polynomial approximation, whose relative error is expected to be l24. The one-term averaged approximation behaves worse in the center of the bridge (z.0) and near the disks (z.1). In particular, the liquid slipping on the disks is not predicted. However, it can be considered to be near the best one-term polynomial fit to the exact velocity field. In general, the parabolic approximation shows to be better than the previous one. It almost coincides with the exact results, except very near the disks, where the profile is harder to be fitted by a two-term polynomial in r. However it succeeds in describing the cited slipping on the disks, and even the slightly negative axial velocity near the anchors, already noted by Sanz.1 FIG. 5. Relative error in the determination of Vby the Lee ~-••-), Cosserat and averaged ~–––!, and parabolic ~---!models as a function of the slenderness, for the first two modes (m51,2) and C50. The numbers on the curves indicate the mode index. FIG. 6. Radial profile of the axial ~a!and radial ~b!velocities for the first mode, C50, and L53. The numbers indicate the values of z/L. Results given by Sanz ~Ref. 1!~—!, parabolic ~---!, and one-term averaged ~–––!approximations. With respect to the role of the slenderness, the error increases as Ldecreases. The radial profiles of the axial velocity are approximately constant on a slice for slender bridges @e.g. see Fig. 6~a!#. This circumstance relaxes for smaller L, and the error of the 1-D approximations increases accordingly. Taking the amplitude of w ¯ as reference, the relative error in the velocity field for C50 and L53is about 2% for the parabolic, and 6% for the one-term averaged approximations. For L52, the analogous errors are 5% and 13%, respectively. B. Paramount viscosity We now focus our attention on the case when viscous Acomparison of the different approximations of the velocity field against a 3-D solution in this limit would be interesting. Unfortunately, the latter result is not available. The main features predicted by the parabolic approximation are the no slipping on the disks, and that the maximum axial velocity happens at the center of the bridge, and not at the interface as for C50. This is in qualitative agreement with the 3-D results for the stream function,5and cannot be predicted by the one-term averaged approximation, due to its limited dependence on r. C. Weakly viscous case Borkar and Tsamopoulos,3and more recently Higuera et al.,4have studied the linear boundary-layer problem associated with weakly viscous liquid bridges. They have shown that there are two Stokes boundary layers, one at each disk, as well as another one at the interface. These regions have a very small thickness, of order (C/V0)1/2, where V0is the eigenvalue for C50. The most important correction to the inviscid values of both the damping rate and the frequency is associated with viscous dissipation at the disks ~of order C1/2), and at the bulk ~of order C). The influence of the interfacial boundary layer is significantly smaller, of order C3/2. The damping rate for small viscosities is approximately given by the real part of V1C1/21V2C. The first coefficient happens to be much smaller than the second, which justifies the second-order analysis of Higuera et al.4Consequently, the damping rate goes as C1/2 for very small Ohnesorge numbers, while it goes as Cfor more realistic values of C. Due to the higher order of the viscous momentum equations, new modes appear. These viscous modes, found by Higuera et al.,4are purely damped ( v 50 and a ,0), and for small viscosities their growth factor is approximately proportional to C, and only for unrealistically small Ohnesorge numbers (C;1027or less!do these modes become important. For increasing, not necessarily small, values of C, the damping rate of such modes also increases monotonically. Therefore the evolution of the bridge is mainly described by the most significant nearly inviscid modes. The oscillation frequency and the shape of the interface are almost the same as in the inviscid case. Therefore, they are well estimated by the 1-D models. Now we examine whether these models can predict the correct values of the damping rate, the appearance of viscous modes, and the large velocity gradients present at the boundary layers. The answer is negative for the mean-velocity models and partially positive for the parabolic one, as we show below. The damping rate, which is not zero as long as C.0, has been provided by Higuera et al.4In Fig. 8, the damping rate 2 a of the first mode is plotted against L, for C50.002 and different approximations. The single crosses are data obtained numerically by Tsamopoulos et al.2Notice that, for this value of C, the result given by the parabolic model agrees very well with the data from the second-order boundary-layer analysis. The error of the mean-velocity models is clearly greater, the averaged one giving slightly worse results. This is to be expected, since the damping rate for very small Cis related to the dissipation in the Stokes FIG. 7. Relative error in the determination of C a by the Lee and averaged ~–––!, Cosserat ~-•-), and parabolic ~---!models as a function of the slenderness, for the first two modes (m51,2) and paramount viscosity. The numbers on the curves indicate the mode index. effects are paramount. This holds for C@ a l2, which occurs not only for very large Ohnesorge numbers, but also near the stability limit.5 V is real for any positive value of L, and its dependence on the Ohnesorge number tends asymptotically to a constant divided by C. Therefore, the quantity C a is considered. The commentaries previously made about the typical axial length are valid for any value of C, provided that the wavelength of maximum growth factor for the infinite jet lmax is computed for the given value of viscosity. This length tends to infinity as C does,7 which means that l;L even for very slender bridges. The error must therefore tend to zero as L tends to infinity. This is found in Fig. 7, where the relative error of C a for paramount viscosity, computed through the data given by Nicola ´s,5 is plotted against L, for m51 and 2. Such deviations are related to the viscous terms neglected in the derivation of the 1-D models. In this limit, the Cosserat model exhibits a worse agreement with the exact solution, due to the inconsistency in the second-order viscous terms. As also occurs for jets, the better results of the averaged model get worse for L,2, since the approximation of some viscous terms in its derivation is not very good for small slenderness. Concerning the parabolic model, the error in the first two modes does not overcome a 10% for L>1. As expected, the error increases with the mode index, due to the decreasing typical axial length. boundary layers, which is not described by the latter. In general, this error diminishes as Lincreases, but is much larger than the error in the oscillation frequency. The dependence of the damping rate on Cis shown in Fig. 9, for L52 and 3. Again the parabolic model is more accurate than the mean-velocity models. The dependence on Cof the damping rate given by the mean-velocity models is very approximately proportional to C, which gives worse results for very small C. This seems to indicate that these models account for the dissipation in the bulk, and not for the boundary layers at the disks. As a consequence, the relative error is significantly greater in the damping rate than in the oscillation frequency. The relative error of the parabolic model also increases for very small C, but it is not proportional to C, and the error is significantly smaller. Note also that the cited 3-D approximation differs from the 1-D ones for large enough viscosities, for which the boundary-layer analysis is not as accurate. As already pointed out, this is more evident for L53, where viscous effects are greater for the same Ohnesorge number, and it becomes apparent in Fig. 8 as well. The differential order in zof the Lee, Cosserat, and averaged models is four, independent of the value of C.In addition, the boundary conditions at the disks are the same for the inviscid and viscous cases. Therefore, the mentioned mean-velocity models cannot account for the boundary layers, and no viscous modes are found. This is to be expected, since the no-slip condition on the disks does not affect either the axial velocity or its mean value. Besides, the large gradient of the axial velocity near the interface is confined to a very thin layer, whose average influence on a slice is very small. Thus, W ¯ and its derivatives change slowly as Cincreases from zero. The same is true for the shape of the interface, related to the former through the kinematic condition. Since the one-term approximation of the velocity field is in terms of W ¯ and F, it can only describe it in the bulk, and not in the boundary layer. Instead, the differential order in zof the parabolic model changes from six when C50 to eight when C.0. Besides, the number of boundary conditions at the disks changes analogously. Therefore this model accounts for the appearance of viscous modes. Although the relative error is large ~up to 50% for the first two viscous modes!, the importance of these modes is secondary, since they die down quite rapidly. When Cis small enough, a large gradient in the radial velocity is to be expected near the wall. This is in fact found in Fig. 10, where an axial profile of the radial velocity at r50.5 is plotted for several values of the Ohnesorge number, as given by the parabolic model. The 3-D inviscid solution v0, given by Sanz,1is also shown. According to the leading term in the 3-D boundary-layer analysis given by Higuera et al.,4a uniform approximation to the radial velocity vcan be readily found from the exact inviscid solution v0, namely FIG. 8. Damping rate 2 a versus L, for the first ~antisymmetric!mode and C50.002. As obtained by the Lee ~-••-!, averaged ~---!, Cosserat ~-•-!, and parabolic ~---!models; and by the boundary-layer analysis of Higuera et al. ~Ref. 4!~—!. The crosses are data from Tsamopoulos et al. ~Ref. 2!. FIG. 9. Damping rate 2 a versus C, as obtained by the Lee ~-••-!, averaged ~---!, Cosserat ~-•-!, and parabolic ~---!models; and by the boundary-layer analysis of Higuera et al. ~Ref. 4!~—!. For the first ~antisymmetric!mode and L52,3. FIG. 10. Axial profile of the radial velocity near an anchoring disk, as given by the parabolic model ~---!and the inviscid solution of Sanz ~Ref. 1!as well as boundary-layer estimations obtained from the latter ~—!. For the first mode, L53, and r50.5. The numbers on the curves indicate the values of C. v~r,z,t!5v0~r,z,t! $ 12exp@2~11i!~L2z!/ d # % , ~36! where d 5(2C/ v 0)1/2 is a measure of the boundary-layer thickness, and v 0is the oscillation frequency for the inviscid case. This uniform approximation is also plotted in Fig. 10 for the same values of C. The gradients of the velocity field are well predicted by the parabolic model, the difference in the values of the velocity itself coming from the error of the inviscid solution. For Cgreater than 1022the physical boundary layer cannot be considered narrow, and the inviscid solution does not fit the viscous one anymore. This illustrates the fact that the viscous effects may become important even for quite small values of C, provided Vis small enough. This occurs near the stability limit for each mode, as in the preceding example, in which L53. D. Moderately viscous case The scarce data available from a 3-D approach prevent us from doing a careful comparison. Nevertheless, some aspects related to the typical axial length, the eigenvalues, and the approximations to the velocity field are discussed. As commented before, relation ~33!holds for any given value of C, provided lmax is computed for that Ohnesorge number. Since lmax increases as Cdoes,7the value of Lfor which new nodes appear in the shape of the interface increases. This fact is put in evidence in Fig. 11, where the shape of the interface fis plotted for the first mode, L510 and several values of C, as given by the parabolic model, the other models giving similar results. Notice that the distance between nodes increases with C, and it becomes Lfor C510. For the latter value of C, the viscous effects can be considered paramount. Meseguer10 has carried out a study of the dependence of Von Lfor several values of C, by means of the Cosserat model. Quite similar results are obtained with the other models derived above. A comparison to the exact linear results of Tsamopoulos et al.2is hardly interesting, since most of these are for quite thick bridges. However, relative errors can be calculated for L510 p /9 and 2 p /3, with C50.1. For the first data (L.3.50), the errors of the Lee, Cosserat, averaged and parabolic models are respectively 0.6%, 0.3%, 0.09%, and 0.095%; and for the second value (L.2.09) we have 3.1%, 0.5%, 0.6%, and 0.6%. The former data show the inconsistency in the viscous terms of the Cosserat model, even for such a small value of C, and are in qualitative accordance with the corresponding results for paramount viscosities. Note that in this case the weight of the viscous terms is large, in spite of the relatively small Ohnesorge number, because u V u is small enough. For the first mode, this happens near its stability limit, at L5 p . V. CONCLUSIONS A set of 1-D models, previously derived for jets, has been formulated in the field of axisymmetric viscous liquid bridges, under the hypothesis of large slenderness. The solution of the linear eigenproblem associated with each model has been compared to the available known results. It has been shown that the relative error of the cited models in determining the frequency v and growth factor a is related to the magnitude of the typical axial length l. For not very slender columns, l.L/m, and the index of the mode m gives the number of nodes of the shape of the interface between the disks. As Lincreases new nodes appear at the eigenfunctions, lremains approximately constant, and the growth factor of the corresponding eigenmode tends to a max . Since lmax increases as Cdoes, this explains why the cited relative errors tend to be constant for increasing values of L—this constant decreasing to zero as Cincreases from zero to infinity. The different 1-D approximations to the velocity field have been discussed and compared in the inviscid and viscous cases. In the central region of the bridge, far from the disks, the results are similar to those for jets. The well known one-term approximation, given by Meseguer, gives an error of order l22. An inviscid and a very viscous two-term averaged approximations, that improve the former for weakly and highly viscous liquid jets, respectively, are more accurate in the same limits of Cin the central part of the bridge. However, it is worse near the disks, since they fulfill the boundary conditions for the velocity field at the disks only approximately. In contrast, the accuracy of the velocity field given by the parabolic model is satisfactory, even near the disks. It accounts for the nonzero radial velocity at the disks when C50 and also predicts the large velocity gradient near the disks when Cis very small. Further work is necessary to determine the value of these models to study the breakup of liquid bridges. A common objection is that the axial length scales near the rupture might not be large, which would affect the error of these approximations.10,18 In this sense, the work of Eggers and Dupont14 with the viscous Lee model is very encouraging, since they predict a good agreement with a 3-D approach in the pinch region. To conclude this study, we would like to stress the small effort necessary to find each eigenvalue and its corresponding eigenfunctions. The 1-D approach is at least four orders of magnitude more economical in time of computation than the 3-D one in the viscous case.2Contrarily to the latter, the former has no problems of convergence, for any value of C or L, and all the modes can be easily distinguished. Besides, the method of solution is formally the same for any value of the Ohnesorge number. FIG. 11. Shape of the interface of the first ~antisymmetric!mode versus z/L, for C50, 0.1, 1, 10, and L510, as given by the parabolic model.