Pierre Gaillard
email: pierre.gaillard@u-bourgogne.fr

Deformations of the seventh order Peregrine breather solutions of the NLS equation with twelve parameters

We study the solutions of the one dimensional focusing NLS equation. Here we construct new deformations of the Peregrine breather of order 7 with 12 real parameters. We obtain new families of quasirational solutions of the NLS equation. With this method, we construct new patterns of different types of rogue waves. We recover triangular configurations as well as rings isolated. As already seen in the previous studies, one sees appearing for certain values of the parameters, new configurations of concentric rings.

Introduction

The nonlinear Schrödinger equation was first solved by Zakharov and Shabat [1] in 1972 by the inverse scattering method. The first expressions of the quasi-rational solutions were given by Peregrine [2] in 1983. From this time, a considerable number of studies were carried out. Eleonski, Akhmediev and Kulagin obtained the first higher order analogue of the Peregrine breather [3] in 1986. Akhmediev et al. [4,5], constructed other analogues of order 3 and 4, using Darboux transformations. Rational solutions of the NLS equation have been written in 2010, as a quotient of two Wronskians in [6]. An other representation of the solutions of the NLS equation has been constructed in [7] in 2011, also in terms of a ratio of two Wronskians determinants of order 2N . In 2012, Guo, Ling and Liu constructed an other representation of the solutions of the focusing NLS equation, as a ratio of two determinants has been given in [8] using generalized Darboux transform. In the same year, Ohta and Yang [9] have given a new approach where solutions of the focusing NLS equation by means of a determinant representation, obtained from Hirota bilinear method. A the beginning of the year 2012, one obtained a representation in terms of determinants which does not involve limits [10]. The two formulations given in [7,10] did depend in fact only on two parameters; this remark was first made by V.B. Matveev. Then we found for the order N (for determinants of order 2N), solutions depending on 2N -2 real parameters. In this article, we restrict ourself the study to the case of the solutions of NLS of order 7; because of the constraints of the publication, we do not have the space to publish all the deformations. With this new method, we construct news deformations at order 7 with 12 real parameters. The explicit representation in terms of polynomials is found, but is too monstrous to be published. One constructs various drawings to illustrate the evolution of the solutions according to the parameters. One obtains at the same time triangular configurations and ring structures with a maximum of 28 peaks. These deformations are completely new and gives by new patterns a better understanding of the NLS equation.

Determinant representation of solutions of NLS equation

We recall the results obtained in [7] and [10]. We consider the focusing NLS equation

iv t + v xx + 2|v| 2 v = 0. (1)
In the following, we consider 2N parameters λ ν , ν = 1, . . . , 2N satisfying the relations

0 < λ j < 1, λ N +j = -λ j , 1 ≤ j ≤ N. (2)
We define the terms κ ν , δ ν , γ ν by the following equations,

κ ν = 2 1 -λ 2 ν , δ ν = κ ν λ ν , γ ν = 1 -λ ν 1 + λ ν , (3)
and

κ N +j = κ j , δ N +j = -δ j , γ N +j = 1/γ j , j = 1 . . . N. (4)
The terms x r,ν (r = 3, 1) are defined by

x r,ν = (r -1) ln γ ν -i γ ν + i , 1 ≤ j ≤ 2N. ( 5
)
The parameters e ν are defined by

e j = ia j -b j , e N +j = ia j + b j , 1 ≤ j ≤ N, (6)
where a j and b j , for 1 ≤ j ≤ N are arbitrary real numbers. We use the following notations :

A ν = κ ν x/2 + iδ ν t -ix 3,ν /2 -ie ν /2, B ν = κ ν x/2 + iδ ν t -ix 1,ν /2 -ie ν /2, (7)
for 1 ≤ ν ≤ 2N , with κ ν , δ ν , x r,ν defined in (3), ( 4) and ( 5).

The parameters e ν are defined by (6).

Here, the parameters a j and b j , for 1 ≤ N are chosen in the form

a j = N -1 k=1 ãk ǫ 2k+1 j 2k+1 , b j = N -1 k=1 bk ǫ 2k+1 j 2k+1 , 1 ≤ j ≤ N. (8)
We consider the following functions :

f 4j+1,k = γ 4j-1 k sin A k , f 4j+2,k = γ 4j k cos A k , f 4j+3,k = -γ 4j+1 k sin A k , f 4j+4,k = -γ 4j+2 k cos A k , (9)
for 1 ≤ k ≤ N , and

f 4j+1,N +k = γ 2N -4j-2 k cos A N +k , f 4j+2,N +k = -γ 2N -4j-3 k sin A N +k , f 4j+3,N +k = -γ 2N -4j-4 k cos A N +k , f 4j+4,N +k = γ 2N -4j-5 k sin A N +k , (10)
for 1 ≤ k ≤ N . We define the functions g j,k for 1 ≤ j ≤ 2N , 1 ≤ k ≤ 2N in the same way, we replace only the term A k by B k .

g 4j+1,k = γ 4j-1 k sin B k , g 4j+2,k = γ 4j k cos B k , g 4j+3,k = -γ 4j+1 k sin B k , g 4j+4,k = -γ 4j+2 k cos B k , (11)
for 1 ≤ k ≤ N , and

g 4j+1,N +k = γ 2N -4j-2 k cos B N +k , g 4j+2,N +k = -γ 2N -4j-3 k sin B N +k , g 4j+3,N +k = -γ 2N -4j-4 k cos B N +k , g 4j+4,N +k = γ 2N -4j-5 k sin B N +k , (12)
for 1 ≤ k ≤ N .

Then we get the following result :

Theorem 2.1 The function v defined by v(x, t) = det((n jk) j,k∈[1,2N ] ) det((d jk) j,k∈[1,2N ] ) e 2it-iϕ (13)
is a quasi-rational solution of the NLS equation ( 1)

iv t + v xx + 2|v| 2 v = 0, depending on 2N -2 parameters ãj , bj , 1 ≤ j ≤ N -1,
where

