F Faugeras

J Coulouvrat

G Laurent

C Marrelec

W Contino-Pépin

C Urbach

Tribet

†, ‡, C. Thimonier, †, ¶, ‡, S R Ramesh

Desgranges

V Faugeras

F Coulouvrat

J Laurent

G Marrelec

C Contino-Pépin

W Urbach

# C Tribet

N Taulier

Acoustic droplet vaporization of perfluorohexane emulsions is induced by heterogeneous nucleation at an ultrasonic frequency of 1.1 MHz

Acoustic droplet vaporization of perfluorohexane emulsions is induced by heterogeneous nucleation at an ultrasonic frequency of 1.1 MHz

Introduction

Liquid PFCs possess a phase change ability: subjecting PFC droplets to sufficient acoustic pressure triggers a liquid-to-vapor transition, called Acoustic Droplet Vaporization (ADV). [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] The therapeutic potentials of ADV have been exploited for various applications including embolotherapy, 2 phase aberration correction, [START_REF] Fabiilli | The role of inertial cavitation in acoustic droplet vaporization[END_REF][START_REF] Aliabouzar | Effects of droplet size and perfluorocarbon boiling point on the frequency dependence of acoustic vaporization threshold[END_REF] and high-intensity focused ultrasound (HIFU) therapy. [START_REF] Kripfgans | In vivo droplet vaporization for occlusion therapy and phase aberration correction[END_REF][START_REF] Sokka | MRI-guided gas bubble enhanced ultrasound heating in in vivo rabbit thigh[END_REF][START_REF] Zhang | An in vitro Study of a Phase-Shift Nanoemulsion: A Potential Nucleation Agent for Bubble-Enhanced HIFU Tumor Ablation[END_REF][START_REF] Zhang | Acoustic Droplet Vaporization for Enhancement of Thermal Ablation by High Intensity Focused Ultrasound[END_REF][START_REF] Chen | Targeted drug delivery with focused ultrasound-induced blood-brain barrier opening using acoustically-activated nanodroplets[END_REF] However , the physics underlying ADV is complex and the understanding of the ADV mechanism is currently an active research domain. It is admitted that ADV is the result of the appearance of a vapor nucleus, during the rarefactional phase of an acoustic pulse, and of its subsequent growth. Literature investigations have reported that the threshold is a function of the driving ultrasound frequency and pressure, with decreasing pressures needed at higher frequencies. [START_REF] Kripfgans | In vivo droplet vaporization for occlusion therapy and phase aberration correction[END_REF][START_REF] Giesecke | Ultrasound-mediated cavitation thresholds of liquid perfluorocarbon droplets in vitro[END_REF][START_REF] Schad | In vitro characterization of perfluorocarbon droplets for focused ultrasound therapy[END_REF][START_REF] Williams | Characterization of Submicron Phase-change Perfluorocarbon Droplets for Extravascular Ultrasound Imaging of Cancer[END_REF] In addition, droplet size also plays a role, with lower pressures being required for larger droplets. [START_REF] Kripfgans | In vivo droplet vaporization for occlusion therapy and phase aberration correction[END_REF][START_REF] Williams | Characterization of Submicron Phase-change Perfluorocarbon Droplets for Extravascular Ultrasound Imaging of Cancer[END_REF][START_REF] Martz | Precision Manufacture of Phase-Change Perfluorocarbon Droplets Using Microfluidics[END_REF][START_REF] Sheeran | Design of ultrasonically-activatable nanoparticles using low boiling point perfluorocarbons[END_REF] There is no concrete model capable of predicting all these behaviors. However, a model developed by Shpak et al. [START_REF] Shpak | Ultrafast dynamics of the acoustic vaporization of phasechange microdroplets[END_REF] is able to predict the location of the initial vaporization nucleus appearance when high ultrasonic harmonics are present with high enough intensity. Indeed, the focusing of these harmonics at a region inside the droplet volume, induces a high negative pressure that in turn creates the nucleus whose growth leads to ADV. This effect is all the more efficient when the droplet is large as it requires lower ultrasonic harmonics. Similar observations reported in related studies [START_REF] Kripfgans | Acoustic droplet vaporization for therapeutic and diagnostic applications[END_REF][START_REF] Sheeran | Phase change events of volatile liquid perfluorocarbon contrast agents produce unique acoustic signatures[END_REF][START_REF] Lacour | Ultimate fate of a dynamical bubble/droplet system following acoustic vaporization[END_REF] have confirmed the validity of this model. However, this effect is not systematically observed for micron-sized droplets as Shpak et al. [START_REF] Shpak | Ultrafast dynamics of the acoustic vaporization of phasechange microdroplets[END_REF] determined a nucleation map constructed from ultra-high speed imaging, displaying bubbles with nucleation spots located at random positions throughout the droplet, when its radius was R = 5.6 ± 2.1 µm.

According to this model, superharmonic focusing effect is not expected at low frequency (f < 1.5 MHz) or for small droplet (R < 3 µm). [START_REF] Shpak | Ultrafast dynamics of the acoustic vaporization of phasechange microdroplets[END_REF] Consequently, outside the range of application of the superharmonic focusing effect, another nucleation mechanism is responsible for ADV, in particular for nanometer-sized droplets which are often used in ADV studies involving medical applications. The search for a mechanism is of course hindered by the fact that individual nanodroplets cannot be detected optically, thus preventing evidences from direct observations. However, two mechanisms can be expected: homogeneous and heterogeneous nucleations. The appearance of a nucleus of a critical size is responsible for ADV in both these mechanisms. In the first case (i.e., homogeneous nucleation), the nucleus appears inside the droplet volume, while in the second case (i.e., heterogeneous nucleation) it appears at the droplet surface. Current nucleation models show that heterogeneous nucleation can require less energy that homogeneous nucleation. [START_REF] Fletcher | Size Effect in Heterogeneous Nucleation[END_REF] Thus, in some conditions a nucleus is most likely to appear at the droplet surface than in its volume. Several investigations has already showed the importance of the properties of surface properties in vaporization, suggesting that heterogeneous nucleation is taken place in their system. [START_REF] Chattaraj | Phase behavior of mixed lipid monolayers on perfluorocarbon nanoemulsions and its effect on acoustic contrast[END_REF][START_REF] Shakya | Vaporizable endoskeletal droplets via tunable interfacial melting transitions[END_REF] In order to gain a better knowledge of the mechanism that initiates ADV in the absence of superharmonic focusing, we have investigated three types of PFC droplets which allow to vary the area and curvature of the surfaces in contact with perfluorocarbon:

1. Plain PFH droplets surrounded by a monolayer of a homemade fluorinated surfactant, 2. Water multi-core PFH droplets, similar to (1) but containing in their core nanodroplets of water surrounded by a shell made of a commercial surfactant, ensuring their dispersion in the PFH, 3. Water single-core PFH droplet, similar to (2) but with a single water microdroplet instead of several water nanodroplets.

For these three types of emulsions, we have determined, at an ultrasonic frequency of 1.1 MHz, the variation of ADV pressure threshold with the droplet number and size, as well as with the volume fraction of water trapped in the PFH.

Materials & methods

Materials

Unless otherwise stated, all chemicals were used without further purification and all solutions were filtered using Acrodisc Syringe Filters (from PALL) with a pore size of 0.2 µm.

All aqueous solutions were made using Milli-Q IQ 7000 Type-1 water Purification System.

Synthesis of F-TAC surfactant

F-TAC surfactants are amphiphilic molecules composed of two structural parts (Fig. 1A).

