Iaroslava Golovkova

Lorraine Montel

Elie Wandersman

Thibault Bertrand
email: t.bertrand@imperial.ac.uk

Alexis Michel Prevost

Lea-Laetitia Pontani
email: lea-laetitia.pontani@sorbonne-universite.fr

Depletion attraction impairs the plasticity of emulsions flowing in a constriction †

We study the elasto-plastic behavior of dense attractive emulsions under mechanical perturbation. The attraction is introduced through non-specific depletion interactions between the droplets and is controlled by changing the concentration of surfactant micelles in the continuous phase. We find that such attractive forces are not sufficient to induce any measurable modification on the scalings between the local packing fraction and the deformation of the droplets. However, when the emulsions are flowed through 2D microfluidic constrictions, we uncover a measurable effect of attraction on their elasto-plastic response. Indeed, we measure higher levels of deformation inside the constriction for attractive droplets. In addition, we show that these measurements correlate with droplet rearrangements that are spatially delayed in the constriction for higher attraction forces.

Introduction

The flow of particulate systems is a problem of great importance both theoretically and practically, with direct applications to the industry. It is relevant for a wide range of soft materials, from granular packings to foams and emulsions. While these materials present obvious differences, they share universal features, e.g. they generically undergo what is known as a jamming transition 1,2 . As the particle or droplet volume fraction φ increases, this rigidity transition between liquid and amorphous-solid states controls the phase behavior of these disordered solids. At a critical volume fraction φ c (random close packing), the system jams and develops a yield stress [3][4][5][6] . The mechanical and rheological properties, such as the elastic modulus or the local pressure, of these systems are known to display a power law dependence with the distance to the jamming onset (φφ c ) 4,5,[7][8][9][10][11][12] .

Jammed solids are characterised by a spatially heterogeneous network of interparticle contacts, with a broad distribution of contact forces exhibiting an exponential tail 4,11,13,14 in which only a small subset of the particles sustain most of the mechanical load [15][16][17][18] . Below the yield stress, these systems responds elastically, while above it, they deform and flow plastically 19 . In these soft glassy flows, it was shown that stress and strain rates are coupled nonlocally 6,20,21 . In two-dimensional materials, the flow properties can easily be probed both at the microscopic and macroscopic scales [22][23][24][25][26][27][28][29][30][31] . As a consequence, previous experimental studies examined the microscopic rearrangements in a variety of two-dimensional model systems under stress 32,33 . This plastic flow is generically governed by local structural rearrangements which relieve stresses and dissipate energy 6,22,23,34 . Local plastic rearrangements have been connected to the fluctuating macroscopic flow in both simulations [35][36][37][38][39] and theoretical studies 6,20,21,40,41 of model systems. Nevertheless, the intimate link between the microscopic dynamics of an amorphous material and its macroscopic elasto-plastic response is still an open question for a broad class of more realistic materials.

In emulsions, the use of surfactants prevents the coalescence of the droplets and leads to short-range purely repulsive dropletdroplet interactions 22,23,42 . As such, dense stable emulsions are examples of jammed solids. In the last decades, a number of experimental works studied the structural, mechanical and rheological properties of purely repulsive emulsions 12,[42][43][44][45][46][47] . In particular, as in other soft materials 29,31,[48][49][50][51][52] , recent studies in quasi-2D flowing emulsions have also highlighted the importance of T1 events for local rearrangements and stress redistribution 22,23 . Monodisperse emulsions allow one to study material properties such as grain boundaries, dislocations and plasticity [53][54][55][56][57] ; in particular, a recent study showed the existence of a spatiotemporal periodicity in the dislocation dynamics of these emulsions 33 . However, none of these studies have so far adressed the question of how interdroplet attractive forces modify the flow response of these emulsions. Indeed, in a variety of natural settings and industrial applications, emulsion droplets do display additional attractive interactions that have been shown to change the nature of the jamming transition [58][59][60] . In contrast with the purely repulsive case, droplets in attractive emulsions can form bonds and thus a soft gel-like elastic structure which can sustain shear stresses below isostaticity 45,[61][START_REF] Poulin | The European Physical Journal B -Condensed Matter and Complex Systems[END_REF][START_REF] Datta | [END_REF] . However, the microscopic dynamics of the material, i.e. at the scale of the particles, was not explored. As a consequence, it is of particular importance to ask how the response to stress and in particular, the structural and mechanical properties of emulsions are modified by the presence of attractive interactions. Despite their broad applicability, our understanding of the influence of particle-particle interactions on the macroscopic properties of soft matter systems with attractive interactions is currently hindered by a crucial lack of controlled experimental settings.

In this article, we propose a first step towards completing our understanding of the microscopic origin for the macroscopic properties of adhesive emulsions. In particular, we study emulsions in which droplets interact through depletion attraction. First, we find that the static structure of 2D polydisperse emulsions remains unchanged by the introduction of depletion forces. However, the response of 2D monodisperse emulsions under mechanical constraint is impacted by the presence of depletion forces. Indeed, we flow the droplets through a microfluidic constriction in which they have to undergo elasto-plastic remodelling in order to go from a wide channel to a narrow one. In particular, we find that attractive droplets deform more inside the constriction, which we correlate to a shift in the positions of rearrangements. These findings show that depletion attraction forces are sufficient to modify the elasto-plastic response of dense emulsions under a mechanical perturbation. This attraction, even though it is not evidenced in static conditions, impairs rearrangements and in turn promotes an enhanced elastic response under flow.

Materials and methods

Emulsion preparation

Polydisperse emulsions are prepared using a pressure emulsifier (Internal Pressure Type, SPG Technology co.). Silicone oil (viscosity 50mPa.s, Sigma Aldrich) is pushed through a cylindrical Shirasu Porous Glass membrane decorated with 10 µm pores, direcly into a 10mM SDS solution that is maintained under vigorous agitation. The resulting droplets display an average diameter of 42 µm (polydispersity 21%). In order to prepare the emulsion with both SDS concentrations, we use the same droplets and only replace their continuous phase. To do so, the emulsion is washed in a separating funnel in order to replace the continuous phase by solutions of 10 or 45mM SDS in a water/glycerol mixture (60:40 in volume). This enhances the optical quality of the oil/water interface visualization through bright field and confocal microscopy.

For experiments in the constriction and in static packings, we use monodisperse emulsions with an average droplet diameter of 45 µm (polydispersity 3.9%). These emulsions are obtained with a custom made flow-focusing microfluidic set-up (channel size 60 µm × 60 µm, width at the flow-focusing junction 30 µm). We use the same oil and continuous phases for polydisperse and monodisperse emulsions.

Observation and image analysis of 2D static packings

When studying 2D static packings, we consider emulsions that are fluorescently labelled with Nile Red (Sigma Aldrich). To label the emulsion, we incubate it overnight in a SDS buffer (with [SDS] = 10 or 45 mM) saturated in Nile Red allowing the dye to partition between the oil and water phases over time. A 10 µL drop of emulsion is placed between a microscope glass slide (76 x 26 mm, Objekttrager) and a cover slip (24 x 60mm, Knittel Glaser) separated by spacers (50 µm or 30 µm polymethylmethacrilate -PMMA-film, Goodfellow). Droplets are imaged through confocal microscopy (Spinning Disc Xlight V2, Gataca systems) using a 20x objective.

To study the local structure of these static packings, we use a custom Matlab (MathWorks) routine that works as follows. We first threshold the images and perform a watershed tessellation, we then measure the perimeter p and area a of each droplet as well as the area a c of the associated watershed tesselation cell (see Fig 1D). Following Boromand et al. 64 , we study the relation between the deformation of the droplets and their local packing fraction. To do so, we compute their shape factor A = p2 /4πa and determine the local packing fraction φ l = a/a c . Note that we only consider droplets in the center of the packing, i.e. we exclude those that are partially cut by the edge of the image frame. The shape parameter A equals 1 for circular disks and is greater than 1 for all nonspherical particles 65 .

Experimental set-up for emulsion flow