n j1 = f j,1 (x, t, 0), n jk = ∂ 2k-2 f j,1 ∂ǫ 2k-2 (x, t, 0), n jN +1 = f j,N +1 (x, t, 0), n jN +k = ∂ 2k-2 f j,N +1 ∂ǫ 2k-2 (x, t, 0), d j1 = g j,1 (x, t, 0), d jk = ∂ 2k-2 g j,1 ∂ǫ 2k-2 (x, t, 0), d jN +1 = g j,N +1 (x, t, 0), d jN +k = ∂ 2k-2 g j,N +1 ∂ǫ 2k-2 (x, t, 0), 2 ≤ k ≤ N, 1 ≤ j ≤ 2N (14)
The functions f and g are defined in ( 9),( 10), (11), (12).

We don't have the space to give the proof in this publication. We will give it in an other forthcoming paper. The solutions of the NLS equation can also be written in the form :

v(x, t) = exp(2it -iϕ) × Q(x, t)
where Q(x, t) is defined by :

Q(x, t) := f 1,1 [0] . . . f 1,1 [N -1] f 1,N +1 [0] . . . f 1,N +1 [N -1] f 2,1 [0] . . . f 2,1 [N -1] f 2,N +1 [0] . . . f 2,N +1 [N -1] . . . . . . . . . . . . . . . . . . f 2N,1 [0] . . . f 2N,1 [N -1] f 2N,N +1 [0] . . . f 2N,N +1 [N -1] g 1,1 [0] . . . g 1,1 [N -1] g 1,N +1 [0] . . . g 1,N +1 [N -1] g 2,1 [0] . . . g 2,1 [N -1] g 2,N +1 [0] . . . g 2,N +1 [N -1] . . . . . . . . . . . . . . . . . . g 2N,1 [0] . . . g 2N,1 [N -1] g 2N,N +1 [0] . . . g 2N,N +1 [N -1] ( 15
)
3 Quasi-rational solutions of order 7 with twelve parameters

Wa have already constructed in [7] solutions for the cases from N = 1 until N = 6, and in [10] with two parameters.

Because of the length of the expression v of the solution of NLS equation with 12 parameters, we can't give here. We only construct figures to show deformations of the analogue of the seventh Peregrine breather; in the following we will call it for simplicity, the seventh Peregrine breather.

Conversely to the study with two parameters given in preceding works [7,10], we get other type of symmetries in the plots in the (x, t) plane. We give some examples of this fact in the following.

It is important to note the similar role played by the parameters ãj and bj for a same j; the same configuration of the peaks is obtained. For this reason one will give the figures only for a parameter ãj or bj . On the other hand, to understand the configuration for a value of the parameter, one will give two sights to see the distribution of the peaks. With different choices of parameters, we obtain all types of configurations : triangles, rings and concentric rings with a maximum of 28 peaks.

Conclusion

We have constructed explicitly solutions of the NLS equation of order N with 2N -2 real parameters. The expressions in terms of polynomials in x and t are too monstrous to be published in this paper.

It is important to note the symmetrical role played by the parameters ãj and bj ; the configurations obtained for one of these two parameters ãj or bj for the index j are the same ones. Thus for each couple (ã j ; bj ) we have only built one associated figure for only for one parameter, ãj = 0, or bj = 0.

In the cases a 1 = 0 or b 1 = 0 we obtain triangles with a maximum of 28 peaks; for a 2 = 0 or b 2 = 0, we have 3 concentric rings with two of them with 10 peaks and an other with 5 peaks with in the center the Peregrine P 2 with 3 peaks. For a 3 = 0 or b 3 = 0, we obtain 4 concentric rings without central peak with 7 peaks on each of them. For a 4 = 0 or b 4 = 0, we have 3 concentric rings with 9 peaks, with a in the center one peak. For a 5 = 0 or b 5 = 0, we obtain 2 concentric rings without central peak with 11 peaks on each of them and the appearance in the center of the Peregrine breather P 3 with 6 peaks. For a 6 = 0 or b 6 = we have only one ring with 13 peaks with inside the appearance of the Peregrine breather of order 5 with 15 peaks. We obtained new patterns in the (x; t) plane, by different choices of these parameters; we recognized rings configurations as already observed in the case of deformations depending on two parameters [7,10]. We get news triangular shapes and multi-concentric rings.

Figure 1 :

 1 Figure 1: Solution of NLS, N=7; all parameters equal to 0, the Peregrine breather of order 7, P 7 .

Figure 2 :

 2 Figure 2: Solution of NLS, N=7; b1 = 10 4 ; we obtain a regular triangle with 28 peaks; on the right, sight of top.

Figure 3 :

 3 Figure 3: Solution of NLS, N=7; ã2 = 10 6 , 3 rings with respectively 5, 10, 10 peaks with in the center the Peregrine of order 2, P 2 ; on the right, sight of top.

Figure 4 :

 4 Figure 4: Solution of NLS, N=7; b3 = 10 10 , 4 rings with 7 peaks on each of them without central peak; on the right, sight of top.

Figure 5 :

 5 Figure 5: Solution of NLS, N=7; ã4 = 10 10 , 3 rings with 9 peaks on each of them with in the center one peak; on the right, sight of top.

Figure 6 :

 6 Figure 6: Solution of NLS, N=7; b5 = 10 15 , 2 rings of 11 peaks with in the center the Peregrine breather of order 3, P 3 ; on the right, sight of top.

Figure 7 :

 7 Figure 7: Solution of NLS, N=7; ã6 = 10 12 , b6 = -10 8 , a ring with 13 peaks and in the center the Peregrine breather of order 5, P 5 ; on the right, sight of top.

This study at the order 7 was never still carried out; it is completely new and makes it possible to provide a better understanding of the phenomena of rogue waves.


Alban Sauret

Ho Cheung Shum

Beating the jetting regime Beating the jetting regime

Keywords: microfluidics, droplet breakup, numerical simulation