One part is the polar head, made of a water-soluble oligomer of Tris(hydroxymethyl) aminomethane (Tris). The second is the fluorinated tail. F-TAC surfactants are readily available at the gram scale in one step by radical telomerization of the monomer tris(hydroxymethyl)acrylamidomethane (THAM) in presence of a telogen reagent (perfluoroalkanethiol) and a radical initiator (α,α'-azobisisobutyronile or AIBN) according to the reaction shown in Fig. 1B. [START_REF] Desgranges | Micron-sized PFOB liquid core droplets stabilized with tailored-made perfluorinated surfactants as a new class of endovascular sono-sensitizers for focused ultrasound thermotherapy[END_REF] The length of the polar head can be tuned according to the concentration of starting reagents and conditions carried out for their synthesis. Within a perfluorocarbon emulsion, the polar head of the F-TAC surfactant is located at the surface of perfluorocarbon droplets, in contact with water, while the fluorinated tail acts as a fluorophilic anchorage ensuring the stabilization of PFH droplets. In order to ensure an optimal droplet stabilization and water solubility, the F-TAC used in the current study is made of a fluorinated tail with eight fluorinated carbons and endowed with thirteen Tris units, and will be referred as "F 8 -TAC 13 " surfactant.

Interfacial tension

The interfacial tensions were measured with an error of 2 mN/m using either a Tracker tensiometer (Teclis, France) or a SDT tensiometer (Krüss) at 20 ℃ (for samples with glycerol).

The first instrument analyzes the shape of either a pendant PFH droplet in water, or a bubble in water or liquid PFH, while the second analyses the shape of a spinning droplet or bubble. All droplet or bubble shapes were fitted by the Young-Laplace equation in order to derive the surface or interfacial tension. Measurements were made after the value of interfacial tension reached an equilibrium value. Some measurements were performed in the presence of 0.1%wt F 8 TAC 13 (solubilized in water or glycerol) or 5% Krytox 157 FSL (solubilized in PFH).

Droplet Generation and Characterization

Microfluidic chip fabrication

We manufactured two types of chips. The first one, suitable for the production of either plain PFH droplets or water multi-core PFH droplets, uses a single classic flow-focusing chip (top and middle figures in Fig. 2). In these chips, all channels are rectangular (50 µm deep and 100 µm wide), while the nozzle size is 10 µm and 40 µm depending of the droplet radius (see Table 1).

The second type of chips was used for the production of water single-core PFH droplets and is made of two flow-focusing junctions in series inside one chip as shown in bottom figure of Fig. 2. All channels are rectangular (50 µm deep and 100 µm wide), the size of the first nozzle is 35 µm, while for the second nozzle has a size of 70 µm.

Both chip types were first designed on AutoCAD, then printed using a two-photon polymerization printer, a Nanoscribe GT Photonic Professional device, with a negative-tone photoresist IP-S (Nanoscribe GmbH, Germany) and 25× objective, directly on silicon substrates after a nitrogen plasma cleaning preparation (to increase the resin adhesion on the substrate). To reduce printing time, a shell writing strategy was applied. It consists of fabricating a dense shell delimiting the structure, the inner being only partly polymerized in the form of a scaffold. After developing in propylene glycol methyl ether acetate (for 30 min) and isopropanol (for 5 min), a batch polymerization is performed with UV-exposure. Soft lithography was then used to replicate the channels. PDMS and its cross-linking agent were mixed at a ratio 10:1 (wt/wt), degassed and poured onto the mold. The polymer was cured for 2 hours at 70 ℃. PDMS chips were then peeled off carefully from the mold, and 350 µm holes were punched into the chip with biopsy punchers to create the inlets and outlets. Once the design was etched onto the PDMS, the patterned surface of the chip and the surface of a glass slide were activated in a plasma cleaner filled with air, for 1 min at a power of 18 W.

This plasma activation ensures the secure bonding of the chip on the glass slide. The walls of the single flow-focusing junction chip were treated hydrophilic. To do this, the chip was plasma activated for 1 minute and water was passed through the channels.

Production of plain PFH droplets

We used two distinct methods to generate perfluorohexane droplets depending on their size.

In both methods, the mass ratio of F 8 -TAC 13 to continuous aqueous phase was 0.1%.

Plain PFH droplets with a nanometer-sized radius were produced using a high-pressure homogenizer (model LV1 from Microfluidics). For this technique, a coarse emulsion was first prepared by vortexing. The coarse emulsion was then used with the high-pressure homogenizer, where the emulsion went at a high velocity (thanks to a high pressure pump) through two capillaries which met inside a Y shape chamber (model F12Y, internal dimension of 75 µm). The high shear forces occurring in this chamber between the two emulsion flows induce a reduction in droplet size. The passage through the chamber was repeated 8 times at a pressure of 138 MPa. The resulting emulsion was centrifuged at 6000 rpm (i.e., at a relative centrifugal force of 2000 g) for 30 s using a bench mini-centrifuge (Mini Star from VWR) to force any droplets with diameter larger than 1 µm to sediment and only the supernatant was kept.

Plain PFH droplets with a micron-sized radius were generated by microfluidics using a conventional single flow focusing chip where the flow rates of the dispersed phase (channel filled with perfluorohexane, white arrow in the top figure of Fig. 2) and of the continuous phase (channel filled with the aqueous solution solubilizing F 8 -TAC 13 (orange arrows in the top figure of Fig. 2) are given in Table 1.

Production of water single-core PFH droplets

These droplets were produced using a two-nozzle device (see Fig. 2). Following the procedure of Bodin-Thomazo et al. [START_REF] Bodin-Thomazo | Marker patterning: a spatially resolved method for tuning the wettability of PDMS[END_REF] , the channel making the junction between the PFH and water channels (i.e., channel C between channels B and D in Fig. 2), was selectively patterned using a black marker to avoid the wetting of the walls by PFH. We used a flat PDMS layer as a substrate for the chip instead of a glass slide. This was followed by the surface treatment procedure previously explained (air plasma for 1 min). After the treatment, the chip and the flat PDMS were aligned to fit the marker pattern. In this chip, the walls of the unpatterned channels remain hydrophilic, while the walls of the patterned channel are kept hydrophobic.

The chip was first flushed with methanol, then flushed with water for 20 minutes before starting the production of droplets. Then, the different fluid phases (PFH and aqueous solutions) were injected into the microfluidic device by applying a pressure of 20 mbar to the headspace of their respective inlets. The flow rates were adjusted to encapsulate exactly one water droplet in each double emulsion within these ranges: 0.8-2.3 µL/min for the dispersed aqueous phase (orange arrow in channel A on the bottom figure of Fig. 2), 1-5.5 µL/min for the intermediate PFH phase (white arrows in channels B), and 2-8 µL/min for the continuous aqueous phase (orange arrows in channels D).

Production of water multi-core PFH droplets

A coarse emulsion of water in PFH was first produced inside a vial from a mixture of normal saline (0.9% wt NaCl), PFH, and Krytox (5% wt). The ratio of normal saline to PFH was adjusted depending on the desired water volume fraction inside a multi-core droplet. The coarse emulsion was obtained by sonication using a Branson digital Sonifier (model 450) with a double step microtip at a frequency of 20 kHz, 30% of maximum amplitude for a pulse of 5 seconds, and 10 seconds pause, for 2 minutes. During sonication, the vial containing the fluids was immersed in an ice-water bath to avoid the vaporization of PFH due to the heating generated during the process. An emulsion made of normal saline nanodroplets surrounded by Krytox 157 FSL and dispersed in the PFH was thus produced.

In the next step the multi-core droplets were generated in the same flow-focusing chip as the one used for a plain droplet generation, with the above primary nanoemulsions as the dispersed phase and water as the continuous phase. The resulting multi-core droplet therefore consists of multiple aqueous nanodroplets inside a larger PFH droplet surrounded by water. We define ϕ w as the volume fraction of water composing a droplet and R w the mean radius of the water droplets encapsulated inside the PFH droplets.

The flow rates used to produce these emulsions are in the same ranges as the ones used for plain PFH droplets.

Size and polydispersity measurements

For micron-sized droplets, the size and polydispersity were estimated from a series of videos taken with an ultra-fast camera (Model SC1, Edgertronic, USA) during their production.

The collected recordings were treated with an in-house MATLAB program. The code uses a Circular Hough Transform based algorithm for locating the droplets and then estimating their diameter along with the polydispersity index (PDI). This approach was chosen owing to its robustness in the presence of noise, occlusion and varying illumination. An example is given in Fig. S1 for water single-core droplets where the inner radius (R w ) and outer radius (R) diameter were extracted using MATLAB. The water volume fraction ϕ w was derived from R w and R for water single-core PFH droplets.

For PFH nanodroplets dispersed in water (and stabilized by F 8 -TAC 13 ) or for water nanodroplets dispersed in PFH (and stabilized by Krytox) and before their encapsulation into PFH droplets, the size and polydispersity were determined by dynamic light scattering. The measurements were performed on emulsions diluted 100 times, using an ALV/CGS-3 platform based goniometer system (from ALV GmbH), at room temperature. An experiment consists of measurements of the same solution at scattering angles, θ, ranging from 50°to 160°, with a step of 10°. At each angle θ, the device provides the decay rate

Γ θ = q 2 k B T 6πηR ,
where k B is the Boltzmann constant, T is the temperature in K, η is the viscosity of the solvent, and q(θ) = (4πn s /λ) sin(θ/2) is the magnitude of the scattering vector. The refractive index of the solvent is n s = 1.33 and λ = 633 nm is the laser wavelength. A fit of the curve by the cumulant method made it possible to determine the hydrodynamic droplet mean radius R along with the polydispersity index PDI [START_REF] Mailer | Particle sizing by dynamic light scattering: non-linear cumulant analysis[END_REF] . A typical example of Γ θ versus q 2 plot is presented in Fig. S2 for plain PFH droplets which exhibit a diameter of 149 nm. was made to emit long and low intensity acoustic pulses and the amplitude of the output signal was monitored and recorded by the hydrophone. The tube/hydrophone position was varied thanks to x-y-z motored stages and we selected the position corresponding to the maximum amplitude recorded by the hydrophone to place the tube at the focal point.

Acoustic droplet vaporization experiments

Acoustic setup

In addition, we used a Fabry-Pérot fiber-optic ultrasonic hydrophone (Precision Acoustics, Dorchester, UK) for the measurement of temperature. [START_REF] Morris | A Fabry-Pérot fiber-optic ultrasonic hydrophone for the simultaneous measurement of temperature and acoustic pressure[END_REF] In our case, temperature variations never exceeded 0.5 ℃ at the applied pressures.

Acoustic parameters

The signal consisted of sine-wave bursts at fundamental frequency 1.1 MHz. We employed short bursts of 5-cycles sine wave gated by a rectangular window of 91 µs pulses at an interval of 3 s to allow sufficient time to record the output. To subsidize any thermal effects we fixed the duty cycle at 5% and a pulse repetition frequency at 11 kHz. The total time of an experiment was 7 min. Acoustic peak negative pressures varied from |P | = 0.3 to 7 MPa.

Samples preparation

Due to their high density, the droplets sediment quickly when dispersed in water. In order to guarantee a homogeneous dispersion of the droplets during the duration of the experiment, glycerol was chosen as the bulk phase because of its viscosity of 11.1 mPa s.

ADV signal processing

The method of measuring the ADV threshold is based on the principle that microbubbles (i.e., vaporized droplets) are much more powerful acoustic scatterers than liquid droplets because of the high compressibility of gases compared to liquids. [START_REF] Holland | Thresholds for transient cavitation produced by pulsed ultrasound in a controlled nuclei environment[END_REF][START_REF] Apfel | Activatable infusable dispersions containing drops of a superheated liquid for methods of therapy and diagnosis[END_REF] When excited at a given frequency f 0 , a bubble not only oscillates at this same resonant frequency, but also at the subharmonic f 0 /2, 28 and at higher harmonics if 0 (where i is an integer larger than 1). However, higher harmonics if 0 can also appear due to the non-linear propagation of an acoustic wave. [START_REF] Muir | Prediction of nonlinear acoustic effects at biomedical frequencies and intensities[END_REF] Hence, the appearance and magnitude of the subharmonic peaks at half the fundamental frequency, here f 0 /2 = 0.55 MHz, was used as the most relevant indicator to detect the onset of ADV. For instance, Fig. S3.A displays the signal emitted by a glycerol solution containing plain PFH droplets with a radius of 20 µm that were insonified by a sound wave of amplitude 2 MPa, sufficient enough to trigger ADV. This is evidenced by the magnified level of the subharmonic at 0.55 MHz (grey solid line) in Fig. S3.B compared to pure glycerol (black dotted line).

The determination of ADV threshold can be achieved using two methods that will be compared to each other hereafter. To do so, three different samples were used for each experiment. For each sample, the amplitude of the emitted signal was progressively increased, at 25 mV intervals, at the waveform generator. For each increment, the pressure wave scattered by the droplets was collected by the hydrophone. 100 pulses were sent onto each sample, leading to the recording of 100 scattered signals, which were then stored for each sample at each applied pressure. During acquisition, the oscilloscope performed analog-to-digital conversion of the signal at a sample rate of 1 GHz. Once recorded, the Fast Fourier Transform (FFT) of each acquired signal was calculated using MATLAB. Therefore, the data points in the charts shown in this article are all averages over 300 FFT data.

Method A: We first derived the FFT spectra of the signal emitted from a droplet sample at an applied acoustic pressure P . Then, we calculated the power spectrum V sub around the subharmonic peak f 0 /2 (from 0.4 to 0.6 MHz). According to Maxwell et al., [START_REF] Maxwell | Probability of Cavitation for Single Ultrasound Pulses Applied to Tissues and Tissue-Mimicking Materials[END_REF] the value of V sub depends linearly on pressure in the absence of cavitation (or vaporization). Consequently, we normalized the value V sub by the applied pressure times a constant so that a plateau appeared at low pressures. Any deviation from this plateau would be characteristic of the occurrence of cavitation. The data points were fitted using a piecewise linear function (see Fig. 4). The pressure at the intersection point of the two linear fits indicated the ADV threshold [START_REF] Aliabouzar | Effects of droplet size and perfluorocarbon boiling point on the frequency dependence of acoustic vaporization threshold[END_REF][START_REF] Aliabouzar | Acoustic Droplet Vaporization in Acoustically Responsive Scaffolds: Effects of Frequency of Excitation, Volume Fraction and Threshold Determination Method[END_REF][START_REF] Radhakrishnan | Relationship between cavitation and loss of echogenicity from ultrasound contrast agents[END_REF] .

Method B:

The differentiation between the presence and absence of bubble formation is made by taking the integral I of the Fourier Transform (IFT) in the subharmonic range between 0.4 and 0.6 MHz 3,10 of the droplet sample I vap (P ) and a control sample of glycerol I gly (P ) containing no droplets, at each applied pressure. Vaporization occurrence at a particular acoustic pressure P was determined according to the criterion,

I vap (P ) ≥ I gly (P ) + ζσ, (1)
where ζ is an integer and σ is the standard deviation value of the values of I gly (P ), recorded from glycerol at pressure P , compared to its means value. Each time the condition in Eq. 1 was met, a vaporization event was counted. As the phenomenon is stochastic, its probability, p (n) , was obtained by repeating the experiment 100 times at the same pressure value and for the same sample. The probability, p (n) , of bubbles appearing within a sample, for a given pressure P , was defined as the frequency of vaporization, i.e., the sum of vaporization events observed divided by 100. This value was then averaged over three independent experiments using three different samples.

We assumed that, in our suspension of droplets, each droplet vaporizes independently of the other with a probability p [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] . The number of vaporization events from n droplets is then given by a binomial distribution with parameters n and p [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] , and the probability to observe at least one vaporization by

p (n) ≥1 = 1 -1 -p (1) n . (2)
When n becomes large, p

≥1 is approximated by a minimal Gumbel distribution with location parameters µ (n) and scale parameters

β (n) p (n) ≥1 ≈ 1 -exp -e P -µ (n) β (n) . (3)
The median of this distribution is P

(n) 0.5 and can be expressed as

P (n) 0.5 = µ (n) + β (n) ln(ln(2)). (4)
The probability p (n) to observe a vaporization event when n droplets are insonified is fitted by Eq. 3 to obtain the value of P (n) 0.5 (using Eq. 4). However, the determination of the P

Results and Discussion

Interfacial tensions

Several interfaces are at play in our system, each of them characterized by an interfacial tension γ. Since nucleation models depend on these interfacial tensions, we measured their values (given in Table 2) at equilibrium for all the configurations met in our droplet systems. This includes the interfacial tensions between liquid PFH and glycerol (γ lc ), liquid PFH and water (γ lw ), gaseous PFH and glycerol (γ gc ), gaseous PFH and water (γ gw ), and finally gaseous and liquid PFH (γ gl ). Note that our value for γ gl is very close to the one measured by Luís et al. [START_REF] Luís | Influence of Nanosegregation on the Surface Tension of Fluorinated Ionic Liquids[END_REF] (that is 11.7 mN/m). Values were measured in the presence of F 8 -TAC 13 when using glycerol or in the presence of Krytox 157 FSL when using water. Since droplets are of millimetric size in these experiments, the measured values are close to those of a flat surface (i.e., R ≈ ∞) and we consider that γ lc = γ lc (∞), γ lw = γ lw (∞), γ gc = γ gc (∞), γ gw = γ gw (∞), and γ gl = γ gl (∞).

However, the effective interfacial tension are significantly different from these values for droplets with a radius of several hundreds of Angström according to Tolman: 34

γ(r) = γ(∞) 1 + 2δ r , (5)
where r is the radius of the droplet, γ(∞) is the interfacial tension for a flat interface (i.e., r = ∞), δ is the Tolman length, whose value is on the order of a few Angström. [START_REF] Holten | Homogeneous nucleation of water between 200 and 240 K: New wave tube data and estimation of the Tolman length[END_REF] Consequently, we expect γ gw and γ gl to strongly depend on the nucleus size, while γ lw should not significantly vary because δ/R 1 in the experimentally explored range of droplet radii.

Cavitation pressure threshold in glycerol

The pressure threshold that induces cavitation in pure glycerol was evaluated in the absence of droplets. In these experiments, all glycerol solutions were vortexed at a speed of 2000 rpm for 30 s before measurements. Using the Method A described in section "ADV signal processing", we measured the normalized subharmonic peak amplitude of the acoustic signal of pure glycerol for an acoustic peak negative pressure P whose absolute value increases from

Assessment of the absence of superharmonic focusing

We followed the procedure described by Shpak et al. [START_REF] Shpak | Acoustic droplet vaporization is initiated by superharmonic focusing[END_REF] We first used the software HIFUbeam Simulator 37 to evaluate the shape of the acoustic wave produced by our transducer in the absence of droplets. In this Matlab script, we used as input parameters, the frequency (f = 1.1 MHz), the geometric focus (63.2 mm) and radius (64.00 mm) of our transducer, as

given by the manufacturer Sonic Concepts. The acoustic pressure wave was evaluated at the transducer focus and then fitted by a Fourier series 36

P 0 (t) = ∞ n=0
a n e i(nωt+φn) ,

where ω = 2π/f , φ n is a phase term, and a summation up to 10 was more than enough to fit the acoustic wave. We used the values a n and φ n to evaluate the contribution to the acoustic pressure P d due to the presence of the droplet 36

P d (r, θ, t) = ∞ n=0 ∞ m=0
a n e i(nωt+φn) α mn j m (nk 1 r)P m (cos θ),

where j m is the spherical Bessel function of the first kind of order m, P m is the Legendre polynomial of order m, and

α mn = χ m j m (x 0 )h (2) m (x 0 ) -h (2) m (x 0 )j m (x 0 ) j m (x 1 )h (2) m (x 0 ) -k 1 ρ 0 k 0 ρ 1 h (2) m (x 0 )j m (x 1 ) , (8)
where χ m = (-i) m (2m + 1), x 0 = nk 0 R, x 1 = nk 1 R, and h

m is the spherical Hankel function of the second kind of order m. Using Eq. 7, we looked at the moment where the pressure reached its minimum value inside the droplet. For all combinations of droplet radius and acoustic pressure tested in our experiments, we observed no superharmonic focusing. An example of curve snapshot for the largest droplet radius (of 30 µm) is given in Fig. S4.

Methods to determine ADV pressure threshold P 0.5

We performed measurement on solutions of plain PFH droplets with a radius of 20 µm. We used Method A and B to determine the value P events associated with the droplets is expected to be much smaller than the number of cavitation events occurring in glycerol. Thus, while the precision of Method A is adequate to detect cavitation in glycerol, the precision of Method B is more adapted for the detection of vaporization event.

P 0.5 as a function of plain PFH droplet number, n

We performed measurements on solutions of plain PFH droplets, whose radius were 20 µm, for droplet volume fractions ϕ d varying from 2.5×10 -4 up to 50×10 -4 . From these measurements, we derived the ADV pressure threshold P (n) 0.5 using method B. In these experiments, only the droplets located in the volume v focus of the transducer focus may be vaporized. This volume is equal to 0.38 µL (determined at -6 dB according to the manufacturer specifications) and contained n = 3, 6, 11, 28, and 57 droplets for ϕ d equal to 2.5 × 10 -4 , 5 × 10 -4 , 10 × 10 -4 , 25 × 10 -4 , and 50 × 10 -4 , respectively. We observe in Fig. 6A that the probability curves are shifted upon the addition of droplets. The fit of these probability curves provide the values of P (n) 0.5 , displayed in Fig. 6B, which represents the values of pressure at which the probability to observe an ADV event is p (n) = 0.5. The value of |P n 0.5 | decreases as the number of insonified droplets n increases. Table 3: Results from experiments performed on a solution of PFH plain droplets dispersed in glycerol, at a droplet volume fraction ϕ d = 5 × 10 -4 . R is the radius of the droplets, n tot is the total number of droplets in the sample tube, while n is the number of droplets inside the acoustic focus volume (determined at -6 dB). P

(1) 0.5 and P (n) 0.5 are the derived pressure at which the probability to observe a vaporization event is 0.5, when 1 or n droplets are respectively present in the acoustic focus volume. σ 1 is the standard deviation on the gaussian distribution having P (1) 0.5 as a median For each radius, we determined the probability p (n) to observe a vaporization event as a function of acoustic peak negative pressure P . From these probability curves, we derived P (n) 0.5 , whose values are displayed in Table 3 and Fig. 7. The values slightly decrease as the droplet radius increases, while the volume droplet fraction is constant (i.e., the corresponding number of droplets increases). P 0.5 as a function of volume fraction of encapsulated water, ϕ w , in multi-core

R n tot n |P (n) 0.5 | |P (1) 0.5 | σ 1 (µm) (MPa) (MPa) (MPa) 0.1 2.4 ×

PFH droplets

We performed experiments with suspensions of PFH droplets (of radius R = 20 µm and of volume fraction ϕ d = 5 × 10 -4 ) comprising many water droplets (of radius R w = 0.215 µm).

The concentration of water droplets was varied so that the volume fraction of water in a droplet ϕ w ranges from 0.2 to 0.8. The ADV pressure thresholds, determined for the various values of ϕ w , are displayed in Fig. 8 along with the value for plain PFH droplets which correspond to ϕ w = 0. We observe that the values of |P P 0.5 between water single-core and multi-core PFH droplets

In the following experiments, we consider suspensions of droplets (of radius R = 20 µm and of volume concentration ϕ d = 5×10 -4 ) made of a single water droplet (of radius R w = 10.5 µm).

The water volume fraction is ϕ w ≈ 0.15. We derived a value of |P (n) 0.5 | = 1.6 ± 0.1 MPa from our probability curve and the value is plotted in Fig. 8 in order to compare it with values determined for multi-core water PFH droplets. We observe that the value |P (n) 0.5 | is much smaller for water single-core droplet than for multi-core ones.

Extrapolating ADV pressure threshold from a suspension of insonified droplets to a single insonified droplet Since experiments were performed on droplet suspensions, the value P (n) 0.5 should depend on the number of insonified droplets and therefore be different from the value P (1) 0.5 that would be obtained if only one droplet was insonified. Our experiments on an increasing concentration of plain PFH droplets of radius 20 µm confirmed the dependence of P (n) 0.5 on n (see Fig. 6B). This result means that comparisons between different measurements can be drawn only if they were performed at the same number of insonified droplets. This is a critical point as most experiments in the literature do not specify or control this number.

To understand the ADV mechanism, it is crucial to determine the value of P (1) 0.5 instead of P (n) 0.5 , as all nucleation models to date, were developed for the case of one droplet. It is possible to achieve this either by diluting the droplet solution until we reach P (1) 0.5 or by directly calculating the value P (1) 0.5 from the probability curves p (n) . Indeed, it can be shown that if the probability p [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] , which is a function of acoustic pressure P , is given by an integral Gaussian distribution with median P (1) 0.5 and standard deviation σ [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] , then µ (n) and β (n) (given in Eq. 4) can be expressed as

µ (n) = P (1) 0.5 -σ (1) d n , (9)
and

β (n) = σ (1) c n , (10)
with

c n = 1 2 ln(n) , (11)
and

d n = 2 ln(n) - ln(ln(n)) + ln(4π) 2 2 ln(n) . ( 12
)
By replacing the value of µ (n) given by Eq. 9 in the Eq. 3 we also obtain a relation between

P (n)
0.5 and P

(1) 0.5

P (n) 0.5 P (1)
0.5 -σ [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] [d n + c n ln(ln( 2))] .

Eq. 3 can be used to fit probability curves in order to determine either P

0.5 or P

(n) 0.5 , when using respectively Eq. 9-12 or 4. While Eq. 13 can be used to extrapolate the value of P [START_REF] Crouse | Opacification and border delineation improvement in patients with suboptimal endocardial border definition in routine echocardio-graphy: Results of the phase III albunex multicenter trial[END_REF] = 1.0 ± 0.3 MPa) is obtained from the fit of all the probability curves depicted in Fig. 6A, except for n = 57 (i.e., magenta curve) where the probability curve could not be correctly fitted. Both these methods are therefore equivalent in providing a value for P (1) 0.5 . However, the second method offers the advantage of requiring the estimation of just a single probability curve for only one value of n.

We used results collected from the experiments performed at a constant volume fraction of droplet but with various droplet radii (see Fig. 7), to derive |P

Failure of the homogeneous nucleation model to predict our results

In a homogeneous nucleation, a nucleus appears randomly inside the droplet volume. The probability q ≥1 that at least 1 nucleus appears is given by

q ≥1 = 1 -e -N (τ ) , (14)
where N (τ ) is the average number of nucleus during a time τ and within a volume V . This can be expressed as N (τ ) = JV τ , J being the volumic rate of nucleation J = J 0 e

- W hom (q ≥1 ) k B T
and J 0 a constant. [START_REF] Blander | Bubble nucleation in liquids[END_REF] Consequently, we have

W hom (q ≥1 ) = k B T ln   J 0 V τ ln 1 1-q ≥1   . ( 15
)
We suppose that the vaporization occurs during the half period where the acoustic pressure is negative, thus τ = 1 2f . In addition, the energy required to create a spherical gas volume, so-called nucleus, of radius r is [START_REF] Blander | Bubble nucleation in liquids[END_REF][START_REF] Katz | Condensation and boiling: Corrections to homogeneous nucleation theory for nonideal gases[END_REF] W hom (r) = 4πr 2 γ gl (r) -

4 3 πr 3 (P g -P l ) + n PFH (µ g -µ l ), (16)
where γ gl is the surface tension of the nuclei surface (interface between gas and liquid), P g and P l are the pressures in the gaseous and liquid PFH, n PFH is the number of PFH molecules, and µ g are µ l the chemical potential of the gaseous and liquid PFH. The fate of a nucleus depends on its radius r compared to the critical radius r * for which the nucleus has equal chance to grow or to shrink (with µ g = µ l ):

19 r * = 2γ gl (r * ) P v -P * l , (17)
where P v and P * l are the vapor pressure of PFH and the critical pressure leading to the nucleus appearance. The nucleation energy W hom (r) required to obtain a nucleus of critical radius r = r * is then

W hom (r * ) = 16πγ 3 gl (r * ) 3(P v -P * l ) 2 , (18)
Since at r * the probability of the nucleus to grow is 1/2, the probability to obtain a vaporization of the droplet is also 1/2, thus W hom (r = r * ) = W hom (q ≥1 = 0.5). Taking into account the Laplace equation between the pressure inside and outside the droplet (respectively P * l and P

(1) 0.5 )

P * l = P (1) 0
.5 + 2γ lc R . (19)
After equalizing Eq. 18 and 15 and replacing V = 4 3 πR 3 , we obtain

P (1) 0.5 = P v - 16πγ 3 gl (r * ) 3k B T ln 2πJ 0 R 3 3f ln(2) - 2γ lc R . (20)
Since all parameters in Eq. 20 are known, this equation can be used to estimate the values

P (1)
0.5 for various droplet radii. On doing so, the equation predicts values of |P

(1) 0.5 | greater than 12 MPa whatever the droplet radius may be, which are much larger than our experimental values as shown in Fig. 9. In addition, the homogeneous nucleation model predicts a smaller variation of P (1) 0.5 compared to the experimental values in the range of investigated droplet radius. The homogeneous model can also be applied to the case of water single-core and multi-core PFH droplets, by considering the total volume V of PFH inside the droplet in Eq. 15, thus excluding the water volume which cannot participate in the vaporization. By doing so, the model predicts that (i) the ADV pressure threshold should not change between water single-core and multi-core PFH droplets as the volume of PFH is identical, (ii) the threshold should diminish when the volume of PFH increase in water multi-core PFH droplets. These two predictions are opposite to what was experimentally observed in Fig. 8. In summary, the model of homogeneous nucleation fails to predict the results of all of our experiments and consequently it is not the mechanism responsible for acoustic droplet vaporization in our experimental conditions.

Development of a model of heterogeneous nucleation and comparison with experimental data

In heterogeneous nucleation, a nucleus appears at an interface present in the system. In the case of plain PFH droplet, it will be the single interface between PFH and glycerol. Similar to the homogeneous nucleation, the energy of heterogeneous nucleation is given by,

W het (q ≥1 ) = k B T ln   Π 0 Aτ ln 1 1-q ≥1   , (21)
where Π 0 is the rate of heterogeneous nucleation 39 and A is the surface area of the interface on which the nucleus can appear. The work needed to create a nucleus is written as 39 (see Fig. 10)

W het (r 1 , r 2 ) =γ gl (r 1 )a(ψ, r 1 ) + γ gc (r 2 )a(ϕ, r 2 ) -γ lc a(χ, R) -(P g -P l )v + n PFH (µ g -µ l ), (22)
where a(ψ, r 1 ) is the surface area of the gaseous PFH nucleus in contact with liquid PFH, a(ϕ, r 2 ) is the surface area of the nucleus in contact with the surfactant droplet shell that was previously in contact with liquid PFH with a surface area a(χ, R). We have

a(ν, r) = 2πr 2 (1 -cos ν), (23)
where ν = ψ, ϕ, or χ, and r = r 1 , r 2 or R. The nucleus volume is

v = π 3 r 3 1 (2 -3 cos ψ + cos 3 ψ) + r 3 2 (2 -3 cos ϕ + cos 3 ϕ) . (24)
If the critical radii r * 1 and r * 2 are known, then the pressure P

(1) 0.5 could be predicted by solving the equation

W het (r * 1 , r * 2 ) -W (q ≥1 = 0.5) = 0. (25)
Since the nucleus is at equilibrium at r 1 = r * 1 and r 2 = r * 2 (where we have µ g = µ l and P g = P v ), the Young-Laplace equations can be written for the three interfaces: liquid PFH/glycerol, gaseous PFH/liquid PFH, and gaseous PFH/glycerol :

P v = P * l + 2γ gl (r * 1 ) r * 1 , (26)
P v = P * A + 2γ gc (r * 2 ) r * 2 , (27)
P * l = P * A + 2γ lc R . (28)
We can easily determine the critical radii from these equations

r * 1 = 2γ gl (r * 1 ) (P v -P * l ) , (29)
r * 2 = 2γ gc (r * 2 ) (P v -P * l ) + 2γ lc R . ( 30
)
The calculation of the surface area and volume of the nucleus depend on the critical angles χ * , ψ * and ϕ * (see Fig. 10). All the cosine of these angles can be written as a function on the cosine of the angle θ * = ψ * + ϕ * (see supplemental informations) leading to the occurrence of nucleation preferentially on one interface rather than the other one are mainly the droplet radius, the droplet surface area and the interfacial tensions.

cos ϕ * = r 2 + r 1 cos θ * r * 1 + r * 2 + 2r * 1 r * 2 cos θ * , ( 31
)
cos ψ * = r * 2 cos θ * + r * 1 r * 1 + r * 2 + 2r * 1 r * 2 cos θ * , ( 32
) cos χ * = 1 - r * 2 R 2   1 - r * 2 + r * 1 cos θ * r * 1 + r * 2 + 2r * 1 r * 2 cos θ * 2   . (33)
The experimental data displayed in Fig. 8, show no change in the values of ADV threshold as the number of encapsulated water droplets increases in water multi-core PFH droplets.

In addition, these values are similar to the value measured for a droplet devoid of water droplets (i.e., plain PFH droplets). If the nucleation would have occurred at the water/PFH interface of water droplets, a decrease in ADV threshold would have been expected because the total surface area of the interface water/PDH has increased (due to the increase in the number of water droplets) since it would lead to a higher probability of nucleation according to Eq. 21. We conclude that nucleation is preferentially occurring on the surface of the PFH/glycerol interface at a mean ADV pressure threshold of | P (1) 0.5 | = 2.8 ± 0.1 MPa. This value is characteristic of the properties of the PFH/glycerol interfaces (i.e., radius, surface area, and interfacial tensions) which are identical between all the droplet types.

For a droplet encapsulating a micron-sized droplet, the value of ADV threshold is |P The difference observed between water single-core and multi-core PFH droplets comes from the difference in the water droplet radius. Indeed, our model indicates that the energy required to form a gas nucleus at a surface of water micron-sized droplet is lower than at a surface exhibiting a larger curvature. Our results also suggest that the energy required to form a nucleus at the water/PFH interface of water nanometric-sized droplets is larger than the one at the PFH/glycerol interface (of the micron-size droplet encapsulated the water droplets). In this case, the difference arises not only due to a difference in surface curvature, but also because the surfactants (hence the respective interfacial tensions) are different between the two interfaces.

Conclusion

The acoustic vaporization pressure threshold P (n) 0.5 strongly depends on the number of insonified droplets. We showed that the threshold value for only one insonified droplet, P (1) 0.5 , can be obtained either by diluting the droplet solution or by using a statistical model which can derive this value from the probability curves measured on droplet solution.

When considering plain PFH droplets, the P (1) 0.5 value decreases as the radius of the plain PFH droplets increases. The threshold P (1) 0.5 can be reduced by encapsulating a micron-sized water droplet. However, the encapsulation of many nanometer-sized water droplets induces no modification in P (1) 0.5 value. Our experimental conditions prevent the occurrence of superharmonic focusing. The experimental data performed on plain PFH droplets, water single and multi-core PFH droplets can not be explained by homogeneous nucleation and the model of homogeneous nucleation fails to predict the data of plain PFH droplets. However, all experimental results can be explained by heterogeneous nucleation, where the nucleation occurs either on the PFH/glycerol surface of the droplet or on the PFH/water surface of the encapsulated water droplets. We derived equations modeling heterogeneous nucleation that were successful to describe the behavior of P (1) 0.5 as a function of the radius of plain PFH droplet. The model uses as parameters three interfacial tensions, the droplet radii, and the vapor pressure of the perfluorocarbon. Thus, the model can be used for any kind of perfluorocarbons and any temperature as long as the values of these parameters are known.

  Methanol and NaCl were obtained from Sigma Aldrich. PDMS Sylgard 184 was purchased from Neyco. Photoresist SU8 was purchased from Chimie Tech Services. Silicon wafers were obtained from BT Electronics. The pressure controllers (model MFCS-EZ) used to inject the fluids into the microfluidic chips were purchased from Fluigent. The tubing connecting the pressure controllers to the microfluidic chip was bought from VWR (Internal diameter: 0.51 mm, External diameter: 1.52 mm). The channel inlets and outlets were punched with a 0.35 mm biopsy punchers delivered by World Precision Instruments (WPI); while the plasma cleaner was purchased from Harrick Scientific. The perfluorohexane PFH and the fluorinated surfactant, Krytox 157 FSL, were purchased respectively from ABCR GmbH and Chemours.

Figure 1 :

 1 Figure 1: (A) Chemical structure of F 8 -TAC 13 surfactant. (B) General synthetic pathway of F-TAC surfactants (reagents and conditions: AIBN/CH 3 OH, reflux 62-73%).

Figure 2 :

 2 Figure 2: Schematic representation of one-(top and middle schemes) and two-nozzle (bottom scheme) microfluidics chips used to produce respectively plain PFH droplets, water multiple core PFH droplets and water single-core PFH droplets. The green and blue colors filling the channels are for perfluorohexane and water, respectively. The black line around droplets indicates the presence of F 8 -TAC 13 surfactant, while the red line indicates the presence of Krytox surfactant. For each droplet type, a picture of the emulsion taken under microscope is given.

Fig. 3 2 ℃

 32 Fig. 3 sketches the ultrasonic setup triggering the ADV. A waveform generator (Model 33220A from Agilent) generates an electrical signal, that first goes through a radio-frequency power amplifier (Model 150A100C from AR France), and then through a power reflection meter (Model & NRT from Rohde Schwarz) that measures the delivered average electrical power. The electrical signal is converted into an acoustic wave by a focused transducer (Model H-101-G from Sonic concepts Inc., with a f-number of 0.99), whose fundamental mode is at 1.1 MHz. The acoustic wave propagates into a water tank thermostated at 20±0.2 ℃ and connected to a degassing machine (Model WDS-1005 from Sonic Concepts), which reduces the concentration of dissolved oxygen below 2 mg/L. A PCR tube, filled with 200 µL of suspension, is placed at the transducer focus. A hydrophone (model Y-107 from Sonic Concepts), whose focus overlaps with that of the transducer's inside the tube, continuously monitors the scattered signal. The positioning of the PCR tube at the transducer focus was performed using a 0.2 mm needle hydrophone (from Precision acoustics) placed inside an open tube. The transducer

Figure 3 :

 3 Figure 3: Experimental setup for ADV and measuring cavitation. A 1.1 MHz transducer (T) vaporized the droplets while the hydrophone (H) recorded at 0.55 MHz the scattered emissions due to subharmonic response of vapor micro bubbles in the sample (S).

  The samples were prepared by taking, for instance, a volume of v = 0.1 µL from the bottom of a solution of droplets in water (so that the pipetted volume is mostly made of droplets) and by dispersing it in a volume of 200 µL of degassed glycerol in a PCR tube. In this case, the volume fraction of droplets was approximately ϕ d = 5 × 10 -4 . Other values of ϕ d were obtained by varying the volume v of droplets added to the 200 µL of glycerol. To ensure a homogeneous distribution of the droplets before each experimental test, the tube was vortexed at 2000 rpm for 30 s. Table3lists the different sizes and the corresponding numbers of droplets dispersed in 200 µL of glycerol solution. A solution of fresh droplets was prepared before each experiment.

Figure 4 :

 4 Figure 4: Normalized subharmonic component peak amplitude of the scattered response of a glycerol solution devoid of droplets (figure A) and of a suspension of PFH droplets (R = 20 µm, ϕ d = 5 × 10 -4 ) in glycerol (figure B) as a function of applied acoustic pressure P .The two lines intersect gives |P | = 4.5 ± 0.3 MPa for the figure A and |P | = 2.1 ± 0.3 MPa for the figure B. Consequently, the vaporization of PFH droplets appears before glycerol cavitation when pressure increases.

  the chosen value for ζ.[START_REF] Fabiilli | The role of inertial cavitation in acoustic droplet vaporization[END_REF] It is expected that when ζ is low, fluctuations of the background, at f 0 /2, may contribute to the events ascribed to vaporization, which is introducing a bias in P (n) 0.5 determination. In contrast if ζ is large, we may fail to detect vaporization events. Fig. 5, displays the probability curves determined for a suspension of plain PFH droplets (R = 20 µm and ϕ d = 5 × 10 -4 ) for ζ varying from 5 to 10. For ζ = 10, some vaporization events are obviously missed as p does not remain on a plateau after reaching 1. A plateau at p = 1 is rapidly reached when ζ = 3, 5 and 6, compared to larger values of ζ, suggesting that signals not associated with a vaporization event are taking into account in the derivation of p. We thus chose a value of 7 for ζ as it was the value exhibiting no bias when considering our whole experimental data, while at ζ = 9 a decrease in the plateau at 1 could be observed in some of our experiments.

Figure 5 :

 5 Figure 5: The probability p (n) to detect an ADV event was evaluated as a function of the applied peak negative acoustic pressure P using different values for ζ (3, 5, 6, 7, 9, and 10 from the left to the right probability curves). From the respective curves, we determined the following P (n) 0.5 values: 1.0, 1.7, 1.8, 2.1, and 2.4 MPa. The experiments were performed on a suspension of plain PFH droplets having a radius R = 20 µm, and at a droplet volume fraction ϕ d = 5 × 10 -4 .

1 up to 5 . 2

 52 MPa, as shown in Fig. 4A. In these experiments, a different glycerol solution was used for each measurement. We observe from Fig. 4B that cavitation starts to occur for |P | > 4.5 MPa in pure glycerol. Consequently, all of our measurements on determining the ADV pressure threshold were performed at |P | < 4.5 MPa, hence, in the absence of cavitation occurring in glycerol.

  of ADV pressure threshold. Fig.4Bshows that Method A gave the value |P (n) 0.5 | = 2.1 ± 0.3 MPa, which is at the point of intersection of the linear piecewise fit. While Method B yielded |P (n) 0.5 | = 2.12 ± 0.03 MPa as shown in Fig. 5 from the fit of the data to Eq. 3 and the use of Eq. 4. Consequently, both methods gave similar ADV threshold values P (n) 0.5 . The error obtained by Method B was noticeably smaller than Method A and for this reason we used Method B to derive the probability of vaporization in the following experiments. Note that the number of detected vaporization

P 0. 5

 5 as a function of plain PFH droplet radius, R, at a constant droplet volume fraction, ϕ d Next measurements were performed at a constant volume fraction of plain droplets, ϕ d = 5 × 10 -4 , irrespective of the droplet radius. The droplets have a radius of either 0.1, 2,5, 20 or 30 µm. As a consequence, the number of droplets n located in the volume of the transducer focus varied. We calculated that n was respectively equal to 4.5 × 10 7 , 2.9 × 10 3 ,

Figure 6 :

 6 Figure 6: (A) The circles represent the probability p to measure a vaporization event as a function of pressure for solutions containing droplets of radius 20 µm and at various droplet volume fraction ϕ d (see text), leading to a number of droplets in the transducer focus of 3 ( ), 6( ), 11 ( ), 28 ( ), and 57 ( ).The lines are fits performed using either Eq. 3 or 13, as there is no difference between the two fits. (B) The values |P (n) 0.5 | obtained from the previous fit are plotted as a function of n, the number of droplets localized in the transducer focus. The line is a fit (R 2 = 0.88) using Eq. 13, 11, and 12, for which |P (1) 0.5 | = 3.3 ± 0.4 MPa and σ (1) = 1.2 ± 0.3 MPa.

5 |

 5 vary around 2.2 MPa (dotted line in Fig.8) for the water multi-core PFH droplets, close to the value |P (n) 0.5 | for plain PFH droplet when ϕ d = 5 × 10 -4 .

Figure 7 :

 7 Figure 7: Values of |P (1) 0.5 | (♦) and |P (n) 0.5 | (•) determined for solutions at constant volume fraction of droplets (ϕ d = 5 × 10 -4 ), with droplet radius varying from 0.1 to 30 µm.

Figure 8 :

 8 Figure 8: The squares represent the values of |P (n) 0.5 | for a suspension of either water multicore PFH droplets ( ) or water single-core PFH droplets ( ). While the lozenges represent the values of |P (1) 0.5 | for a suspension of either water multi-core PFH droplets (♦) or water single-core PFH droplets ( ). All suspensions contain an identical droplet volume fraction of ϕ d = 5 × 10 -4 and an identical radius of R = 20 µm. The water droplets are stabilized by Krytox 157 FSL and possess a radius of either R w = 10.5 µm or 250 nm when a PFH droplet contains, respectively, a single water droplet ( and ) or many of them ( and ♦).

5

 5 values determined at different n. In another words, we can extrapolate the value of P (1) 0.5 through two methods: the first one by fitting P (n) 0.5 values determined at several droplet concentrations, or by fitting a probability curve at a known n. Using the first method, the fit of the P (n) 0.5 values illustrated in Fig. 6B leads to |P (1) 0.5 | = 3.3 ± 0.4 MPa (and σ (1) = 1.1 ± 0.3 MPa). While for the second method, an average value of |P (1) 0.5 | = 3.1 ± 0.5 MPa (and σ

5 |

 5 in order to compare our data with nucleation models. We observe that the values of |P (1) 0.5 | and |P (n) 0.5 | are very different when n is very large, but they converge to the same value where n = 3. It is notable that the values of |P (1) 0.5 | decreases from 8.5 down to 1.5 MPa when the droplet radius increases from 0.1 to 30 µm. A correct model of nucleation should be able to predict this variation.

Figure 9 :

 9 Figure 9: The diamond are the values |P (1) 0.5 | previously given in Fig. 7. The dotted line represents the prediction of the homogeneous nucleation model. While the solid line is the prediction given by the heterogeneous nucleation model where γ gc is the only unknown parameter, whose values are given in Fig. 11B.

Figure 10 :Figure 11 :

 1011 Figure 10: Schema of a nucleus on the surface of a plain PFH droplet surface (represented by the black solid line). The surface of the nucleus are represented by the green (with no surfactants) and orange (covered with surfactants) solid lines.

Figure 12 :

 12 Figure 12: The figure displays the shape of the nucleus having radii r * 1 and r * 2 as showed by the solid lines. The dotted lines show the curvature of the droplet with radius R, where the liquid PFH is located on the top part of the graph. The nucleus size changes with the droplet radius R, where the colors green, red, blue, brown and yellow are for R = 30, 20, 5, 2.5, and 0.1 µm.

1

 1 MPa, thus smaller than 2.8±0.1 MPa. This suggests that nucleation is preferentially occurring at the surface of the PFH/water interface of the encapsulated water droplet rather than of the PFH/glycerol interface.

Table 1 :

 1 Nozzle diameters, flow rates of the dispersed phase Θ d and of the continuous phase Θ c used in the microfluidic chip to produce plain PFH droplets of radius R.

	R	Θ d	Θ c	nozzle size
	(µm) (µL/min) (µL/min)	(µm)
	2.5	0.3	0.8	10
	5	0.9	8	10
	20	0.9	8	40
	30	1	2.5	40

Table 2 :

 2 Values of interfacial tensions, where the subscripted letters g, l, c, w are for gaseous PFH, liquid PFH, glycerol and water, respectively. The first five values are of interest for our configuration. The next five values are controls or intended as a comparison between water and glycerol. The error in the values is ±0.002 N/m.

	notation interfacial tension	surfactant
		(N/m)
	γ gl	0.012	none
	γ lc	0.014	F 8 -TAC 13
	γ gc	0.033	F 8 -TAC 13
	γ lw	0.018	Krytox 157 FSL
	γ gw	0.066	Krytox 157 FSL
	γ gw	0.072	none
	γ lc	0.026	none
	γ gc	0.062	none
	γ lw	0.025	F 8 -TAC 13
	γ gw	0.038	F 8 -TAC 13

  10 10 4.5 × 10 7 2.51 ± 0.02 8.5 ± 0.6 1.08 ± 0.11 2.5 1.5 × 10 6 2.9 × 10 3 2.42 ± 0.01 5.8 ± 0.2 0.86 ± 0.05 5 1.9 × 10 5 3.6 × 10 2 2.39 ± 0.01 4.9 ± 0.2 0.86 ± 0.

				05
	20	2.9 × 10 3	6	2.12 ± 0.01 3.1 ± 0.1 0.75 ± 0.09
	30	8.8 × 10 2	2	1.41 ± 0.01 1.5 ± 0.1 0.21 ± 0.01
	3.6 × 10 2 , 6 and 2 droplets.

Acknowledgement

The authors acknowledge the funding from Plan cancer 2014-2019 (Project BubDrop4Glio).

The Ph.D scholarship of C. T. was funded by ITMO cancer.

In addition, the interfacial tensions should verify the equation (see Fig. 10) γ gl (r * 1 ) cos ψ * -γ lc cos χ * + γ gc (r * 2 ) cos ϕ * = 0.

It can be used it to determine the value of cos θ * as it can be written as a second degree equation a cos 2 θ * + b cos θ * + c = 0, where

cos θ * is thus simply

The value W het (r * 1 , r * 2 ) can be calculated numerically using the values r * 1 , r * 2 , and cos θ * , knowing that at equilibrium µ g = µ l and P * g = P v . However, in the presence of surfactant, an interfacial tension strongly depends on the surface concentration of surfactant. The appearance of a nucleus at the interface glycerol/PFC should lead to a local increase in surface area of the interface, and hence to a decrease in surfactant surface concentration. A variation in the value of γ gc (r *

2 ) is then expected. While no change in the values of γ lc (R) and γ gl (R) are expected since, respectively, the droplet radius R did not change and that no surfactant was located at the interface between gaseous and liquid PFH. Consequently, we numerically solved Eq. 25 by considering γ gc (r *

2 ) as an unknown parameter and using the experimental values P (1) 0.5 (from Fig. 7). By doing so, Eq. 25 can be numerically solved only if γ gc decreases as the nucleus radius r * 2 increases as displayed in Fig. 11A. The values r * 1 and r * 2 are plotted in Fig. 12 to provide the shape of the nucleus for the various plain PFH droplet radii. We observe that both radii decrease as R decreases, while the value |P * l | increases (see Fig. 11B). The nucleus volume is equal to v = 22, 41, 50, 104, and 175 nm 3 when the droplet radius is respectively equal to R = 0.1, 2.5, 5, 20 and 30 µm. Note that we used a Tolman length equal to 4Å. [START_REF] De Aguiar | The Interfacial Tension of Nanoscopic Oil Droplets in Water Is Hardly Affected by SDS Surfactant[END_REF] The system could not be solved using a larger value. Smaller values of δ only lead to smaller variation in nucleus size and γ gl (r * 2 ). We use the fit of γ gl (r *

2 ) values as a function of R (given in Fig. 11B) to predict the values of P (1) 0.5 by numerically solving Eq. 25. The predicted values are given in Fig. 9 as a solid line and it shows that the model of heterogeneous nucleation is able to predict the correct behavior of P (n) 0.5 as a function of the radius of plain PFH droplets. The prediction line is identical whatever the choice of δ, as it will be balanced by a change in γ gc values.

The question remains to understand the behavior of γ gc (r * 2 ). It cannot be the result of a curvature effect as given by the Tolman equation (dotted line in Fig. 11A) or of a decrease in surfactant surface concentration since an increase in γ gc (r * 2 ) would have been expected. But, it can be explained if some PFH molecules, close to the droplet surface, are getting embedded into the surfactant layer during the formation of a gaseous PFH nucleus. This insertion leads to a decrease in interfacial tension γ gc (r *

2 ) as the surface area is more covered. The insertion is possible if there is enough room for the surfactant layer to accommodate perfluorohexane molecules. The available space is decreasing as the curvature (i.e., 1/r * 2 ) is increasing, indeed in such a case the surfactant fluorocarbon chains are getting compact.

Consequently, the lower the curvature, the larger number of PFH can be inserted into the surfactant layer and the smallest is the value of γ gc (r *

2 ) The heterogeneous model can also be applied to the case of water single-core and multicore droplets. In this case, the nucleation can take place either at the surface of the interface in contact with glycerol, as previously, or at the surface of the encapsulated water droplet(s).

Similarly to plain PFH droplets, a model of heterogeneous nucleation can be written for the PFH/water interface of the encapsulated water by replacing R by the water droplet radius R w (with a change of sign in one of the Young-Laplace equation), and using the values γ lw and γ gw in the presence of Krytox (see Table 2)) instead of γ lc and γ gc . The parameters

Supporting Information Available

The following file is available free of charge.

• supplemental_materials: detailed derivations of equations used in the manuscript; figures displaying: size distribution of water single-core droplets, light scattering measurement on nanometric-sized PFH droplets, acoustic signal of samples, curves showing the absence of superharmonic focusing, schema with notation used in the previous equations; and finally a short bibliography.


H Amino

C Flageul

S Benhamadouche

M Ferrand

A SECOND ORDER CONSERVATIVE TIME-STAGGERED SCHEME FOR ALL SPEED FLOW

A time-staggered scheme for variable density flow is presented. The compressible Navier-Stokes equations are used by this pressure correction scheme which is implemented in the collocated finitevolumes open-source computational fluid dynamics solver code saturne. The Helmholtz equation is solved for the pressure increment, taking the thermodynamic pressure into account and avoiding the acoustic time step limitation. The internal energy equation is used to compute the temperature and a numerical analysis providing conditions ensuring the positivity of the thermodynamic variables is proposed. The scheme accuracy is first verified on a 0-D pressure cooker system. Then, a discontinuous shock problem shows its capability to reproduce shocks correctly and finally, the scheme is validated on a natural convection numerical reference case.

WORK OVERVIEW

Designing time schemes for Computational Fluid Dynamics (CFD) that are able to cover both compressible and incompressible regimes is a challenging [START_REF] Guillard | On the behaviour of upwind schemes in the low Mach number limit[END_REF] research topic which has been given attention in the last decades [START_REF] Degond | All speed scheme for the low Mach number limit of the isentropic Euler equations[END_REF][3][4] [START_REF] Boscheri | A second order all Mach number IMEX finite volume solver for the three dimensional Euler equations[END_REF] [START_REF] Park | Multiple pressure variables methods for fluid flow at all Mach numbers[END_REF]. Based on Pierce and Moin [START_REF] Pierce | Progress-variable approach for large-eddy simulation of non-premixed turbulent combustion[END_REF] finite difference time staggered scheme, we present a finite volumes conservative second-order time scheme for variable density flow. The proposed semi implicit pressure correction [START_REF] Guermond | An overview of projection methods for incompressible flows[END_REF] scheme uses the compressible Navier-Stokes equations and is designed for regular and discontinuous solutions. The latter is implemented in the open-source CFD solver code saturne [START_REF] Archambeau | Code Saturne: A finite volume code for the computation of turbulent incompressible flows-Industrial applications[END_REF] and validated against analytical and reference test cases.

Inner iterations are performed for each time step, characterised by the upper-script k and the scheme variable arrangement depends on the parameter θ (θ = 1 for Implicit Euler or 1/2 for Crank Nicolson): for instance, the momentum equation is solved between times n -1 + θ and n + θ, k:

Θ ρ n , ρ n+1,k-1 u k -Θ ρ n-1 , ρ n,k-1 u n-1+θ ∆t + div Θ u n-1+θ , u k ⊗ q n+θ,k-1 n-1+θ = ∇ p| n+θ,k-1 n-1+θ + div (τ k ), (1)
where ρ is the density, u the velocity, q = ρu the momentum and p| n+θ n-1+θ the mechanical pressure. Θ(•) is the temporal interpolation of a given field. If θ = 1/2 the velocity location is at the centre of the time interval [n, n + 1] whereas other fields are evaluated respectively on times n and n + 1 (Fig. 1). The internal energy equation is used to compute the temperature as in [START_REF] Herbin | A cell-centred pressure-correction scheme for the compressible Euler equations[END_REF]:

c v ρ n+1,k-1 T n+1,k -ρ n,k-1 T n ∆t + div Θ T n , T n+1,k q n+1,k-1 n = τ n+θ,k-1 : ∇u n+θ,k-1 + Γ u 2 /2 n+1,k-1 n + ∇ • λ∇Θ T n , T n+1,k -div Θ p n , p n+1,k-1 u n+θ,k-1 + u n+θ,k-1 • ∇ p| n+θ,k-1 n-1+θ (2)
T is the temperature, c v the heat capacity and λ the thermal conductivity. Following [START_REF] Herbin | A cell-centred pressure-correction scheme for the compressible Euler equations[END_REF], Γ u 2 /2 is a corrective source term derived from the discrete kinetic energy equation; for refined grids, the total energy *Corresponding author: martin.ferrand@edf.fr is conserved. The time scheme includes the effect of the thermodynamic pressure in the correction step by solving the following Helmholtz equation for p n+1,k :

ρ n+1,k -ρ n ∆t -θdiv (∆t∇p n+1,k ) = -div Θ ρ n , ρ n+1,k-1 u k + ∆t ∇ p| n+θ,k-1 n-1+θ + (1 -θ)div (∆t∇ p| n-1+θ n-2+θ ) ρ n+1,k = ρ k + p n+1,k -p n+1,k-1 ∂ρ ∂p T T n+1,k , p n+1,k-1
, ũk and ρk where are predicted variables.

(

) 3
The acoustic waves are treated implicitly, thus removing the acoustic CFL restriction on the time step. An upwind spatial discretization of the equations convective terms lead to a numerical analysis on the positivity of the temperature, pressure and density of a cell c; if the following conditions are below the unity, the positivity is preserved among time for all the aforementioned variables:

CFL + T 1 = (1 -θ) ∆t M n c ∑ f ∈F c Ṁ+ f n+1,k-1 n , Fo + T = λ c (1 -θ)∆t c v M n c ∑ f ∈F c |S f |, CFL + T 2 = ∆t M n c ∑ f ∈F c θ(γ -1) T n+1,k-1 c T n c + (1 -θ)γ Ṁ+ f n+1,k-1 n + ∆t(γ -1)u n+θ,k-1 c • ∇ c p| n+θ,k-1 n-1+θ p n c , CFL + p = ∆t M n c ∑ f ∈F c Θ ρ n , ρ n+1,k-1 u k + ∆t∇ p| n+θ,k-1 n-1+θ f • S f + (1 -θ)∆t∇ f p| n-1+θ n-2+θ • S f , (4)
where S f is the outward surface vector, Ṁ+

f n+1 n = max(0, Ṁ f n+1 n
) the positive contribution of the face mass flux and M c the cell mass.

x t T n-1 , Y n-1 , ρ n-1 T n ,Y n , ρ n T n+1 ,Y n+1 , ρ n+1 n -1 n n + 1 f 0 f 1 q n n-1 q n+1 n u n-1 , p n-1 u n , p n u n+1 , p n+1 p| n n-1 p| n+1 n x t T n-1 ,Y n-1 , ρ n-1 T n ,Y n , ρ n T n+1 , Y n+1 , ρ n+1 n -1 n n + 1 f 0 f 1 q n n-1 q n+1 n u n-1 2 u n+ 1 2 p| n-1 2 n-3 2 p| n+ 1 2 n-1 2 p n+1 p n p n-1

RESULTS

The present scheme is verified and validated along different verification and validation cases. Note that for all test cases, the time step is constant.

0-D Pressure cooker like cell

A single computational [1 × 1 × 1] m 3 cell with no fluxes is heated by its lateral faces while others are adiabatic. First, the scheme ability to account pressure variations is verified while conserving the mass of the system by imposing a Neumann condition on the walls. Then, its time convergence rate is verified through a Dirichlet condition. The exact temperature expressions for the Neumann (N) and Dirichlet (D) case read:

T N (t) = T 0 1 + t τ 1 , τ 1 = Ω ρ 0 c v T 0 S Q i . T D (t) = (T 0 -T w ) e -t τ 2 + T w , τ 2 = Ω ρ 0 c v S h .
T 0 = 300 K and ρ 0 = 1.177 kg m -3 are the initial temperature and density. Ω is the total cell volume, c v its heat capacity, and S the total heated wall surface. T w = 313 K denotes the wall temperature and the interior air heat transfer coefficient h = 30 W m -2 K -1 remains constant. Simulations are performed with a time step ∆t = 1 s and one inner iteration. For the Neumann case, the pressure, temperature and density L 2 error norm related to the exact solution are studied over time and compared to a simulation using a Poisson equation to correct the pressure (i.e by considering ρ n+1,k = ρ k in Eq. ( 3)). Figure 2, left, shows that the pressure variation is well reproduced by the present scheme and the mass is well conserved. The simulation using the Poisson equation led to more important errors. The increasing behaviour of errors are to be seen as an accumulation of truncation error in the range of the solver precision (10 -8 ). Figure 2, right, shows the pressure L 2 error at t = 0.8τ s for the Dirichlet boundary condition study. It is verified that the first and second order time convergence rate are obtained when θ = 1 and 1/2, respectively.

Sod shock case

The Sod [START_REF] Sod | A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws[END_REF] 1-D Riemann problem is considered to test the scheme accuracy related to singular solutions. Simulations are performed on a tube of length L = 400 m with 2 m × 800 cells, 0 ≤ m ≤ 5. An interface at the middle of the tube is used to initialise the right (R) and left (L) states as following: ρ L = 1.0 kg m -3 , ρ R = 0.125 kg m -3 , p L = 100 000 Pa and p R = 10 000 Pa. Walls are considered as symmetries except the two end faces of the tube, set as outlets. While refining the mesh at two CFL numbers, the variables L 1 error norm at t = 0.3 s are studied for both θ values. For each time step, three inner iterations are performed. Figure 3 shows, for a 3200 mesh and θ = 1, the different fields at t = 0.3 s. The exact solution, in black lines, is accurately reproduced by the simulation using the source term based on the kinetic energy equation (red lines). Even if the density, pressure and velocity fields are close to the exact solution, the temperature plateau is not achieved without the source term (blue lines). Figure 4 shows for each configuration the variables L 1 error norm related to the exact solution. It is well known that even even second order schemes can not achieve a second order accuracy for irregular solutions. When reducing the cells distance, the time convergence rates tend to 0.5 for the density and temperature and to 1.0 for the pressure and velocity, which are in agreement with these available in the compressible solvers literature [START_REF] Herbin | A cell-centred pressure-correction scheme for the compressible Euler equations[END_REF] [START_REF] Colas | An implicit integral formulation to model inviscid fluid flows in obstructed media[END_REF]. Lastly, the presented new CFL like conditions (Eq.4), represented in Figure 5, are shown to be less restraining than the speed based CFL condition.

Natural convection case

A heated cavity numerical benchmark case presented in [START_REF] Le Quéré | Modelling of natural convection flows with large temperature differences: a benchmark problem for low Mach number solvers. Part 1. Reference solutions[END_REF] is studied. The system of characteristic length L = 0.4603 m presents a steady flow driven by buoyancy effects at Ra = 10 6 . The side lateral walls difference of temperature (T h = 960 K and T c = 240 K) is responsible for the cavity motion and the fluid properties are constant. Other walls are considered adiabatic. A 740 × 740 cells mesh is used values. λ 0 = µ 0 γR/((1-γ)Pr), with Pr = 0.71, µ 0 = 1.68×10 -5 kg m -1 s -1 and R = 287 J kg -1 K -1 . The time step is fixed to ∆t = 0.0025 s and three inner iterations are performed. To accelerate the convergence to the steady state, the cavity temperature is initialised as T 0 = 600 K. Figure 6, left, shows that the steady state is reached after 10 s of simulation. Moreover, the normalised density stays constant and equal to the unity, proving the scheme's system mass conservation. The mean pressure reaches as well its reference value of 0.856. Furthermore, the Nusselt number profiles are in agreement to the reference (Figure 6, right), which is reflected by the correct reproduction of the centrelines velocity profiles (Figure 7).

CONCLUSIONS

A staggered time scheme has been proposed, validated and verified against analytical and benchmark problems. A numerical analysis has provided new conditions to insure the positivity of the thermodynamic variables. The latter were shown to be less constraining than the material CFL condition. Further- more, the presented scheme had good accuracy on accounting pressure variations and reproducing shock and natural convection solutions. An extension of this time scheme to moist-air is in development and further work shall focus on applying the presented scheme to more complex industrial studies.

Figure 1 :

 1 Figure 1: Time discretized variables localisation. Left frame : θ = 1. Right frame : θ = 1/2.

Figure 2 :

 2 Figure 2: Left: L 2 error norm of the density, temperature and pressure for the Neumann boundary condition case. (-) Poisson scheme. (-) presented scheme. Right: Dirichlet boundary condition case. L2 error for the pressure at t = 0.8τ 2 for a first-order implicit Euler (θ = 1) and Crank Nicolson (θ = 1/2) time scheme.

Figure 3 :

 3 Figure 3: Sod case results at t = 0.3 s for θ = 1 using a 3200 cells mesh. (-) exact solution (--) simulation using the source term Γ u 2 /2 (--) simulation without the source term Γ u 2 /2 .

  CFL = 0.04, θ = 1/2.

  CFL = 1.0, θ = 1.

  CFL = 1.0, θ = 1/2.

Figure 4 :

 4 Figure 4: L 1 error convergence for the SOD case for two values of CFL and θ.

Figure 5 :

 5 Figure 5: CFL x numbers for the SOD case for two different velocity based CFL simulations with θ = 1. Left frame: CFL = 0.04. Right frame: CFL = 1.0.

Figure 6 :

 6 Figure 6: Left: relative domain-averaged quantities ψ/ψ 0 over time. Right: Nusselt number profiles for the hot and cold walls compared to the reference [14].

Figure 7 :

 7 Figure7: Distribution of the velocity on the centrelines of the convective cavity compared to the reference[START_REF] Darbandi | Numerical simulation of thermobuoyant flow with large temperature variation[END_REF].

17th UK Heat Transfer Conference (UKHTC2021) 4-6 April 2022, Manchester, UK