We designed the constriction in a microfluidic channel composed of three main sections (Fig. 1): at the entrance, the channel is 50 µm deep and 200 µm wide over a 5 mm length, then at the constriction the width is reduced from 200 to 50 µm over a length of 200 µm, finally the channel remains 50 µm wide over a final 5 mm length. Note that the droplets are not confined in the vertical direction in this geometry since their diameter is slightly smaller than the height of the channel. The channel is made in polydimethylsiloxane using a negative cast micromachined in a block of PMMA (50 × 50 × 5 mm 3 ) using a desktop CNC Mini-Mill machine (Minitech Machinary Corp., USA). After passivating the channel with casein 0.05 mg/ml (β -casein from bovine milk, Sigma Aldrich) for 20 minutes, the emulsion is flowed in the device using a pressure pump (MFCS-8C Fluigent, P = 30 mbar). After droplets fill the constriction area, the pressure is decreased to stop the emulsion flow, and droplets are left to cream in the supply tube overnight, thus compressing the droplets in the microfluidic device in order to reach high values of packing fraction. After this passive compression phase, the emulsion is flowed again in the channel at a constant pressure. The flow of the droplets at the constriction is imaged in bright field microscopy with a 10x objective at a frequency of 20 frames per second (fps).

Image analysis of the emulsions flowing in the constriction

To analyse the videos of flowing emulsions, we first threshold the images to subsequently determine the center and perimeter of each droplet in the channel using a custom made Matlab rou- tine. When studying droplet deformation, we only consider the droplets located in the constriction region. We define this area along the channel as a window that includes the 200 µm of the constriction itself, plus 50 µm before and after the constriction (Fig. 1). To quantify the deformation of each droplet, we use the approach proposed by Chen et al. 23 . The perimeter of the droplet is interpolated and discretized at 1024 evenly spaced angles θ and the deformation d is calculated as a standard deviation of the radii r(θ ) for each of these angles divided by the mean value of r:

d = r 2 -r 2 r (1)
We also determine the global packing fraction of the emulsion in each video frame. To this end, we calculate the ratio between the sum of all droplets area and the area of the channel within the window of 200 × 200 µm located before the constriction region. Finally, frames are sorted according to the emulsion packing fraction, and the distributions of droplet deformations for each packing fraction are computed.

For rearrangements and flow analysis, the droplets were tracked using a custom Python routine and the FastTrack software (http://www.fasttrack.sh/). All droplets are sorted according to the lane they belong to in the channel ahead of the constriction. In our experiments, they are thus sorted into four lanes. The instantaneous velocity of the droplets was computed as the distance travelled between two consecutive frames acquired at a fixed frame rate

v(t) = √ (x(t+dt)-x(t)) 2 +(y(t+dt)-y(t)) 2 dt
. In our images the spatial resolution yields 300nm per pixel.

Results

Analysis of static packings

We first study 2D static packings of monodisperse and polydisperse emulsions with two distinct depletion interactions. Using silicon oil droplets stabilized with two different concentrations of SDS (10mM and 45mM) allows us to change the depletion forces between the droplets. In our experiment, the continuous aqueous phase is supplemented in glycerol (40 % in volume of glycerol). Note that in addition to allowing for a better imaging of the droplets, it also shifts the critical micellar concentration (CMC) of SDS. However, the CMC is only raised from 8mM (in pure water) to about 9mM in our experimental conditions 66,67 , which ensures that the system is still above the CMC under both SDS concentrations and that the surface tension remains the same when the concentration of SDS is increased from 10 to 45 mM. Above the CMC, depletion attraction forces increase linearly with the concentration of micelles 68 , which itself grows with increasing concentrations of SDS. Given the aggregation numbers of SDS (i.e. the number of SDS molecule per micelle at a given concentration) found in the literature [69][70][71] , we estimate that there is approximately 30 times more micelles at 45mM SDS than at 10mM SDS (see ESI † ). Depletion forces at 45mM SDS are thus expected to be 30 times larger than at 10mM SDS.

To study the impact of depletion forces on static 2D packings, we first quantify the deformation of the droplets as a function of their local packing fraction in these monodisperse and polydisperse emulsions. Thus, we measure the asphericity and local packing fraction of each droplet in several images of 2D monodisperse and polydisperse packings for both 10 and 45mM SDS concentrations (see Materials and Methods). In Fig. 2A, we first plot φ l vs A -1 for monodisperse emulsions for both SDS concentra- tions. Surprisingly, we find that data for both SDS concentrations are collapsing on the same master curve showing that the depletion-induced attractive interactions between the droplets have little effect on the structure of static packings.

1-9 | 3
To confirm the results we obtained for monodisperse emulsions, we also performed experiments on disordered packings (see Fig. 1C-D for an example of polydisperse emulsion). In Fig. 2B, we plot φ l vs A -1 for polydisperse emulsions for both SDS concentrations. Recent studies 64,65 developed a new numerical model to study the structural and mechanical properties of disordered 2D packings of bubbles and emulsions, including at high compressions. In the so-called deformable particle (DP) model, particles deform in response to mechanical constraints to minimize their perimeter while keeping their area fixed. This leads to a model of deformable disks with potential energies that includes an energy term associated to the line tension and a penalization energy term quadratic in the change of area of the droplets, thus associated to their compressibility. Further, the deformable particles interact via a purely repulsive potential energy. Within the framework of this DP model and in our range of deformations, it was predicted that for disordered packings the distance to jamming onset φ lφ c scales with asphericity A -1 as

φ l -φ c = α(A -1) ω (2)
with ω ≈ 0.3.

In the case of the disordered packings, we can compare our experimental data to this theoretical prediction. We find that the data for both SDS concentrations are well-fitted by Eq. ( 2) with scaling exponent ω = 0.33 and critical volume fraction φ c = 0.842, i.e. the scaling form obtained in the DP model. A summary of the fitting procedure and a table of all performed fits for both SDS concentrations are given in ESI † . While the results of the DP model were obtained for disordered systems, we also find that the change in local packing fraction as a function of droplet deformation in monodisperse emulsions is well described by the power law scaling in Eq. (2) with, in this case, a scaling exponent ω ≈ 0.25 and a critical packing fraction φ c = 0.91 which can be explained by a high degree of crystallization.

These scalings are shown on a log-log scale in Fig. 2C, where we show that, in monodisperse emulsions (respectively polydisperse emulsions), data points corresponding to both depletion forces overlap and are captured by the same scaling function with φ c = 0.91 and ω = 0.25 (respectively φ c = 0.842 and ω = 0.33). This indicates that depletion induced attractive interactions do not affect significantly the scaling φ lφ c versus A -1, i.e. changing SDS concentration does not induce any measurable modification in the static packings of droplets. This might seem counterintuitive. Indeed, while it was shown that purely repulsive polydisperse emulsions become fluid-like below random close packing (i.e. only respond elastically above random close packing), experimental studies on the rheology of attractive emulsions showed that attractive emulsions are elastic both below and above random close packing hinting at the fact that loose emulsions can be stabilized by attraction [START_REF] Datta | [END_REF] . Here, we do not observe any significant change in φ c for static packings with depletion interactions. In our SDS stabilized emulsions, the droplets are essentially frictionless and thus are free to roll. In such a system, the slightest amount of compression will lead to a rearrangement of the structure without noticeable deformation of the droplets until their purely repulsive jamming packing fraction is reached. This intuition is confirmed by the measured values of φ c which correspond to random close packing φ RCP ≈ 0.842 for the polydisperse (disordered) emulsions and to hexagonal close packing φ HCP = 0.91 for the monodisperse emulsions.

Despite the fact that static packings cannot be distinguished as a function of depletion forces, we reveal in what follows that significantly distinct behaviors can be evidenced in the context of a dynamic flow.

Emulsion flow in a constriction

In order to study their response under mechanical perturbations, monodisperse emulsions are flowed in microfluidic channels exhibiting a single physical constriction (Fig. 1). In particular, we use monodisperse droplets whose diameter is comparable to the channel height, constraining the system to a 2D monolayer of droplets. We focus our analysis on the area of the constriction in which droplets have to rearrange and deform in order to go from a large channel into a narrower one. The width of the narrow channel is chosen such that it only allows for the passage of one droplet diameter (Fig. 1) in order to maximize the number of rearrangements.

A typical experiment is carried out in two phases. The channel is first filled with the emulsion using a pressure pump. After a waiting time (see Materials and Methods), the pressure is increased again so that this packed emulsion can flow in the channel. We usually require a typical pressure of the order of 30 mbar to establish a continuous flow. For each experiment, we image the droplets upstream, in order to evaluate their packing fraction, as well as inside the constriction to measure their deformation. We choose to quantify the deformation d of each droplet in the channel through the standard deviation of droplet radii as previously done 23 (see Materials and Methods).

Deformation along the channel

We first study the deformation of the droplets inside the channel. To do so, we measure the packing fraction of the emulsion in a window located upstream of the constriction (on the left of the image) and that encompasses 200µm of the channel length (Fig. 3A). We show in Fig. 3 the average deformation d along the channel for both SDS concentrations.

The obtained curves differ for the two SDS concentrations both in the constriction region and in the thinner channel. For both conditions (Fig. 3B-C) the deformation builds up in the constriction to a first maximum average deformation until it is released to a lower value of d at x ≈ 450µm. Then the deformation builds up again to a second maximum and is decreased to a lower deformation. Qualitatively, this behavior can be explained as the signature of a local stress release after a rearrangement. Indeed, Chen et al. 23 showed that in compressed emulsions, T1 events were immediately followed by a local decrease of deformation ). The deformation is averaged in bins that are 25 µm wide along the x-axis. The average deformation peaks inside the area of the constriction for both conditions. The error bars correspond to the standard error of the mean for the distributions of d obtained in each bin. The total number of droplets, combining all packing fractions, is N = 27219 for 10mM SDS and N = 91391 for 45 mM SDS. inside compressed emulsions. Here the localized peaks indicate that droplet rearrangements indeed occur at positions that are set by the topology of the packing in the channel 33 .