We study numerically the dynamics of jets and drops in a microcapillary co-flow device. The co-flow stream encounters different flow regimes, including dripping, jetting. Using a level-set/finite element axysimmetric numerical simulation, we study the dynamics of breakup of a jet subject to flow rate perturbations. A scaling law for the width of the unperturbed jet is presented and compared to existing experimental results as well as numerical measurements. Then, we show that the introduction of a sinusoidal perturbation of the inner fluid flow rate can facilitate breakup of the fluid in regimes where a jet is usually observed. Moreover, the flow rate perturbation leads to a good control over the size and the frequency of the resulting droplets. Using theoretical modelling, we provide a criterion to determine the optimal frequency to break up the jet. We also derive scaling laws to determine the volume of inner fluid encapsulated in the emulsion droplets as a function of the frequency and to estimate the distance for the jet to break up as a function of the amplitude of perturbation. These scaling laws are in good agreement with results of numerical simulations. Our work suggests a novel approach and offers guiding principles to break up liquid jets in cases where dripping is difficult to achieve.

I. INTRODUCTION

The dynamics of a jet of one dispersed phase injected in a co-flowing continuous phase constitutes an interesting problem in fluid mechanics and is of fundamental importance in industrial applications that involve emulsions and bubbles. Indeed, a particular use of this system is to generate emulsions, which consist of droplets of one inner phase dispersed in a second immiscible fluid, called outer or continuous phase. The resulting multiphase flow and droplets are widely used in the industries of food, pharmaceutics, cosmetics and bio-medicine. Recent advances in microfluidics have enabled generation of emulsions with more sophisticated geometries. For instance, multiple emulsion can now be formed with a high degree of control over the size uniformity and a perfect encapsulation efficiency [START_REF] Utada | Monodisperse double emulsions generated from a microcapillary device[END_REF] .

To generate monodisperse emulsion, a proven microcapillary device consists of two capillaries aligned coaxially. Improvements in the generation of simple emulsion using this method cannot be done without a good understanding of fluid dynamics at a microscopic scale. We can separate the dynamic of co-flow streams between dripping and jetting. In the dripping regime, drops are directly formed at the tip of the inner capillary and the resulting emulsion drops are monodisperse in size. In the jetting regime, a jet is formed at the tip of the inner capillary; the jet often breaks up into droplets at certain distance from the inner capillary tip. The size of the resulting droplets is usually non-uniform and different from the size of the inner tip [START_REF] Utada | Dripping to jetting transitions in coflowing liquid streams[END_REF][START_REF] Utada | Absolute instability of a liquid jet in a coflowing stream[END_REF] . The mechanism of the dripping and jetting regimes are quite well understood; in particular the transition between these two regimes is well described by a linear stability analysis of the Rayleigh-Plateau instability in a co-flow stream [START_REF] Guillot | Stability of a jet in confined pressure-driven biphasic flows at low reynolds number[END_REF][START_REF] Guillot | Stability of a jet in confined pressure-driven biphasic flows at low reynolds number in various geometries[END_REF] . The most unstable wavelength predicted by this analysis is in agreement with experimental results and in particular the transition from dripping to jetting is wellpredicted. A number of groups attempt to produce droplets of a controlled size in the jetting regime. A useful tool to study this flow dynamics is numerical simulation, which enables simulation of the different regimes encountered in co-flow device with arbitrary geometries [START_REF] Hua | Numerical simulation of microdroplet formation in coflowing immiscible liquids[END_REF] . Moreover, different physical values can be extracted from the numerical simulations, such as the pressure at a given point of the flow or the velocity field, which are often difficult to measure in the experiments.

A jet of a dispersed fluid can break up into droplets due to a Rayleigh-Plateau instability [START_REF] Rayleigh | On the capillary phenomena of jets[END_REF][START_REF] Plateau | Sur les figures d'équilibre d'une masse liquide sans pesanteur[END_REF][START_REF] Eggers | Physics of liquid jets[END_REF] . The cylindrical jet is rendered unstable by capillarity effects; this instability is influenced by parameters such as viscosity of the outer fluid or inertia. A dynamic description of the system suggests a dependence of the destabilization of the jet on the frequency. Depending on the nature of the instability, convective or absolute, a dripping or a jetting regime prevails [START_REF] Guillot | Stability of a jet in confined pressure-driven biphasic flows at low reynolds number[END_REF] . The dynamics of co-flow jet has been studied experimentally [START_REF] Utada | Dripping to jetting transitions in coflowing liquid streams[END_REF][START_REF] Utada | Absolute instability of a liquid jet in a coflowing stream[END_REF] and the transition between dripping and jetting has been shown to depend on two dimensionless numbers. The Weber number of the inner phase is defined as:

We in = ρ in r jet u in 2 γ (1)
where ρ in and u in denote, respectively, the density and the average velocity of the dispersed phase, r jet denotes the radius of the jet and γ is the interfacial tension between the two phases. We in represents the ratio between the inertial and interfacial tension effects. The dynamic of the flow also depends on the capillary number of the continuous phase

Ca out = η out u out γ , (2)
where η out and u out are, respectively, the dynamic viscosity and the average velocity of the continuous phase. It represents the balance of viscous and interfacial tension effects. Experimentally, it has been observed that for We in < 2 and Ca out < 0.2, dripping occurs 10 . However, numerical results suggest that the transition also depends on the viscosity ratio of the two phases η in /η out and occurs at smaller dimensionless numbers, typically We in < 0.7 and Ca out < 0.1 [START_REF] Homma | Numerical simulation of the formation of a jet and droplets from a capillary in a co-flowing ambiant fluid[END_REF] .

In this work, we study numerically the formation of drops in devices with the co-flow geometry under conditions where a jet is normally observed. The typical parameters used in the systematic study of the flow rate perturbation are We in ∼ 3 and Ca out ∼ 0.4. This paper constitutes a first step towards the use of induced perturbation to modify the flow. A good understanding of the fluid dynamics governing the differents regimes of flow in a coflow device is essential to increase the efficiency and the modularity of the microfluidics technology for future industrial applications. It is also instrumental to some recent applications, such as generation of water-in-water emulsions, [START_REF] Ziemecka | Monodisperse hydrogel microspheres by forced droplet formation in aqueous two-phase systems[END_REF][START_REF] Shum | Microfluidic fabrication of water-in-water (w-w) jets and emulsions[END_REF] which have provided a new platform for studying the generation of droplets. As the surface tension of these water-water systems is very low, typically 1000 times lower than in oil-water emulsion, the generation of droplet is impossible by classical methods. The values of the dimensionless numbers We in and Ca out observed in these systems typically lie in the range used in this paper. At even lower values of these dimensionless numbers, the initial perturbation on the interface sometimes lead to a corrugated jet as observed in three-phases flow [START_REF] Shum | Corrugated interfaces in multiphase core-annular flow[END_REF] . A method to break the inner jet into droplets relies on pressure perturbation of the dispersed phase, leading to a monodisperse emulsion. However, the dynamics of this mechanism is not clearly understood from a fluid mechanical point of view. Thus, it is of interest to investigate the mechanism of droplet formation and the effects of the different parameters on the resulting size of the droplets as well as the flow rates that can be used for industrial applications.