The other difference between the two conditions can be observed in the thinner channel region, after the constriction, where droplets enter one by one and release their deformation. In the case of low depletion forces ([SDS]=10mM), droplets relax to a deformation value that is close to the initial one at the entry of the channel ( d outd in ≈ 0.0025). However, with high depletion forces ([SDS]=45mM), droplets relax to a plateau at higher values of deformation than at the entry ( d outd in ≈ 0.01). This impaired relaxation could be a signature of long range effects that could also explain why droplets enter the constriction with a slightly higher value of deformation in the high attraction case.

Deformation as a function of packing fraction

To further confirm these observations, we study the distribution of deformation of all droplets at all positions inside the constriction (taken in a window whose length spans 50µm before and after the constriction -see Materials and Methods). Since the global packing fraction can evolve over the course of one experiment, we separate each experiment into stacks according to their upstream packing fraction. We then pool together the image sequences corresponding to the same packing fraction throughout all performed experiments, for each concentration. Note that we also checked that the deformation in the constriction does not depend on the instantaneous droplet velocity within the investigated range (from 120 to 360 µm/s, see ESI † ).

We compare the distributions of the deformations observed for different packing fractions and for each SDS concentration (Fig. 4). The distributions peak at smaller values of deformation in the low attraction case than in the case of strongly attractive droplets (Fig. 4A). This shift can also be clearly evidenced when plotting the cumulative distributions for each condition at various packing fractions (see Fig. 4B). As expected, for low depletion forces (10mM SDS) we find that the distributions exhibit lower values of deformation in all conditions. When attraction is introduced between droplets, all curves are shifted to higher values of deformation.

In the previous section we showed that depletion alone was not sufficient to induce significant additional deformations in static packings of droplets. The shift observed in these deformation distributions must thus originate from differences in the local topological changes of the emulsions. Hence, we next examine the spatial location of rearrangements in the constriction as a function of SDS concentration.

Rearrangements and velocity distributions in the constriction

We here test the hypothesis that rearrangements are impaired by the attraction between the droplets, which would in turn force the droplets to deform more to overcome the constriction.

Since the size of the channel as well as the diameter of the droplets are fixed, there are always four lines of droplets flowing in the channel, ahead of the constriction, and one line after the constriction (Fig. 5A). In this framework, droplets will exchange neighbors to do the necessary rearrangements in given areas of the channel that are defined by geometry. We can thus estimate the distance ∆x between two rearrangements by calculating the distance between the points where one can accomodate for 4 droplets and 3 droplets of diameter D in the constriction. Knowing the slope of the constriction α this leads to ∆x = D/(2α) = 60µm.

To study these rearrangements, we measure both the velocity and number of neighbors of the droplets along the channel. In the constriction, droplets are stalled transiently until they perform a rearrangement. This effect is evidenced by local minima in their individual velocity profiles as shown in Fig. 5B, that are separated by the expected distance between rearrangements ∆x. We extract the positions of these minima in two specific areas of the channel, corresponding to zones of rearrangements, for each droplet in lines 1 and 4 for both SDS concentrations. We compare both conditions by plotting the cumulative distributions of the minimum velocity location in each zone (Fig. 5C). We observe that the distributions for attractive droplets are shifted by 7µm and 4µm in zones 1 and 2 respectively (measured shift at 50%), indicating that rearrangements are indeed delayed in the channel compared to the low depletion case. Similarly to the deformation, this shift does not seem to depend on the flow speed in the channel (See ESI † ). Furthermore, we analyzed the average number of neighbors per droplet along the constriction (Fig. 5D). Droplets in lines 1 and 4 enter the constriction with 4 neighbors in a hexagonal lat- tice, and exit with 2 neighbors in the smaller channel. We measured the positions along the channel where droplets transition on average from 4 to 3 and from 3 to 2 neighbors respectively. We find that, for the more attractive emulsions, the transitions from 4 to 3 and from 3 to 2 neighbors are delayed by 11µm and 8.5µm respectively.

Discussion

Attractive interactions between particles is expected to affect their packing topology as well as their rheological and mechanical response to local mechanical perturbations. Below the jamming transition, previous work showed that attraction induced by depletion forces tuned significantly the structure of 3D packings and could mechanically stabilize them below the isostatic limit 60 . Above the jamming transition, one expects adhesive forces in packings of deformable spheres to change how droplet deformation and coordination numbers scale with the packing fraction 64,65 . To the best of our knowledge, this issue has been addressed neither in theoretical models nor in experimental systems.

In our experimental study, we provide a first step towards the understanding of the mechanical response of adhesive emulsions by introducing attractive interactions induced by depletion between oil droplets. We first evidence that such attraction forces are too low to induce any measurable effect in 2D static packings of droplets. Indeed, for both attraction forces, we recover the scaling laws predicted by Boromand et al. 64 for purely repulsive packings, with a critical packing fraction φ c ≈ 0.842. However, using monodisperse emulsions, we uncovered distinct changes in their elasto-plastic response when the droplets are flowed through a 2D physical constriction. Note that recent numerical studies dealing with attractive soft sphere packings showed similarly that attractive packings with structures very similar to their repulsive counterpart (e.g. in the fact that the critical packing fraction is not significantly different) could have qualitatively different mechanical properties 72 . In our experiments, the first man-ifestation of attraction is an increase of the average deformation of the droplets in the constriction. The second one is the delay of topological rearrangements inside the constriction as attraction forces are increased. Depletion forces thus appear adequate to redmodulate the elasto-plastic response of emulsions in our system.

Such findings could be relevant for biological tissues in which adhesion controls to a large extent remodelling events that occur on timescales that are beyond those of cytoskeletal activity. In order to isolate the role of adhesion in biological processes, cellular tissues can indeed be mimicked with droplet assemblies connected by specific binders [73][74][75] . Within that framework, emulsions have been shown to exhibit similar mechanical properties and have for this reason been used to measure cellular forces both in vitro 76 and in vivo 77,78 . This reductionist approach could thus shed light on behavioral transitions in developping tissues upon adhesion modulation and will be the focus of future investigations.

Fig. 1

 1 Fig. 1 Experimental set-up and image analysis -(A) The oil in water emulsion is pushed using a pressure pump (P = 30 mbar) through the microfluidic channel that consists of three parts: a 200 µm wide channel, a constriction, and 50 µm wide channel. The depth of the channel is 50 µm over the whole length, and the diameter of the droplets is ≈ 45 µm. (B) Typical image of a monodisperse emulsion flowing in the constriction. In the area of the constriction, the flow of the droplets is imaged in bright field microscopy at 20 fps. The packing fraction of the emulsion is determined within the window of 200x200 µm located before the constriction area. (C) A typical confocal microscopy image of compressed 2D droplets at [SDS]=10mM. (D) Result of the image analysis performed on (C). Droplet contours are shown in red and watershed tessellation cells with the green curves.Based on these measurements, we calculate the local packing fraction φ l as the ratio between the area of the droplet and that of its corresponding watershed tesselation cell, as well as the shape parameter A .

Fig. 2

 2 Fig. 2 Analysis of static 2D packings -Local packing fraction φ l versus A -1 for monodisperse (A) and polydisperse (B) emulsions. In (A), the experimental points (red diamonds for 10mM SDS, blue squares for 45mM SDS) are plotted together with the scaling function in Eq. (2) with exponent ω = 0.25 and φ c = 0.91 (black dot-dashed line). In (B), the experimental points (red open diamonds for 10mM SDS, blue open squares for 45mM SDS) are plotted together with the DP model, i.e. scaling function in Eq. (2) with exponent ω = 1/3 and φ c = 0.842 (black dashed line). (C) Log-log plot of φ lφ c versus A -1 for monodisperse and polydisperse emulsions for both SDS concentrations. We used φ c = 0.91 and φ c = 0.842 for mono and polydisperse packings respectively. The total number of droplets for: (1) polydisperse emulsions are N = 1193 and N = 1735 for 10 and 45mM SDS respectively and (2) monodisperse emulsions are N = 630 and N = 530 for 10 and 45mM SDS respectively.

Fig. 3

 3 Fig. 3 Analyzing the droplet deformation in the constriction -(A) Still snapshot of the image analysis in the channel at a given instant for an attractive emulsion ([SDS]=45mM). The color of the droplets codes for their deformation d calculated for their detected contours displayed on the image. (B-C) Average deformation of the droplets along the x-axis of the channel for different packing fractions in (B) the low attraction case and (C) high attraction case ([SDS]=45mM). The deformation is averaged in bins that are 25 µm wide along the x-axis. The average deformation peaks inside the area of the constriction for both conditions. The error bars correspond to the standard error of the mean for the distributions of d obtained in each bin. The total number of droplets, combining all packing fractions, is N = 27219 for 10mM SDS and N = 91391 for 45 mM SDS.

Fig. 4

 4 Fig. 4 Statistics of deformation under flow -(A) Probability density function of the deformation d calculated in the constriction for different packing fractions in the case of low attraction forces ([SDS] = 10mM, open circles) and high attraction forces ([SDS] = 45mM, stars). (B) Cumulative distributions of the deformation d in the constriction for low attraction forces (open circles) and high depletion forces (stars) for different packing fractions.

Fig. 5

 5 Fig. 5 Rearrangements and velocity distributions in the constriction -(A) Image of droplets in the constriction for 10mM (top, blue) and 45mM (bottom, orange) SDS. The zones where the number of droplets decreases from 3 to 2 and from 2 to 1 are indicated with yellow lines and referred to as zone 1 and 2 respectively. (B) Typical velocity curves of individual droplets in lines 1 and 4 along the channel axis for 10mM (blue) and 45mM SDS (orange) emulsions. The velocity drops to a minimum value (indicated by an arrow) in zones 1 and 2 each time the droplets stall before a rearrangement. For each droplet trajectory we find the location of this minimal velocity to build panel (C). (C) Cumulative distributions of the minimum velocity location for lines 1 and 4. The 45mM SDS droplets slow down farther into the constriction, as evidenced by the shift in distributions at 50% probability, by about 7µm in zone 1 and 4µm in zone 2. (D) Average number of neighbors per droplet along the channel for lines 1 and 4. For more attractive emulsions, the number of neighbors drops from 4 to 3 and then from 3 to 2 respectively 11µm and 8.5µm later in the constriction.

| 1-9

| 1-9

| 1-9

Acknowledgements

Conflicts of interest

There are no conflicts to declare.

The authors thank Adèle Dramé-Maigné and Yannick Rondelez for providing flow focusing devices for the production of the emulsions, as well as Jacques Fattaccioli for allowing us to use his pressure emulsifier, they also thank Georges Debregeas and Zorana Zeravcic for fruitful discussions and acknowledge financial support from Agence Nationale de la Recherche (BOAT, ANR-17-CE30-0001).

Notes and references

1 A. J. Liu


Ismail Alaoui Abdellaoui

Nathan Souviraà-Labastie
email: nathan.souviraa-labastie@a-volute.com

BLENDING THE ATTENTION MECHANISM IN TASNET

Keywords: Audio source separation, Self-attention mechanism, Transformer, TasNet

Audio source separation has seen a rapid progress due to the advances in data-driven approaches. Motivated by the success of the transformer architecture and the self-attention mechanism in many fields, we investigate in this paper how to include these recent advances in an audio source separation algorithm and in particular in the separator of a TasNet architecture. Our main result leads us to think that LSTM layers are required before the attention mechanism. Our other architecture design experiments focus on finding an optimal number of dense nodes within the point-wise neural networks of the transformer block as well as the right weights initialization in relationship with the multiphase gammatone filterbank. Finally, the non-causal version of our proposed algorithm shows promising results on the MUSDB18 test set (SiSEC challenge data).

INTRODUCTION

Data-driven audio source separation greatly benefited from the rapid progress of deep learning techniques. While traditional approaches relied on various priors such as the harmonics of the signal or a linear correlation within the audio signal, recent techniques are data-driven and centered around deep neural networks [START_REF] Rafii | An overview of lead and accompaniment separation in music[END_REF][START_REF] Stöter | The 2018 signal separation evaluation campaign[END_REF]. Within the field of audio source separation, we can distinguish two families of approaches, namely the frequency and time domain approaches. Early approaches to audio source separation were frequency domain approaches and the first data-driven approaches were also frequency based [START_REF] Rafii | An overview of lead and accompaniment separation in music[END_REF].

More specifically, it consists in using a time-frequency representation, e.g. Short-Time Fourier Transform (STFT), of the time domain audio signal as input and predicting the time-frequency representation of the unknown sources before converting back into the time domain. A frequency mask is usually predicted for each source and transcribes the fact that a particular frequency bin corresponds to a given source. Various techniques including ensemble learning as well as deep clustering among others used a time-frequency representation as input, yielding various levels of success [START_REF] Zhang | A deep ensemble learning method for monaural speech separation[END_REF][START_REF] Isik | Single-channel multi-speaker separation using deep clustering[END_REF][START_REF] Kolbaek | Multitalker speech separation with utterance-level permutation invariant training of deep recurrent neural networks[END_REF][START_REF] Chen | Deep attractor network for single-microphone speaker separation[END_REF][START_REF] Luo | Speaker-independent speech separation with deep attractor network[END_REF][START_REF] Jansson | Singing voice separation with deep u-net convolutional networks[END_REF].