The paper is organized as follows. In section II, the governing equations and the mathematical description of the system as well as the relevant different dimensionless parameters are presented. Section III is devoted to the numerical simulation: we introduce the numerical procedure used for solving the two-phase coflow problem and we validate the method by comparison with previous results obtained in the literature. Then, in section IV, we show that the radius of the unperturbed jet can be obtained from a simple scaling law, which is confirmed by our numerical simulation. In section V, we focus on the dynamic of the perturbed jet: a scaling law for the resulting size of the droplets is derived together with a condition on the amplitude of disturbance.

II. GOVERNING EQUATIONS AND MATHEMATICAL DESCRIPTION

A. Description of the system

We consider a microcapillary coflow geometry similar to the one used in previous experimental studies, which consists of two coaxially aligned capillary tubes [START_REF] Utada | Dripping to jetting transitions in coflowing liquid streams[END_REF][START_REF] Utada | Absolute instability of a liquid jet in a coflowing stream[END_REF][START_REF] Guillot | Stability of a jet in confined pressure-driven biphasic flows at low reynolds number[END_REF] . The inner capillary tube is cylindrical with a tip tapered to an inner diameter of typical size d tip = 20 µm and outer diameter d out = 30 µm. The inner tube is placed in a surrounding square capillary with an inner diameter equal to the outer diameter of the untapered portion of the inner capillary; this ensure a good coaxial alignment 2 . We use another cylindrical outer capillary of inner diameter d cap = 200 µm axially aligned in the square capillary as a collecting tube. We assume all the capillaries to be axisymmetrically aligned to simplify the problem and the subsequent numerical simulation.

The experimental device leads to the numerical modelization represented schematically in figure 1. An inner fluid is injected at the center in an axisymmetric capillary tube of inner diameter r in and outer diameter r out at the constant mean velocity U in (r). This inner fluid is surrounded by an outer phase injected through the axisymmetric coaxial annulus of external radius R and internal radius r out at the mean velocity U out (r). Regarding the order of magnitude of the different forces, gravitational effects are negligible due to the small dimensions of the channel in this study and consequently is not taken into account in the numerical simulation. The length of the computational domain is L = 1000 µm which represents 100 r in and we ensure that this is sufficiently long for all the numerical computations presented in this article. Both fluids are assumed to be incompressible with a constant density of ρ in and ρ out for the inner and outer fluids repectively. Moreover, in this study, all the fluids are newtonian and have a constant viscosity η in and η out respectively. Finally, γ is the interfacial tension between the two phases. As the two fluids are not miscible, an axisymmetric interface is formed between the two phases. We solve the axisymmetric, incompressible, unsteady, fully developed Navier-Stokes equations in cylindrico-polar coordinates (r, θ, z). The velocity field in this system of coordinate is written as (u, v, w) where u, v, w denote, respectively, the radial, azimuthal and axial velocity. Moreover, due to the symmetry of the system, we have v = 0 and thus we will use only the couple of coordinates (u, w).

B. Governing equations

Both inner and outer fluids are considered as incompressible. The Navier-Stokes equation and the continuity equation with variable density and viscosity are expressed as follows:

ρ ∂u ∂t + (u • ∇)u = -∇p + ∇ • η ∇u + ∇u T +∇ • γ I -(nn T ) δ (3)
∇ • u = 0 ( 4
)
where u = (u, v, w) is the velocity of the fluid, p the pressure, ρ the density, η the dynamic viscosity, γ the surface tension, I the identity matric, δ is a Dirac delta function that is zero everywhere except at the fluid-fluid interface and n is the unit normal vector to the interface defined by

n = ∇φ |∇φ| ( 5
)
where φ is the level-set function and will be defined more precisly in the next section. We consider a two-phase flow system, and the two fluids are incompressible. It follows that the normal velocity at the interface is continuous and is expressed as

u in • n = u out • n (6)
where the subscript in and out denote, respectively, the inner and outer fluids. Then, the boundary conditions for the stress tensor at the interface writes, for the normal and tangential stress respectively,

-p + η ∇u + ∇u T • n = γ C n (7)
η ∇u + ∇u T • t = 0 (8)
where C represents the curvature of the interface and t the unit tangent vector at the interface.

On the solid boundaries of the capillary walls, denoted as Σ, the no-slip condition implies

u Σ = 0 (9)

C. Inlet and outlet flow

Using the cylindrico-polar coordinate (r, θ, z) where the velocity field is expressed as (u,v,w), the inner fluid is injected through the inner capillary at the mean velocity w m,in , which is related to the inner flow rate Q in by the relation

w m,in = Q in π r in 2 (10)
We implement in the numerical model the mean velocity w m,in as the inner inlet laminar flow. Then, the inner fluid at the outlet of the capillary adopts a Poiseuille profile, which can be expressed as

   u in (r) = 0 w in (r) = 2 Q in π r in 2 1 - r 2 r in 2 for r ∈ [0; r in ] (11)
The inner fluid velocity profile satisfies the no-slip condition at r = r in . The inner flow rate is

Q in = rin 0 2π r w in (r) dr (12)
The average velocity in the outer annular channel w m,out is given by

w m,out = Q out π (R 2 -r out 2 ) ( 13
)
For the outer inlet laminar flow, we implement in the numerical model the mean velocity w m,out . The outer fluid in the outer annular channel has a more complicated profile given by

             u out (r) = 0 w out (r) = 2 Q out π (R 2 -r out 2 ) 1 - r 2 R 2 - 1 -r out 2 /R 2 ln (R/r) ln R r out 1 + r out 2 R 2 - 1 -r out 2 /R 2 ln(R/r out ) for r ∈ [r out ; R] (14)
which satisfies the no-slip boundary conditions at r = r out and r = R. The outer external flow rate is directly obtained by the relation:

Q out = R rout 2π r w out (r) dr (15)
We can compare the analytical velocity profile with the numerical profile obtained at the exit of the capillary by using a mean value of the velocity for the inlet flow. The comparison is plotted in figure 2, where a good agreement confirms the validation of our approach to implement the inlet flow.

Moreover, as the whole problem is axisymmetric, we solve it only in a half part of the domain presented in figure (1) and the axis r = 0 is taken as the axis of axisymmetry. Furthermore, an outflow boundary condition is applied at the exit of the capillary located at the right of the domain by imposing the condition of non-viscous stress together with p = 0.

D. Dimensionless parameters

Regarding the parameters present in this problem, we can define dimensionless parameters to describe the resulting flow. These parameters are used to compare the various forces in the problem. First, the ratio of the inertia force to viscous force is characterized by the Reynolds number Re defined, respectively, for the inner and outer phase by ) FIG. 2. Comparison of theoretical axial velocity profile given by the relations ( 11) and ( 14) (red continuous line) and the numerical profile (black dotted line) for an inner fluid and outer fluid flow rates of Qin = 500 µL•h -1 and Qout = 20 000 µL•h -1 , respectively. The axial velocity is taken directly at the exit of the inner capillary.

Re in = 2 ρ in w m,in r in η in and Re out = 2 ρ out w m,out R η out (16)
where the prefactor of two account for the typical length of the problem, which is the diameter of the capillaries and is equal to 2 r in and 2 R respectively. We can also define the ratio of inertia force to the surface tension effects using the Weber number

We in = 2 ρ in w m,in 2 r in γ and We out = 2 ρ out w m,out 2 R γ (17)
Finally, the capillary number is defined by the ratio of viscous force to surface tension force:

Ca in = ρ in w m,in γ and Ca out = ρ out w m,out γ (18)
We can also define the ratio of the density, ρ in /ρ out , and the ratio of dynamic viscosity, η in /η out .

III. NUMERICAL PROCEDURE

We elucidate the dynamic of the two-phase flow in the microcapillary coflow system using a level set method with a finite-element numerical scheme implemetend in a commerical software, COMSOL Multiphysics c . We solve the laminar, unsteady, two-phase flow with a conservative level set method, which is used to describe the fluid interface between the two phases [START_REF] Olsson | A conservative level set method for two phase flow[END_REF][START_REF] Yue | A diffuse-interface method for simulating two-phase flows of complex fluids[END_REF] . All simulations are perfomed using a numerical grid composed of triangular elements. All elements are of standard Lagrange P 1 -P 2 type (i.e. linear for the pressure field and quadratic for the velocity field). The temporal solver is IDA solver (Implicit Differential-Algrebic solver), based on a backward differencing formulas [START_REF] Hindmarsh | Sundials: Suite of nonlinear and differential/algebraic equation solvers[END_REF] . At each time step, the system is solved with the sparse direct linear solver PARDISO 18 . In this work, we are using relatively high frequency, which ranges from 600 Hz, up to 20 000 Hz. The time steps used in the simulation goes up to T = 5 × 10 -5 s, which is much smaller than the inverse of the applied frequency. This ensures good convergence and high precision of the numerical simulations. As such, the numerical simulation is carried out for a sufficiently long time until a stable droplet breakup time is reached, typically around t f inal = 0.05. An example of the time step obtained in the simulation at a frequency f = 20000 Hz is shown in figure 3.

A. Level-set function φ

The interface between the two fluids is defined by a level-set function, φ. This function is such that in the inner phase its value is φ = 0 whereas in the outer phase, φ = 1. The function φ goes smoothly from 0 to 1 in a transition layer of thickness . In all the simulations presented in this paper, we have used = h c /2 ∼ 10 -6 where h c represents the typical mesh size. The interface between the inner and outer phases is defined by φ = 0.5 and is advected at the fluid velocity u. The convection of the reinitialized level set function is described by the equation: where describes the thickness of the transition layer. The parameter λ determines the amount of reinitialization and its value is of the order of the maximum magnitude occurring in the velocity field.

∂φ ∂t + ∇ • (φ u) + λ ∇ • φ (1 -φ) ∇φ |∇φ| -∇φ = 0 ( 19
)
The aim of the level set function, apart from describing the interface, is to give the density and the viscosity of the fluid at each points through the relations:

ρ = ρ in + (ρ out -ρ in ) φ (20) η = η in + (η out -η in ) φ (21)

B. Meshgrid

In all simulation, we have used an Eulerian meshgrid which consists of triangular elements (see figure 4).

We have performed a simulation with different number of degrees of freedom (DoF) to test the convergence of the numerical results for dimensionless parameters We in = 1.5323 and Ca out = 0.0904. The results are similar for 143 982 DoF and higher, as shown in figure 5. This ensures the convergence of our results at 143 982 DoF. To compare quantitatively the influence of the number of DoF, we measure the distance for the jet to break up into drops as a function of the number of degree of freedom (see figure 6). Using 143 982 DoF for the numerical simulation ensure a well-converged simulation (around 1% of precision for the breakup distance). This number of DoF will be used in the following simulations; this provides the needed accuracy while keeping a reasonable computational time (i.e. about 4 h on a standard workstation).

C. Numerical validation

Even though our numerical simulation is well-converged, the results given by our numerical model should also agree with experimental data. To check this, we have performed a systematic study to determine the dripping/jetting transition in our system and compare these results to theoretical and experimental results obtained by Utada et al. for a similar system 2 . The results are summarized in figure 7.