On the other hand, time domain approaches quickly became prevalent as they are able to address both the magnitude and the phase of the signal in an end-to-end fashion [START_REF] Fu | End-to-end waveform utterance enhancement for direct evaluation metrics optimization by fully convolutional neural networks[END_REF]. These architectures typically incorporate three subparts: an encoder, a separator, and a decoder. Similarly to the time-frequency approaches, time domain approaches can also be mask-based, and it is the separator that learns to predict relevant masks for each source. Conversely to frequency approaches where the representation is fixed, the representation is learned by the encoder, which converts the raw audio data into a weight matrix. Finally, the decoder applies the masks produced by the separator on the weight matrix and reconstructs each time domain signal for each predicted source. One of the foundational time-domain approaches, TasNet, was first applied to address real time speech separation [START_REF] Luo | Tasnet: time-domain audio separation network for real-time, single-channel speech separation[END_REF], and led to a great state-of-the-art improvement. A noncausal version of TasNet using bidirectional long short-term memory (LSTM) layers in the separator block was also proposed, yielding a better performance than the causal counterpart on the WSJ0-2mix dataset [START_REF] Garofalo | Csr-i (wsj0) complete[END_REF]. Other recent approaches include models using convolutional layers, recurrent layers, and U-Net based architectures [START_REF] Luo | Conv-tasnet: Surpassing ideal timefrequency magnitude masking for speech separation[END_REF][START_REF] Stoller | Wave-u-net: A multi-scale neural network for end-to-end audio source separation[END_REF][START_REF] Défossez | Music source separation in the waveform domain[END_REF][START_REF] Luo | Dual-path rnn: efficient long sequence modeling for time-domain single-channel speech separation[END_REF].

An important subpart of the time domain approaches is the separator block. Since it is the one with the most learnable parameters, it is also most responsible for the performance of the separation. While the separator of the original TasNet is based on LSTMs, there has been an urge to use the multihead self-attention mechanism instead of the LSTM layers since the introduction of the transformer architecure [START_REF] Vaswani | Attention is all you need[END_REF]. This mechanism has proven to be effective for sequence transduction thanks to its parallelizable scaled dot product. Self-attention was successfully applied in sequence transduction for speech recognition [START_REF] Chiu | State-of-the-art speech recognition with sequence-to-sequence models[END_REF], language modeling [START_REF] Devlin | Bert: Pre-training of deep bidirectional transformers for language understanding[END_REF], or recommender systems [START_REF] Zhang | Deep learning based recommender system: A survey and new perspectives[END_REF]. Concerning the usage of the self-attention mechanism for audio separation, a few pioneer approaches have been proposed [START_REF] Li | Sams-net: A sliced attention-based neural network for music source separation[END_REF][START_REF] Li | Atss-net: Target speaker separation via attention-based neural network[END_REF]. Through depthwise separable CNNs and multi-head self-attention, Sam-sNet [START_REF] Li | Sams-net: A sliced attention-based neural network for music source separation[END_REF] obtains results at the level of the state of the art on the MUSDB18 dataset, which is quite promising.

In this paper, we first investigate the incorporation of the self-attention mechanism inside the separator block of TasNet. More specifically, we demonstrate which kind of configuration works well for this particular task, as well as show the impact of the weight initialization with respect to a deterministic encoder. In addition, we propose a noncausal transformer-based separator which equals the original TasNet on MUSDB18 test set with additional benefits.

PROPOSED LTL ARCHITECTURE

In this section, we present the proposed LTL architecture (LSTM-Transformer-LSTM) applied to TasNet's separator and formalize the self-attention mechanism incorporated in it.

Baseline system

The backbone of our model is TasNet [START_REF] Luo | Tasnet: time-domain audio separation network for real-time, single-channel speech separation[END_REF]. Starting from the input signal, a segment is selected to be processed. This segment then passes through two 1D convolutional layers in parallel with different activation functions. One of the layer uses the rectified linear unit (ReLu) while the other one uses the sigmoid function. An element-wise multiplication (Hadamard product) is then performed between the outputs of each layer. This produces the mixture weight matrix w ∈ R L×C , where C is the number of filters used in the 1D convolutional layer and L is the length of the kernel. The decoder block uses a dense layer on the result of the product from the output of the separator and w. This acts as a learned deconvolutional operation, mimicking an inverse operation of the 1D convolutional layers [START_REF] Luo | Tasnet: time-domain audio separation network for real-time, single-channel speech separation[END_REF]. The last step is to concatenate the recovered signals across all segments to obtain each estimated source s i and compute the loss function.

Multi-headed self-attention

Self-attention was introduced by Vaswani et al. [START_REF] Vaswani | Attention is all you need[END_REF] to capture long-range dependencies in sequences through an attention mechanism. Other key aspects of this mechanism are a reduction of the complexity per layer, as well as an increase in the amount of parallel processing during the training. Given a sequence of inputs I ∈ R S×E where S is the sequence length and E is the embedding dimension of each input feature, we first compute the query Q, key K, and value V matrices are computed as follows:

Q = IW q , K = IW k , and V = IW v , (1)
where W q , W k , and W v are learnable weights through a linear function. The attention output, also called the head output, is then computed as follows:

Attention(Q, K, V ) = sof tmax( QK T √ d k )V, (2)
where d k is the dimension of the key vector K ∈ R 1×d k . The operation inside the softmax function is called a scaled dot product, and is implemented to refrain the dot products from having big values [START_REF] Vaswani | Attention is all you need[END_REF]. In addition to this attention mechanism, multiple heads using multiple projections W q , W k , and W v were used to allow the network to attend various areas of the sequence. Given a number of heads H, multi-head attention is computed as:

M ultiHead(Q, K, V ) = Concat(OutHead 1 , ..., OutHead H )W O (3)
where OutHead i is computed through equation 2. The original transformer mentioned is based on an encoder-decoder structure. The encoder part of the transformer contains 2 sub-layers: the multi-head self-attention and the positionwise feed-forward neural network. While the multi-head self-attention part is explained above, the feed-forward neural network uses the output of the multi-head attention as an input. Given the output of the multi-head attention M ultiHeadOut, we feed this output to a 2 layer fullyconnected network with a ReLu activation function after the first layer:

F F N (M ultiHeadOut) = ReLu(M ultiHeadOutW 1 + b 1 )W 2 + b 2 (4)
Residual connections are also included around each of the 2 sub-layers followed by a layer normalization. In this work, we are only interested in the encoder part of the transformer since we focus into the latent space of TasNet's separator.

Proposed LTL architecture

As mentioned previously, we investigate the TasNet separator by incorporating both LSTM layers and encoders of the transformer architecture. More specifically, we try zero or more LSTM layers before and/or after the encoder block of the transformer. This transformer block is made of at least one transformer encoder. Fig. 1 shows the LTL model. For the sake of simplicity, we refer in the legend and later a transformer encoder as simply a transformer layer.

EXPERIMENTAL PROCEDURES

We present in this section a quick description of the music dataset as well as the modalities of the experiments. More specifically, we define the data used, show some hyperparameters used, and formalize the objective function.

Dataset and training procedure

Each input data is a mixture track m i ∈ R l where l is the length in samples of the track, while an output is the set of source tracks S = {s i } C i=1 , where C is the number of target source signals. The entire dataset is a collection of songs divided into training and validation files (80%/20%) also described in [START_REF] Emery | A frugal approach to music source separation[END_REF]. Unless otherwise stated, 10% of all the dataset is used. This choice has been made to shorten experimental duration. Each data sample is made of 3 files: the input mixture song, the target vocal ground truth, and the target accompaniment ground truth. We use audio excerpts of 20 seconds as recommended in [START_REF] Emery | A frugal approach to music source separation[END_REF]. The sampling rate is 44.1 kHz.

Similarly to TasNet, we use the scale-invariant source-tonoise ratio (SI-SNR) as training objective because this loss function proved to be effective for this particular use case [START_REF] Luo | Conv-tasnet: Surpassing ideal timefrequency magnitude masking for speech separation[END_REF][START_REF] Luo | Dual-path rnn: efficient long sequence modeling for time-domain single-channel speech separation[END_REF][START_REF] Emery | A frugal approach to music source separation[END_REF][START_REF] Nachmani | Voice separation with an unknown number of multiple speakers[END_REF][START_REF] Pariente | Asteroid: the PyTorch-based audio source separation toolkit for researchers[END_REF].

Adam algorithm [START_REF] Diederik | Adam: A method for stochastic optimization[END_REF] was used to optimize the SI-SNR, with a learning rate of 1e -4 . Gradient clipping by norm was also used, with a threshold of 5. In addition, a basic learning rate scheduler was implemented with a halving of the learning rate if the validation loss does not improve after three consecutive epochs. An early stopping mechanism is also used after ten consecutive epochs. This is why the curves plotted in Fig. 2 and 3 are of different length. In order to make a fair comparison between the experiments, some hyperparameters were kept constant across all the trainings. These include the hidden size of the LSTMs (500), the number of filters in the 1D convolution of the encoder block (500), the length of these filters (220), as well as the number of the sources [START_REF] Stöter | The 2018 signal separation evaluation campaign[END_REF].