We observe that the dripping/jetting transition occurs for inner Weber number, We in > 0.75 and outer capillary number, Ca out > 0.11. These two values are lower than the one obtained by Utada et al. However, this discrepancy has already been observed by Homma et al [START_REF] Homma | Numerical simulation of the formation of a jet and droplets from a capillary in a co-flowing ambiant fluid[END_REF] and is attributed to the difference in the viscosity ratio. Therefore, the transition obtained with our numerical model is in good agreement with previous values obtained in the literature.

IV. SIZE OF THE UNPERTURBED JET

We consider the axisymmetric co-flow device presented in figure 1. If no perturbation is applied and the Weber number of the inner fluid and capillary number of the outer fluid are sufficiently large, We in > 2 and Ca out > 0.2 respectively, we do not observe formation of drops in a dripping regime and instead a jet is formed. Here, we evaluate the width of the undisturbed jet as a function of the different parameters including the inner and outer fluid flow rates, the viscosity of both phases, the interfacial tension and the characteristic size of the device.

The typical Reynolds number, which describes the ratio of inertial forces to viscous forces is quite low, Re = O(1). Therefore, the Navier-Stokes equation is reduced in the stationary regime to

η ∇ 2 u = ∇p (22)
where u and p are, respectively, the velocity and the pressure of the fluid in the considered phase. The symmetry of the problem imply that the velocity field can be expressed as

u out = w out (r) e z and u in = w in (r) e z (23)
and the pressure depends only on the axial coordinate:

p out (z) and p in (z) (24)
On the outer boundaries r = R, we apply the no-slip boundary condition as it remains an excellent approximation for flows at scales above micrometers [START_REF] Stone | Engineering flows in small devices: Microfluidics toward a lab-on-a-chip[END_REF] :

w out (r = R) = 0 (25)
At the interface between the two fluids, the continuity of the velocity field is required:

w in (r = r j ) = w out (r = r j ) (26)
together with the continuity of the tangential shear stress

η in dw in dr r=rj = η out dw out dr r=rj (27)
From these boundary conditions, the equation ( 22) leads to the following expressions for the velocity of the dispersed phase:

w in (r) = - 1 4 η in dp in dz (r 2 -r j 2 ) + (R 2 -r j 2 ) 4 η out dp out dz - r j 2 2 η out dp in dz - dp out dz ln r j R ( 28
)
and for that of the continuous phase:

w out (r) = - 1 4 η out dp out dz (r 2 -R 2 ) - r j 2 2 η out dp in dz - dp out dz ln r R (29)
Therefore, the fluid flow rates of both phases, respectively, writes

Q in = 2 π r j 4 16 η in dp in dz - r j 4 4 η out dp in dz - dp out dz ln r j R + R 2 -r j 2 8 η out r j 2 dp out dz ( 30
)
Q out = π 8η out dp in dz - dp out dz 4r j 4 ln r j R + 2r j 2 (R 2 -r j 2 ) +(R 2 -r j 2 ) 2 dp out dz (31)
The two pressure gardients can be expressed using the continuity of the normal stress at the interface

p in -p out = 2 γ r j (32)
In the unperturbed case, the radius of the jet, r j , remains constant and a differentiation with respect to z leads to

dp in dz = dp out dz = dp dz (33)
Thus, the expressions (30) and (31) becomes

Q in = 2 π r j 4 16 η in + R 2 -r j 2 4 η out r j 2 dp dz (34) Q out = π 8η out (R 2 -r j 2 ) 2 dp dz ( 35
)
The ratio of the inner and outer fluid flow rates writes

Q in Q out = 2 r j 2 (R 2 -r j 2 ) + η out η in r j 4 (R 2 -r j 2 ) 2 (36)
Thus, we arrive at the radius of the inner jet

r j = R      1 + Q in Q out η in η out 1/2 -1 + Q in Q out η in η out -2 - Q in Q out      1/2
(37) Surprisingly, we see that the scaling law suggested by Utada et al. nearly collapses on the present scaling law. Most of the experimental measurements were done at a low ratio of the inner and outer fluid flow rates, Q in /Q out 1. Considering the expression (37) in this limit lead to an asymptotic scaling law:

r j = R Q in 2 Q out 1/2 (38)
which is exactly the scaling law derived by Utada et al. In this case, the influence of the viscosity of both phases is negligible.

We have also performed numerical simulations in the regime where an inner jet is produced for a different viscosity ratio, η in /η out = 1.25. The resulting measurements collapse well with the scaling law with no adjustable parameters as shown in figure (9).

The scaling law (37) covers the range of viscosity ratio observed by Utada et al. [START_REF] Utada | Dripping to jetting transitions in coflowing liquid streams[END_REF] and is more general as it does not require a low ratio of the inner and outer fluid flow rates. The scaling law is expected to be valid for any flow rates as well as any ratio of viscosity.

V. DROP FORMATION

To apply a perturbation on the dispersed flow, we introduce the inlet velocity as follows:

w in (t) = w m,in [1 + cos(2πf t)] (39)
where and f are respectively the amplitude and the frequency of perturbation and w m,in the average inlet velocity. Depending on and f , different flow morphologies are obtained.

A. Influence of the frequency

It is known since the work of Rayleigh and Plateau that a jet breaks up into drops due to the interfacial tension effects for a given wavelength of the perturbation [START_REF] Rayleigh | On the capillary phenomena of jets[END_REF][START_REF] Plateau | Sur les figures d'équilibre d'une masse liquide sans pesanteur[END_REF][START_REF] Eggers | Physics of liquid jets[END_REF] . In addition to the spontaneous breakup of the jet at a wavelength chosen by the system, we can facilitate the generation of droplets using an external forcing of the flow at a given frequency. Here, we study the influence of the frequency at a constant amplitude of perturbation, = 0.1 and tune the frequency f in the range [6 000 Hz; 20 000 Hz]. Outside range of frequencies, no visible effects are observed. Indeed, it is well-known, from the mechanism responsible for the Rayleigh-Plateau instability, that for perturbation smaller than the radius of the jet, no perturbation can grow and thus the jet does not break up.

For a perturbation frequency of f = 16 000 Hz, a jetting regime is observed with the production of droplets. We have observed in our numerical simulation that in this case, the generated droplets are monodisperse in size. The morphology of the jets at three consecutive time steps is represented in figure 10(a),(b),(c). The perturbation on the interface grows along the jet and the radius of the jet becomes periodically thinner with a wavelength λ = 2 π/k = 2π v/f , as shown in figure 10(a). Then, when the radius becomes very thin, a droplet is generated as exhibited in figure 10(c).

For a lower frequency, f = 10 000 Hz, the dynamics of the breakup remains similar except for the longer pinch-off length due to the larger wavelength λ = 2 π v/f . Therefore, the initial breakup occurs for a mother droplet (see figure 11(a),(b)). An elongated thinner fluid stream later forms a smaller drop called satellite droplet (figure 11(c),(d)).

Further decreasing the frequency of perturbation leads to the formation of a mother drop followed by two satellites droplets (figure 12(a),(b),(c),(d)).

For given inner and outer fluid flow rates, as well as surface tension and viscosities, the dynamics of the jet is controlled by two key parameters: the frequency and the amplitude of perturbation. We have performed a systematic numerical study to determine the distance at which the jet breakup into droplets as a function of these two parameters and f . The results are plotted in figure (13). Increasing the amplitude of the initial perturbation leads to a smaller distance for the breakup. An optimal frequency is around f = 16 000 Hz in this study but depends on the velocity of the jet. Different theories have been developed to describe the break up of a jet into droplets. In the earlier studies of Rayleigh [START_REF] Rayleigh | On the capillary phenomena of jets[END_REF] and Plateau 8 , an ideal case of an inviscid jet in air is considered; thus the effect of the surrounding fluid is neglected. Tomotika has shown that a surrounding viscous fluid has a non trivial effect on the stability of the jet and consequently, the optimal frequency, where the growth rate of the instability is maximum 20 , is modified. In the present study of a coflowing stream in a microcapillary device, both the effect of a surrounding outer fluid and a non-zero velocity of both phases have to be considered. Therefore, we use the theory developed by Guillot et al., [START_REF] Herrada | Spatiotemporal instability of a confined capillary jet[END_REF] which considers the dependence of the breakup on the wavenumber of the perturbation, k, for a co-flow jet. They found that for a perturbation of the form,

r = r 0 (1 + e iωt ) (40)
a dispersion relation between ω and k can be derived where

ω = -i K a E(x) x 2 k + F (x) x 5 (k 2 -k 4 ) (41)
E(x) = d(x) b(x) -a(x) e(x) a(x) + d(x) F (x) = c(x) d(x) -f (x) e(x) a(x) + d(x) with a(x) = -(1 -x 2 ) 2 b(x) = 4 x (1 -x 2 ) c(x) = 2 x 4 -2 x 2 -4 x 4 ln(x) d(x) = 2 x 2 (x 2 -1) - η out η in x 4 e(x) = 8 x 3 -4 η out η in x 3 -4 x f (x) = 4 x 4 ln(x) - η out η in x 4
The capillary number based on the capillary radius and the ratio of radii are defined, respectively, by

K a = -δ z P R 2 γ (42) x = r j R (43)
These two dimensionless numbers are related to the inner and outer fluid flow rates through the relations

Q out = π γ R 2 8 η out 1 -x 2 2 K a ( 44
)
Q in = Q out η out η in x 4 (1 -x 2 ) 2 + 2 x 2 1 -x 2 (45)
Therefore, for an imposed wavenumber k related to f through the relation

k = 2 π v f (46)
The growth rate of the perturbation is given by the real part of ω = ω r : The growth rate ω r as a function of the perturbation frequency is plotted in figure (14). It has been shown using energetic consideration that modes whose wavelength 2π/k are smaller than the perimeter of the jet have a negative difference in energy between the perturbed and the unperturbed jet. Therefore, no perturbation with wavenumber k > 1/r 0 can grow [START_REF] Eggers | Physics of liquid jets[END_REF] . This is in agreement with our numerical simulation where no growth of perturbation for f > 23 000 Hz can be observed; this corresponds to k r jet > 1 (see also figure (14)). The optimal frequency of perturbation to break up the jet into drops at the present flow rate (Q in = 1 800 µL • h -1 and Q out = 80 000 µL • h -1 ) is f 16 000 Hz, in agreement with the results from numerical simulation (see figure 13). This theoretical description is thus welladapted to describe the destabilization of a co-flow jet and therefore provides a criterion to choose experimentally the optimal frequency to break up a co-flow jet into droplets using induced perturbations.

ω r = F (x) x 5 (k 2 -k 4 ) (47)

B. Scaling of the droplets size

For given inner and outer fluid flow rates, the generation of monodisperse droplets is possible in a range of perturbation frequency. For Q in = 1 800 µL • h -1 and Q out = 80 000 µL • h -1 , this is achievable for f ∈ [6 000; 18 000] Hz as exhibited in figure (13). However, the size of the resulting drops together with the number of satellite droplets varies with the applied frequency (see figures 10, 11, 12). The inability to predict the number and volumes of the satellite droplets limit the ability to estimate the exact radius of the main emulsion droplets. Therefore, we consider the total volume of fluid contained in all drops as

V d = 4 3 π r drop 3 + i r sat,i 3
where the first term accounts for the mother drops and the second term for the satellite droplets. In some cases, more than one satellite droplets are present, the symbol i denotes the particular satellite droplet and its radius r sat,i . The volume of the jet, which is a sum of the volume of all satellite droplets, writes

V jet = π r jet 2 w jet f = V d (49)
where π r jet 2 is the transverse area of the jet and w jet /f is the periodic distance over one perturbation. Moreover, the inner flow rate can be expressed as:

Q in = π r jet 2 w jet (50)
leading to a simple relation for the volume of one main drop and its satellite droplets: This relation shows that the encapsulated volume of fluid does not depend on the surface tension γ or on the viscosity. We measure the radius of the different drops obtained in the numerical simulation for inner and outer flow rates of Q in = 1 800 µL • h -1 and Q out = 80 000 µL • h -1 respectively, with perturbation with an amplitude of = 0.1 and frequencies in the range f ∈ [6000; 18000]. The resulting measurements show a very good agreement with the scaling law derived without any adjustable parameters, as shown in figure 15.