Concerning the used hardware, most experiments used a Nvidia GTX 1080 GPU, while all experiments involving more than 10% of the dataset used a Nvidia Titan X GPU.

MP-GTF

Another recent work in [START_REF] Ditter | A multi-phase gammatone filterbank for speech separation via tasnet[END_REF] proposed a deterministic encoder block instead of a learned one, which is based on the auditory system and which showed to yield an improvement of the scale-invariant source-to-noise ratio (SI-SNR) of 0.7 dB. We decided to incorporate such a deterministic encoder introduced in [START_REF] Ditter | A multi-phase gammatone filterbank for speech separation via tasnet[END_REF] in our experiment. This approach has been adapted to audio source separation, allowing real-time processing and enforcing the non-negative constraint.

Weight initialization

The initialization of the weights in the transformer block at the beginning of the training revealed to be a crucial hyperparameter, because sometimes allowing the gradient descent to converge faster. Here we present the 4 distribution functions used to initialize the weights. The first two functions are the normal and uniform distributions K n and K u of the Kaiming function defined in [START_REF] He | Delving deep into rectifiers: Surpassing human-level performance on imagenet classification[END_REF]. We also use the normal and uniform versions of the Glorot initialization [START_REF] Glorot | Understanding the difficulty of training deep feedforward neural networks[END_REF] X n and X u .

RESULTS AND DISCUSSION

In this section, the findings of a set of experiments are depicted. For the sake of simplicity, we write X trans Y to denote X LSTM layers before the encoder part of the transformer, and Y LSTM layers after it.

The right LTL configuration

A first investigation was based on the right number of LSTM layers before and after the transformer block. Several configurations were tested (i.e., LSTMs before only, after only, both before and after, no LSTMs), and our research showed that using LSTM layers only before the transformer block yielded the best results. More specifically, the best separator architectures were 1 trans 0 and 2 trans 0, performing similarly to our TasNet baseline model with less parameters, or outperfoming the baseline by 0.25 dB SI-SNR, respectively.

Low impact of the fully connected layer

The number of fully connected nodes within the position-wise feed-forward layer of the transformer layer was also subject to analysis. Based on the previous results, both the 1 trans 0 and 2 trans 0 architectures were examined. The number of fully connected nodes tested were 128, 256, 512, 1024, and 2048 for both of these configurations. From a general standing, decreasing the number of fully connected nodes is beneficial to the models because it allows simpler networks without being detrimental in terms of performance. Moreover, 256 dense nodes is the optimal number of dense nodes. Indeed, using 256 nodes allowed the 1 trans 0 model to perform comparably to the baseline model using 11% less parameters, and the 2 trans 0 model to slightly outperform the baseline model. Fig. 2 shows the impact of these experiments with respect to the average SI-SNR of both the loss of the vocals and the accompaniment. An interesting behavior is the suboptimal loss then a sudden drop when using 128 or 2048 dense nodes, while the other numbers of dense nodes make the network easily escape local minima. Therefore, it seems that the network architectures need neither too few nor too many dense nodes.

Impact of the transformer weight initialization

As mentioned in section 3.3, the initial weight distribution plays a critical role in the trainings performed. Therefore, these experiments used different function distributions to initialize the weights of the transformer. Additionally, we also use in some of the experiments the deterministic encoder presented in section 3.2 (i.e. MP-GTF). The results are shown in Fig. 3.

Concerning the trainings involving the 1 trans 0 model, if we analyze the results by pairs of experiments where each pair is made of an experiment using a specific distribution and the other experiment is made of the same distribution jointly used with MP-GTF, we can observe an interesting pattern. The experiment "X u + M P -GT F " gives worse results than using X u only. Similarly, the experiment "K u + M P -GT F " gives worse results than using K u only. In contrast to this observation, using "K n + M P -GT F " increases the performance by approximately 1.5 dB compared to the experiment using the distribution K n only. In addition, the experiment "X n + M P -GT F " leads to an improvement of around 1.25 dB compared to using X n only, even if the training seemed stuck in a local minimum during the first 70 epochs.

To sum up the results observed for 1 trans 0, we can see there is a good synergy between the normal distribution on the transformer and the MP-GTF on the encoder of TasNet when used on the 1 trans 0 architecture since the filterbank consistently improved the performance of the experiments where the normal distribution is used. It would be hard to come up with an interpretation of such result mainly because of the position of the LSTM layer between the encoder (potentially set to MP-GTF) and the transformer block. Furthermore, the opposite effect can be seen when using a uniform distribution. Similar experiments for the 2 trans 0 architecture leads us to choose the Xavier uniform distribution without using the MP-GTF for later experiments.

Bidirectionnal version scores on MUSDB18 test set

Motivated by the superior results of the bidirectional version of TasNet, a LTL architecture using bidirectional LSTMs was also investigated. In particular, a 3 trans 0 model using 2048 dense nodes within the transformer yielded comparable separation performances to our TasNet baseline model (here with bidirectional LSTMs) when using 10% or 100% of the dataset, hence showing good scaling capability. Despite the 2048 dense nodes being unfavorable for the previous causal models, using this number of dense nodes for the bidirectional LTL leads to stable and efficient learning of the separation task. We assume that using three bidirectional LSTM layers induces a transformer block with more paramaters in order to learn the appropriate masks. Moreover, a Xavier uniform was used to instantiate the weights of the transformer block, and no MP-GTF was used.

In terms of SDR on the MUSDB18 test set, the proposed 3 trans 0 achieved an encouraging 7.30 SDR (dB) on the vocals source and 13.46 SDR (dB) on the accompaniment source which is also closed but slightly inferior to the results obtained by our TasNet baseline [START_REF] Emery | A frugal approach to music source separation[END_REF]. However, this proposed 3 trans 0 architecture brings additional benefits compared to our TasNet baseline: firstly, a faster convergence (58% of the time necessary to TasNet to end its training, i.e., 29 days instead of 50) with around 3 times less epochs but with a double epoch duration and secondly a decrease by 53 % of the inference time compared to TasNet, but the number of learnable parameters increased (+43%, i.e., 33M ). One can argue that the faster convergence is due to a higher number of parameters but similar behaviors are observed during experiments of Section 4.2 without the bias of different number of parameters.

CONCLUSION

This paper explores the use of the self-attention mechanism for the task of music source separation within the state of the art framework of TasNet. In particular, we described a series of experiments that lead us to our proposed new architecture scheme for the separator module of the TasNet framework/architecture: 1 to 3 LSTM layers followed by layers composed of the encoder part of a transformer. Experiments on the combination of a non-learned TasNet encoder (using Gammatone) and different transformer weight initialization are also described with intriguing behaviours. Finally, a bidirectional version of the proposed model is evaluated on the MUSDB18 test set (SiSEC challenge data) and reaches comparable separation performances as the state of the art while significantly reducing training, convergence and in-use time. Future experiments will focus on augmenting the length of the audio excerpts to conclude on the capacity of our architecture to better capture long dependencies.

Fig. 1 .

 1 Fig. 1. Schema of the LTL architecture.

Fig. 2 .Fig. 3 .

 23 Fig. 2. The impact of the number of fully connected nodes in the transformer part of the 2 trans 0 architecture.


Yannick Petit
email: yannick.petit@u-bordeaux.fr

Alexandra Pena

Patricia Segonds

Jérôme Debray

Simon Joly

Andrius Zukauskas

Fredrik Laurell

Valdas Pasiskevicius

Carlota Canalias

Benoit Boulanger
email: benoit.boulanger@neel.cnrs.fr

Alexandra Peña

Benoît Boulanger