V d = 4 3 π r drop 3 + r sat 3 = Q in f (

C. Influence of the amplitude of perturbation

For a given perturbation frequency, the breakup length show a clear dependence on the initial perturbation amplitude (see figure 13). Considering an initial perturbation of the jet radius such that: r jet (t) = r jet 1 -0 e ω t (52) the jet can be considered as broken up into droplets at t 0 for r jet (t 0 ) ∼ 0. This leads to a condition on 0 0 e ω t0 ∼ 1 (53)

The distance d break for the jet to break up is related to the time t 0 and to the velocity of the jet w jet through the relation t 0 = d break /w jet . It follows a relation between the initial disturbance amplitude 0 and the distance d break :

d break ∼ - v jet ω r ln 0 (54)
For given inner and outer fluid flow rates together with fixed interfacial tension and viscosities, v jet and ω r remain constant and thus a simple scaling law is obtained. We measure, in the numerical simulation, the distance d break for a perturbation frequency f = 16 000 Hz, Q in = 1 800 µL • h -1 and Q out = 80 000 µL • h -1 . The results from the numerical simulation agree well with the scaling law, as plotted in figure (16).

VI. CONCLUSION

In this paper, we study the dynamics of a co-flow jet in the regime where usually the jet does not break up into monodisperse droplets. Using an axisymmetric numerical simulation, we show that the introduction of an initial perturbation to the velocity of the dispersed phase can lead to a faster breakup of the jet with a good control over the size of the droplets. A scaling law is derived for the radius of the unperturbed jet in a co-flow stream; the prediction from the scaling law agrees well with previous experimental results in the literature and numerical measurements. The two controlable parameters, the frequency f and the amplitude of the perturbation, are shown to have a strong effect on the flow. A scaling law for the resulting volume of dispersed phase contained in the droplets is provided and the distance to break the jet is characterized with respect to the frequency. The new results lead to a better understanding of the fluid dynamics of a perturbed jet; this is particulary helpful in the new field of water-in-water emulsion where a perturbation in the pressure or velocity field must be introduced to form monodisperse drops. We show that the optimal frequency can be estimated using theoretical consideration and can be directly applied to the numerical simulations or in experiments.

By providing guidelines on how to break up jets into monodisperse droplets with controlled sizes, our work suggests new strategies to extend the applicability of droplet microfluidic techniques to a wider variety of fluidic systems.

FIG. 1 .

 1 FIG. 1. (Color online) Schematic of the problem. The inner fluid is injected from a capillary at the center of the domain in a coflowing outer fluid. This results in a jet that may break up into droplets downstream in the microchannel.

FIG. 3 .FIG. 4 .

 34 FIG. 3. Timestep for a simulation at f = 20 000 Hz, Wein = 2.5, Caout = 0.317. The time step is always smaller than the period of the applied forcing

FIG. 5 .FIG. 6 .

 56 FIG. 5. Comparison of the shape of the interface for Wein = 1.5323 and Caout = 0.0904 and a viscosity ratio ηin/ηout = 0.1 at the same time t = 0.082 showing a jetting regime (a) for 17 333 DoF, (b) for 36 919 DoF, (c) for 74 587 DoF, (d) for 99 977 DoF, (e) for 143 982 DoF and (f) 218 559 DoF. The numerical results look similar beyond 74 587 DoF, which is sufficient for ensuring a good convergence of our numerical simulation using 143 982 DoF.

FIG. 7 .

 7 FIG. 7. State diagram showing the dripping/jetting transition as a function of Caout and Wein for a viscosity ratio ηin/ηout = 0.1. Red circles indicate conditions leading to the jetting regime and blue squares indicate the dripping regime.

FIG. 8 .FIG. 9 .

 89 FIG.8. Diameter of the inner jet rescaled by the diameter of the outer capillary. The black line indicates the scaling law (37) with a viscosity ratio ηin/ηout = 0.1. The blue circles are the experimental measurements obtained by Utada et al.[START_REF] Utada | Dripping to jetting transitions in coflowing liquid streams[END_REF] together with their suggested scaling law (red dotted line).

FIG. 10 .

 10 FIG. 10. Morphology of the jet at three successive times for Wein = 2.533, Caout = 0.3168, ηin/ηout = 1.25 and γ = 20 mN•m -3 . The amplitude and the frequency of perturbation are = 0.1 and f = 16 000 Hz, respectively.

FIG. 11 .FIG. 12 .

 1112 FIG. 11. Morphology of the jet at three successive times for Wein = 2.533, Caout = 0.3168, ηin/ηout = 1.25 and γ = 20 mN•m -3 . The amplitude and the frequency of perturbation are = 0.1 and f = 10 000 Hz, respectively.

FIG. 13 .

 13 FIG.13. Length of the jet before breakup (in µm) as a function of the perturbation frequency for different amplitude: = 0.1 (in blue), = 0.2 (in green), = 0.4 (in red), = 0.6 (in magenta), = 0.8 (in cyan). The inner and outer fluid flow rates are Qin = 1 800 µL • h -1 and Qout = 80 000 µL • h -1 , respectively.

! (s - 1 )FIG. 14 .

 114 FIG.14. Growth rate of the Rayleigh-Plateau instability ωr (in s -1 ) for an imposed perturbation frequency f (in Hz). The dashed-dotted line located around f = 23 000 Hz indicates the critical frequency beyond which no perturbation can grow. The black dotted line around f = 16 000 Hz indicates the frequency where the growth rate is maximum.

)FIG. 15 .

 15 FIG.15. Volume encapsulated in the mother and satellites drops V d (in µL 3 ) as a function of the perturbation frequency f (in Hz). The red squares are measurements using numerical simulation together with the uncertainties and the dashed black line is the scaling law (51) without any adjustment parameters. The inner and outer fluid flow rates are Qin = 1 800 µL • h -1 and Qout = 80 000 µL • h -1 , respectively.

0 FIG. 16 .

 016 FIG.16. Distance for the jet to break up into drops, d break (in µm), as a function of the initial amplitude of perturbation 0 (red squares) together with the scaling law (54) (dotted line).