Development and Experimental Demonstration of Negative First-Order Quasi-Phase Matching in a periodically-poled Rb-doped KTiOPO4 crystal

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Phase-matching is the privileged way to optimize the performance of nonlinear optical frequency conversion. Usually, it can be obtained by using the birefringence of anisotropic crystals, which corresponds to birefringence phase-matching (BPM) [1,2]. However, BPM conditions are strongly limited by the refractive indices of the material and they usually forbid the exploitation of the highest coefficient of the second-order electric susceptibility tensor. To overcome this limitation, quasi-phase matching (QPM) is a powerful alternative. It is based on a periodical modulation of the sign of the second-order optical nonlinear coefficient of anisotropic or isotropic non-centrosymmetric crystals [3,4]. We propose here the development and the experimental validation of a new scheme of QPM that allows us to enlarge the spectral range of phase-matched generation and to get access to giant spectral acceptances. For that purpose, we studied Second-Harmonic Generation (SHG) in the ferroelectric periodically-poled Rubidium-doped KTiOPO4 crystal (PPRKTP) [5].

The key parameter of QPM is the periodicity Λ of the spatial modulation of the second-order nonlinearity along the direction of propagation of three interacting waves such as their wavelengths fulfil: 1/λ ! -1/λ " -1/λ # = 0. Actually, the goal is to compensate the spatial phase-mismatch Δk ****⃗ = k ! ****⃗ -k " ****⃗ -k # ****⃗ , with ,k $ ****⃗ , =

#% & ! n $
where n $ is the refractive index at the wavelength λ $ (p = 1, 2, 3). Then the QPM relation writes [3]:

Δk ****⃗ = mk ' ****⃗ (1)
where

k ' ****⃗ , with ,k ' ****⃗ , = #%
' , is the grating vector oriented along the propagation direction, and m is the order of the Fourier series. Usually, m is chosen as a positive integer (m ≥ 1), as it is the case for instance while targeting the largest nonlinear coefficient χ !! (#) . It is shown in Fig. 1 (left) in the optimal case, i.e. for the lowest QPM positive order m = 1. This is the classical scheme of QPM that we called QPM-A. However, it is also possible to consider propagation directions where m ≤ -1. This configuration has been previously theoretically considered, but without any phenomenological interpretation or experimental verification [6]. We called here this new scheme QPM-B. As for the case of positive orders, the lowest QPM-B order (namely m = -1) shown in Fig. 1 (right) leads to the best efficiency compared to higher QPM-B orders (namely m < -1), which requires the shortest needed QPM periodicity. New efficient directions for parametric frequency conversion are expected while using QPM-B since the corresponding phase-matching relation is different from that of QPM-A. Fig. 1. Schemes of collinear first-order Quasi-Phase-Matching (QPM). k ' """"⃗ , where $k ' """"⃗ $ =

There is also no difference between the two cases (m = 1, m = -1) in the initial phases φ1, φ2 and φ3 of the three interacting waves when Eq. ( 1) is fulfilled, i.e. for a perfect QPM. Actually, in this case the phases are locked in such a way that the relation φ3 -φ1 -φ2 = ± π/2, as in the case of BPM [1]. . This is similar to another four-particle interaction corresponding to four-wave mixing (FWM). Another instructive analogy deals with light diffraction by a periodic acoustic wave under the Bragg condition in a thick acoustic grating configuration. In such a quasi-elastic diffraction condition, it is well established that both positive and negative first-order Bragg interactions exist and are associated either to the annihilation or creation of an acoustic phonon, respectively [8].

From the thermodynamics point of view, efficient QPM can occur whatever the sign of the evolution of the phase mismatch between the nonlinear polarization and the generated optical field. As a consequence, although we generally deal with optical materials with normal spectral dispersion (namely ∂n ∂λ ⁄ < 0), the generalized QPM approach with positive or negative orders would also apply to materials with an anomalous spectral dispersion (namely ∂n ∂λ ⁄ > 0). Because QPM most generally deals with materials with normal dispersion, situations where Δk ****⃗ = m k ' ****⃗ ≠ 0 *⃗ with m > 0 mostly occur, while Δk ****⃗ = 0 *⃗ is reached for BPM situations. In this sense, the existence of BPM solutions is expected to be a necessary condition to access QPM-B solutions. As a consequence, the QPM-B condition can be fulfilled only for types that can exhibit BPM, so that only Types I, II and III are likely to be addressed by QPM-B schemes in materials with normal dispersion of their refractive indices.

We used an experiment of Angular-Quasi-Phase-Matching (AQPM) for the validation of QPM-B. AQPM is a generalization of QPM since the interacting waves are allowed to propagate at any angle with respect to the grating vector k ' ****⃗ , giving then access to new spectral ranges and giant spectral acceptances [9][10][11]. Then, it can be applied to QPM-A as well as QPM-B, leading to the so-called AQPM-A and AQPM-B used in the following. In these cases, the grating periodicity Λ has to be replaced by an effective periodicity defined along the collinear propagation direction as follows:

Λ /00 (θ, ϕ) = Λ |sin(θ) cos (ϕ)| ⁄ (2)
Equation ( 2) gives the grating period as seen in the direction of propagation that is considered, with the angles of spherical coordinates θ and ϕ in the dielectric frame (X,Y,Z) [5]. Λ /00 (θ, ϕ) ranges from a minimal value corresponding to a propagation in the x-axis that is orthogonal to the domains plane, i.e.

Leff(θ = 90°, ϕ = 0°) = Λ, to a maximal value when the propagation occurs in the YZ plane that is parallel to the domains plane, i.e. Leff (θ, ϕ =90°) ® ¥. The last situation corresponds to BPM, so that the QPM solutions converge for continuity reasons to the BPM solutions, while approaching the YZ plane. It is important to notice that the wave vectors of the three interacting waves remain collinear during an AQPM-B, as it is the case for AQPM-A [9].

The measurements were carried out in the PPRKTP biaxial crystal with a grating period of Λ = 38.52 µm, as it has been done in our previous study to determine the conditions of AQPM-A for SHG (1/λ 1 + 1/λ 1 → 1/λ #1 ). The corresponding curves are plotted in Fig.

2 [10]. We considered the same fundamental wavelength, i.e. λω = 2.15 µm for the present study of AQPM-B. At this wavelength, using Eqs. ( 1) and ( 2) and the Sellmeier equations of reference [12], AQPM-B is expected to be possible for both type I and type II (the

𝑘 " -𝑘 % = 𝑘 ' + 𝑘 ) 𝑘 " + 𝑘 % = 𝑘 ' + 𝑘 ) 3 k !! " 1 k !" 2 k !! " k L !!" 3 k !! " 1 k !" 2 k !! " k L !!" 𝑘 " -𝑘 ' -𝑘 ) ± 𝑘 % = 0 3 1 2 k k k > + !! " !" !! " 3 1 2 k k k < + !! " !" !! " QPM-A QPM-B
latter being equivalent to type III in the case of SHG), as shown in Fig. 2. We experimentally focused on type II because the conversion efficiency of Type I AQPM-B is more than one order of magnitude smaller than that of type II, as it is also the case for BPM. Moreover, the tuning curve of Type II AQPM-B does not cross any other SHG BPM or AQPM curve as shown in Fig. 2, leading to an unambiguous identification of these newly proposed solutions.

Fig. 2. Second-Harmonic Generation for a fundamental wavelength of λω = 2.15 µm in a PPRKTP crystal with a grating period Λ = 38.52 µm and shaped as a sphere: calculated Types I and II BPM angles (black dashed lines), calculated Types I, II, IV and V AQPM-A angles (blue dashed lines), Types I and II AQPM-B angles (red continuous lines for the calculations and red stars for the present measurements). Black and blue squares are the experimental results obtained in a previous work [11]. θ and ϕ are the angles of spherical coordinates in the dielectric frame (X,Y,Z): ϕ = 0° corresponds to the XZ plane, ϕ = 90° to the YZ plane, and θ = 90° to the XY plane.

The PPRKTP crystal with a grating period of Λ = 38.52 µm was cut as a sphere in order to be able to access any direction of propagation, with the goal of being able to determine all AQPM-B directions. The sphere, with a diameter of 4.76 mm, was oriented by the Laue method, mounted on a goniometric head, and placed at the center of a motorized Kappa circle allowing the sphere to rotate on itself, as shown in Fig. 3 [10]. A homemade interfacing program enables the correspondence between the angles of the spherical coordinates (q, ϕ) and the Kappa circle angles (k, Fk, Wk) depicted in Fig. 3. In order to get the whole type II AQPM-B curve, the sphere was successively oriented along the x-axis and the z-axis for orientations verifying 42° ≤ q ≤ 65° and 65° ≤ q ≤ 90°, respectively. A pump beam with a fundamental wavelength of 2.15 µm was emitted by a nanosecond BBO Optical Parametric Oscillator. It was carefully focused at the center of the PPRKTP sphere, which enabled a collimated propagation of the incident beam along the diameter of the crystal sphere [9]. The polarization of the incident beam was controlled with a half-wave plate. The SHG phase-matching angles corresponded to the orientations leading to a maximal value of the associated conversion efficiency that was detected using an amplified Si Hamamatsu C2719 photodiode. A NIRquest 512 Ocean Optics spectrometer, with an accuracy of ±3 nm, was used to control and confirm the same phase-matching wavelength for each of the 12 experimental points shown in Fig. 2. The agreement between the measured and calculated phasematching curves is excellent, hence providing a clear demonstration of type II AQPM-B in PPRKTP, and thus the proof-of-concept for negative order QPM. Based on the existence of these new AQPM schemes, we investigated numerically the associated benefit in terms of spectral tunability and associated spectral acceptances. Calculations of SHG solutions were performed in the XY plane, as this plane is associated with high effective nonlinear coefficients for distinct AQPM types, such as types II and V.

As shown in Fig. 4(a), Type II BPM led to solutions for λω spanning from 0.994 µm to 1.079 µm, and from 3.180 µm to 3.308 µm. Type II AQPM-A showed different solutions in distinct spectral ranges, with λω spanning from 0.952 µm to 0.994 µm, and from 3.308 µm to 4.377 µm. Type II AQPM-B completed the spectral range of the SHG interaction in the XY plane, with λω spanning from 0.994 µm to 3.308 µm. Type V AQPM-A covers the spectral range from 2.057 µm to 3.099 µm. Non-critical spectral conditions are evidenced at λω = 1.567 µm for type II AQPM-B and at λω = 2.503 µm for Type V AQPM-A, when f = 25.53° and 23.22°, respectively. The associated spectral acceptances were determined from the full width at 0.405 of the sinc 2 interference function. It comes: L.δλ = 76 nm.cm at λω = 1.567 µm for Type II AQPM-B, and 161 nm.cm at λω = 2.503 µm for Type V AQPM-A, as shown in Fig. 4(b). Thus, it appears that the proposed AQPM-B scheme opens new potentialities, especially as it provides giant acceptances in the telecom spectral range for PPRKTP.

To ensure that the calculated situations considered in Fig. 4 are practically relevant, it is pertinent to consider the figure of merit (FOM), which is classically defined as:

FOM (θ, ϕ) = 2 "## $ (3,5) 6 $% (3,5)•6 &,% (3,5)•6 $,% (3,5) (3)
The quantity d /00 (θ, ϕ) is the effective coefficient, and n 8,1 (θ, ϕ) with i = 1 or 2, and n #1 (θ, ϕ) are respectively the refractive indices of the fundamental and second-harmonic waves in the direction of propagation that is considered. In Table 1 are given the FOM for the two non-critical situations depicted in Fig. 4 where the wavelength acceptances are the broadest, i.e. Type V AQPM-A at (f = 23.22°, λω = 2.502 µm) and Type II AQPM-B at (f = 25.53°, λω = 1.567 µm). The values of the refractive indices and nonlinear coefficients d15 and d24 that are used come from reference [12]. The walk-off angle that is negligible in the XY plane of PPRKTP has been omitted. Note that the effective nonlinear coefficient is defined by the same expression for both the positive and negative first-order QPM conditions given a type and the associated polarization configuration. As a consequence, the effective nonlinear coefficient contains a Fourier factor of 2/π since we consider the first-order modulation of the sign of the nonlinearity, as in the case of AQPM-A. However, the effective coefficient shows different values, as it is estimated along distinct propagation directions, respectively corresponding to the AQPM-A and AQPM-B solutions of Eq. ( 1). Even if the FOM of Type V AQPM-A in the XY plane is 14 times higher than that of Type II AQPM-B, the FOM value for Type II AQPM-B is not negligible and it shows real practical potential as this interaction occurs in the telecom C-band close to 1.55 µm.

In conclusion, we performed a comprehensive interpretation and experimental validation of a generalized description of the QPM interaction by incorporating negative QPM orders, in addition to the well-known positive ones. We have labelled QPM with negative orders as the QPM-B cases, contrarily to the QPM-A cases with positive QPM orders. We have demonstrated the existence of such QPM-B situations in the case of type II SHG with a fundamental beam at 2.15 µm in a PPRKTP crystal exhibiting a period Λ = 38.52 µm perpendicularly to the ferroelectric domains. We have highlighted that AQPM-B gives access to extended spectral ranges and non-critical spectral tolerances. A giant spectral tolerance of 76 nm.cm at λω = 1.567 µm was determined for type II AQPM-B in the XY plane, which interestingly lies in the telecom spectral range. The associated figure of merit was estimated to be 0.36. The generalized QPM description enlarges the palette of QPM interactions and widens the interaction of nonlinear interactions in periodicallypoled media. Knowing the diversity of BPM cones of solutions and the even larger one for the AQPM-A case [2,13], the newly demonstrated negative QPM solutions keep on enlarging the possibilities of favourable configurations for cascaded nonlinear parametric frequency conversion interactions, while considering simultaneous phase-matching and/or positive/negative quasiphase matching conditions. Finally, this approach naturally extends to other periodically poled uniaxial crystals such as the periodicallypoled LiNbO3 (PPLN).

Funding. Centre National de la Recherche Scientifique (CNRS) International Project of Scientific Cooperation France -Sweden (PICS no. 7739).

Disclosures. The authors declare no conflicts of interest.

From

  photon at k ! ****⃗ . For the lowest positive order (m = 1) three particles with k " ****⃗ , k # ****⃗ and k ' ****⃗ are consumed to produce a new particle associated to k ! ****⃗ = k " ****⃗ + k # ****⃗ + k ' ****⃗ . This is somehow similar what occurs while considering the four-particle interaction that corresponds to Third Harmonic Generation (THG) governed by the third-order optical nonlinearity. For the lowest negative order (m = -1) two particles with k " ****⃗ and k # ****⃗ are consumed to produce a new pair of particles, one being the new photon associated to k ! ****⃗ and the

Fig. 3 .

 3 Fig. 3. Schematic diagram of the experiment of Second-Harmonic Generation (𝜔 + ω → 2ω) used for the determination of the angular distribution of QPM-B. (X,Y,Z) is the dielectric frame of the PPRKTP sphere with a grating period Λ. The sphere is placed at the centre of a Kappa circle whose rotation angles are (k, Fk, Wk).

Fig. 4 .

 4 Fig. 4. (a) Spherical AQPM angle ϕ versus fundamental wavelength λω for SHG solutions in the XY dielectric plane of the PPRKTP crystal with the poling period Λ = 38.52 µm: type II BPM (dashed pink line); type II AQPM-A and type II AQPM-B (yellow and blue lines, respectively); and type V AQPM-A (red line). (b) Spectral dependence of the sinc 2 (∆kL/2) interference function for a crystal length L = 1 cm under the undepleted pump approximation (UPA) and at the non-critical wavelengths for type II AQPM-B SHG (λω = 1.567 µm, fixed solution angle at f = 25.53°, blue line), and for type V AQPM-A (λω = 2.503 µm, fixed solution angle at f = 23.22°, red line).

Table 1 .

 1 Figure of merit, FOM (θ= 90°, 𝛟), values in the XY plane for type II AQPM-B and type V AQPM-A at non-critical spectral conditions.