Dolachai Boniface

Cécile Cottin-Bizonne

François Detcheverry

Christophe Ybert
email: christophe.ybert@univ-lyon1.fr

Role of Marangoni forces in the velocity of symmetric interfacial swimmers

Interfacial swimmers are objects that self-propel at an interface by autonomously generating a gradient of surface tension, often through the continuous release of a surfactant. While the case of asymmetric swimmers has long been studied, experiments have shown that spontaneous motion is also possible for symmetric swimmers. The basic mechanism of symmetry-breaking is qualitatively well-established but one key aspect of the phenomenon that has proved particularly difficult to elucidate is the role of Marangoni effects in the self-propulsion. We address this question by numerical methods, which can fully handle the complex interplay between swimmer motion, fluid flow, surfactant distribution, and Marangoni stresses. Our swimmer is a disk releasing a soluble surfactant in a deep-layer fluid. We investigate how the swimming velocity, represented by a Péclet number Pe * depends on its characteristics, as encapsulated in the Marangoni number M. We analyze the properties of the swimming diagram Pe * (M) and compare with approximate models to understand their origin. We find that the low-Pe * regime exhibits a bistability region: spontaneous swimming involves a threshold Marangoni number, a discontinuity in velocity and possibly hysteresis. Those features are present only for a full description of the problem and reveal the subtle but key role of Marangoni flows. The large-Pe * regime features a robust asymptotic scaling law Pe * ∼ M α , whose exponent α 0.72 is close to the 3/4 value predicted by a simplified model, indicating a much weaker influence of Marangoni flows. While our results were obtained assuming a point-source swimmer in the Stokes flow regime, we show that the picture remains very similar when considering a spatially extended source size, finite Reynolds number, or a fixed concentration swimmer. We discuss our findings in relation to experiments.

I. INTRODUCTION

A nonuniform surface tension generates at the liquid-gaz interface a tangential stress, the consequences of which are known as "Marangoni effects" [START_REF] Scriven | The Marangoni effects[END_REF][START_REF] Edwards | Interfacial Transport Processes and Rheology[END_REF]. Be it intentional or not-the contamination of water surfaces is often unavoidable in practice-the presence of surface-active species is the most pervasive cause for change in surface tension but local variation in temperature, solute content or liquid composition may also be at play. Given these various possible origins, Marangoni effects are ubiquitous. They reduce the velocity of a rising bubble [START_REF] Ybert | Ascending air bubbles in protein solutions[END_REF][START_REF] Takagi | Surfactant effects on bubbly flows[END_REF], damp capillary and ocean waves [START_REF] Lucassen-Reynders | Properties of capillary waves[END_REF][START_REF] Alpers | The damping of ocean waves by surface films: A new look at an old problem[END_REF], slow down the drainage of soap films and foams [START_REF] Breward | The drainage of a foam lamella[END_REF], thicken the film deposited on a dipped solid [START_REF] Quéré | Fluid coating on a fiber[END_REF], and thwart the drag reduction of superhydrophobic surfaces [START_REF] Peaudecerf | Traces of surfactants can severely limit the drag reduction of superhydrophobic surfaces[END_REF]. The common theme here is that any tendency to make the surfactant distribution inhomogeneous is counteracted. Besides drop break up or coalescence [START_REF] Eggleton | Tip Streaming from a Drop in the Presence of Surfactants[END_REF][START_REF] Dai | The mechanism of surfactant effects on drop coalescence[END_REF][START_REF] Liu | Coalescence of Bubbles with Mobile Interfaces in Water[END_REF] and tears of wine [START_REF] Venerus | Tears of wine: New insights on an old phenomenon[END_REF], the relevance of Marangoni phenomena also extends to the biological world, in bacterial colonies [START_REF] Trinschek | Modelling of surfactant-driven front instabilities in spreading bacterial colonies[END_REF] or surfactant-producing marine phytoplankton [START_REF] Botte | Numerical modelling of the Marangoni effects induced by plankton-generated surfactants[END_REF], and to technological applications such as desalination devices [START_REF] Morciano | Solar passive distiller with high productivity and Marangoni effect-driven salt rejection[END_REF]. In fact, in the many situations where the role of interfaces can not be neglected, Marangoni effects are potentially essential. It is therefore no surprise that they have been under investigation for more than two centuries.

Within the rich history of Marangoni phenomena [START_REF] Scriven | The Marangoni effects[END_REF], one particularly old strand of research is the one devoted to interfacial swimmers. Also called Marangoni surfers, these objects self-propel at a liquid-air or liquid-liquid interface by autonomously generating a gradient in surface tension. First observations with small camphor scrapings moving spontaneously on water surface date back to the 18th Century [17,[START_REF] Tomlinson | On the motions of camphor on the surface of water[END_REF] and as noted by Lord Rayleigh in 1890, the phenomenon, before being understood, "had puzzled several generations of inquirers" [START_REF] Rayleigh | Measurements of the amount of oil necessary in order to check the motions of camphor upon water[END_REF]. More than one century later, the interest in Marangoni swimmers was revived, this time focusing on the nonlinear and collective phenomena, such as synchronization, pattern formation and self-organization [START_REF] Suematsu | Synchronized intermittent motion induced by the interaction between camphor disks[END_REF][START_REF] Suematsu | Evolution of self-propelled objects: From the viewpoint of nonlinear science[END_REF][START_REF]Self-organized Motion: Physicochemical Design based on Nonlinear Dynamics[END_REF][START_REF] Nakata | Physicochemical design and analysis of self-propelled objects that are characteristically sensitive to environments[END_REF][START_REF] Satterwhite-Warden | Co-operative motion of multiple benzoquinone disks at the air-water interface[END_REF][START_REF] Soh | Swarming in shallow waters[END_REF][START_REF] Grzybowski | Self-assembly: From crystals to cells[END_REF]. The advent of active matter [START_REF] Bechinger | Active Brownian particles in complex and crowded environments[END_REF] gave further impetus to the study of interfacial swimmers and they are by now under intense scrutiny [START_REF]Self-organized Motion: Physicochemical Design based on Nonlinear Dynamics[END_REF][START_REF] Fei | Active colloidal particles at fluid-fluid interfaces[END_REF][START_REF] Grosjean | Surface swimmers, harnessing the interface to self-propel[END_REF]. The recent demonstration of swimming velocity reaching 10 4 body length per second show how efficient Marangoni propulsion can be [START_REF] Dietrich | Microscale Marangoni Surfers[END_REF]. While autonomous motion was first investigated for asymmetric swimmers such as camphor boats [START_REF] Lauga | Viscous Marangoni propulsion[END_REF][START_REF] Würger | Thermally driven Marangoni surfers[END_REF], it was later realized that this is not a necessary condition. Spontaneous swimming is also possible for systems such as camphor disks, that are perfectly symmetric, and thus release surfactant or heat in an isotropic manner. The underlying symmetry-breaking mechanism was first explored around the thermal or solutal hydrodynamic instabilities of active drops [START_REF] Ye Rednikov | On the development of translational subcritical Marangoni instability for a drop with uniform internal heat generation[END_REF][START_REF] Ye Rednikov | Drop motion with surfactant transfer in a homogeneous medium[END_REF][START_REF] Velarde | Drop motions and interfacial instability[END_REF] and later in the context of autophoretic colloids [START_REF] Michelin | Spontaneous autophoretic motion of isotropic particles[END_REF]. Beyond the physics of propulsion, ongoing work has also started to explore possible applications: interfacial swimmers may be designed as biocompatible and multifunctional [START_REF] Pena-Francesch | Multifunctional and biodegradable self-propelled protein motors[END_REF], and exploited for energy conversion and power generation [START_REF] Okawa | Surface tension mediated conversion of light to work[END_REF][START_REF] Frenkel | Mini-generator of electrical power exploiting the Marangoni flow inspired self-propulsion[END_REF], or as a building block for more advanced devices such as beating filaments [START_REF] Tiwari | Periodic oscillations in a string of camphor infused disks[END_REF].

There is in the motion of interfacial swimmers one phenomenon whose coupling with other processes has made it particularly difficult to account for: Marangoni flows. As a consequence, modeling efforts have often avoided the hurdle: they focused on the surfactant distribution at the swimmer edge and on the resulting driving force, but neglected the flow induced all over the surrounding fluid by a nonuniform surfactant distribution. The body of work relying on such a simplification is significant, as seen for instance in Refs. [START_REF] Nagayama | A theoretical and experimental study on the unidirectional motion of a camphor disk[END_REF][START_REF] Mikhailov | Self-motion in physicochemical systems far from thermal equilibrium[END_REF][START_REF] Iida | Theoretical study on the translation and rotation of an elliptic camphor particle[END_REF][START_REF] Kitahata | Spontaneous motion of an elliptic camphor particle[END_REF][START_REF] Koyano | Oscillatory motion of a camphor grain in a one-dimensional finite region[END_REF][START_REF] Shimokawa | Power law observed in the motion of an asymmetric camphor boat under viscous conditions[END_REF] and literature cited therein. We also adopted this approach in our recent study of symmetric Marangoni swimmers [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF]. Our point-source model entirely discarded Marangoni flows but the prediction for velocity and its size dependence seemed to capture surprisingly well the experimental data obtained on camphor disks. Such an agreement was unexpected in at least two ways. First, the existence of Marangoni flows around camphor swimmers has been evidenced since the 19th Century [START_REF] Tomlinson | On the motions of camphor on the surface of water[END_REF][START_REF] Tomlinson | On the motions of camphor on the surface of water[END_REF] and neglecting them is far from straightforward owing to a magnitude comparable with the swimmer velocity. Second, one exact calculation with an asymmetric swimmer [START_REF] Lauga | Viscous Marangoni propulsion[END_REF] indicates that for the system studied, Marangoni flows may significantly decrease the propulsion velocity. There is certainly no obvious reason as why one could simply dispense with Marangoni flows in the modeling of symmetric Marangoni swimmers, but doing so leads to predictions that were nonetheless reasonable. Our motivation for this work was to understand this paradox.

The problem of Marangoni flows in interfacial swimmers is challenging because of the intricate couplings between fluid motion, surfactant advection and diffusion and swimmer displacement. There are, however, less involved situations where Marangoni flows could be analyzed and understood. The first is the transient spreading dynamics of surfactant punctually deposited at the interface [START_REF] Jensen | The spreading of insoluble surfactant at the free surface of a deep fluid layer[END_REF][START_REF] Bickel | Spreading dynamics of reactive surfactants driven by Marangoni convection[END_REF]. The second involves a fixed, steady release of surfactant, a configuration which has received much attention in recent years [START_REF] Roché | Marangoni Flow of Soluble Amphiphiles[END_REF][START_REF] Bandi | Hydrodynamic Signatures of Stationary Marangoni-Driven Surfactant Transport[END_REF][START_REF] Kitahata | Effective diffusion coefficient including the Marangoni effect[END_REF][START_REF] Arangalage | Dual Marangoni effects and detection of traces of surfactants[END_REF][START_REF] Benouaguef | Solutocapillary Marangoni flow induced in a waterbody by a solute source[END_REF]. The flow in the vicinity of the source exhibits some generic features [START_REF] Le Roux | Soluble surfactant spreading: How the amphiphilicity sets the Marangoni hydrodynamics[END_REF] whereas at large distance, the velocity field decays as a power law whose exponent depends on the type of surfactant [START_REF] Mandre | Axisymmetric spreading of a surfactant driven by self-imposed Marangoni stress under simplified transport[END_REF]. Besides, experiments with fixed symmetric and asymmetric interfacial swimmers-camphor disks and boats, respectively-hint at differences in the resulting Marangoni flow [START_REF] Matsuda | Acceleration or deceleration of self-motion by the Marangoni effect[END_REF]. Moving closer to the interfacial swimmer configuration, the case of a mobile source was addressed theoretically in Ref. [START_REF] Williams | Two-dimensional nonlinear advection-diffusion in a model of surfactant spreading on a thin liquid film[END_REF]. The predicted scaling laws, however, are not directly 104006-2 applicable because the shallow water assumption does not hold for most experiments. Finally, the first characterization of Marangoni flows around a swimmer in motion were obtained recently using experimental and numerical methods. They reveal a divergent flow behind the swimmer [START_REF] Sur | Translational and rotational motion of disk-shaped Marangoni surfers[END_REF] and how the direction of motion may be reversed on very shallow films [START_REF] Kang | Forward, reverse, and no motion of Marangoni surfers under confinement[END_REF]. To date though, only asymmetric swimmers have been considered.

The role of Marangoni flows in the self-propulsion of symmetric interfacial swimmers thus remains an open question. The goal of this work is to fill this gap. Using numerical methods, the only approach that can handle the full complexity of the problem, we focus on the swimming velocity, the key experimental observable. By comparing swimming with and without Marangoni flows, we assess their influence on self-propulsion. Depending on the regime considered, we find that Marangoni flows play a very different role. They are essential in the vicinity of the threshold to decide how symmetry breaking occurs. However, in the regime of large velocity, their influence is surprisingly modest. We discuss the implications of these findings in connection with experimental results.

While the writing of this work was being finalized, two studies by Ender and coworkers appeared that address a very similar problem [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF]. As will be discussed below, our results are essentially consistent with their findings on overlapping topics, but differ with their work on a number of elements. Let us point out the main differences. (i) Our swimmer is not a hemisphere but an infinitely thin disk. This geometry is by far the most widespread [START_REF] Satterwhite-Warden | Co-operative motion of multiple benzoquinone disks at the air-water interface[END_REF][START_REF] Nagayama | A theoretical and experimental study on the unidirectional motion of a camphor disk[END_REF][START_REF] Shimokawa | Power law observed in the motion of an asymmetric camphor boat under viscous conditions[END_REF][START_REF] Matsuda | Acceleration or deceleration of self-motion by the Marangoni effect[END_REF][START_REF] Nakata | Self-motion of a camphanic acid disk on water with different types of surfactants[END_REF][START_REF] Soh | Dynamic self-assembly in ensembles of camphor boats[END_REF][START_REF] Suematsu | Self-motion of a Benzoquinone disk coupled with a redox reaction[END_REF][START_REF] Nakata | Mode change in the self-motion of a benzoquinone disk coupled with a NADPH system[END_REF][START_REF] Koyano | Relationship between the size of a camphor-driven rotor and its angular velocity[END_REF][START_REF] Akella | Dynamics of a camphoric acid boat at the air-water interface[END_REF], presumably because it is most convenient to craft efficient-high velocity-swimmers. As we shall see, this seemingly minor difference can modify qualitatively some aspects of the swimming response. (ii) Our study focuses mostly on steady creeping flow, as Refs. [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF], but also briefly considers finite Reynolds number and transient regimes. (iii) We examine in detail the vicinity of the symmetry breaking mechanism, where the influence of Marangoni flows is the strongest. Our study is thus largely complementary to the work of Ender et al. and should help to provide a complete picture.

The remainder of this article is organized as follows. Section II presents the modeling assumptions and governing equations, as well as the implementation of the numerical method. Section III introduces the force balance that characterizes the swimming point and the resulting swimming diagram. We examine in turn the regime of low velocity, i.e., close to the instability where the motions sets in, and the asymptotic regime of large velocity. Several extensions-including finite Reynolds number-are discussed in Sec. IV. A summary and some perspectives are given in Sec. V. Finally, several technical points-among which some experimental results-are gathered in the Appendices.

II. MODEL AND METHOD

A. Physical system and equations

Governing equations

Our swimmer is a disk releasing a soluble surfactant and placed at the interface between air and a liquid that occupies the semi-infinite space below, as illustrated in Fig. 1. We note a the radius of the disk, D its surface, ∂D its perimeter, and S the air-liquid interface. The experimentally observed spontaneous swimming corresponds to such a disk moving in the laboratory reference frame at a stationary velocity U * , set by the system geometry and physicochemistry. To gain more insights, we consider here a more general situation where the disk moves at an arbitrary instantaneous velocity U , which unless specified is assumed constant. For U = U * , this requires that an external force F op is applied by an operator so that the total force vanishes in steady state F = 0. In the following, this problem is treated in the disk reference frame with origin O at the disk center; positions are noted r = (x, y, z) and the fluid flow far away from the swimmer is fixed to -U e x (see Fig. 1). We assume a Newtonian incompressible liquid, whose velocity field u(r, t ) obeys the Navier-Stokes equation

∇ • u = 0, ρ ∂u ∂t + u • ∇u = ∇ • σ = -∇ p + η u. ( 1
)
104006-3 Here ρ is the fluid density, η its dynamic viscosity, p the pressure, and the Laplacian. The stress tensor σ is given by σ = -pI + τ, with I the identity tensor and τ = η(∇u + ∇u T ) the viscous contribution. The concentration field c(r, t ) of the soluble surfactant released by the swimmer is governed by advection and diffusion,

∂c ∂t + ∇ • (uc) = D c, (2)
with D the diffusion coefficient. These equations have to be complemented by boundary conditions (BCs). For the concentration field under the swimmer, we have

-D∂ z c| D = J (r), (3)
where J (r) is the flux of surfactant at position r for a swimmer at the origin. Unless mentioned otherwise, we will consider a simplified configuration where the release occurs only at the swimmer center, that is J (r) = J δ s (r), with J the total release rate of solute and δ s the two-dimensional Dirac delta function. In addition, vanishing of the far-field concentration c(r, t )| r→∞ = 0 and of the normal solute flux on the free-surface -D∂ z c| S = 0 is imposed. While solute transport is coupled to the flow equation through the advective term of Eq. ( 2), reciprocal coupling with the flow field occurs through BCs. Indeed, the presence of the surfactant modifies the surface tension γ through

γ = γ 0 -κc s , ( 4
)
with γ 0 the surface tension of a bare interface, κ a constant coefficient and c s ≡ c(z → 0) the volume concentration of surfactant in the vicinity of the interface. Spatial heterogeneities in concentration thus induce gradients in surface tension which result in a Marangoni stress at the liquid-gas interface S,

σ • n| S = ∇ s γ , ( 5
)
where n is the unit vector normal to the interface and pointing outside the liquid and ∇ s = I s • ∇ is the surface gradient, with I s = In ⊗ n the surface identity tensor. This BC is complemented 104006-4

again by the usual conditions on fluid velocity: no-slip below the swimmer u(r, t )| D = 0, uniform far-field u(r, t )| r→∞ = -U e x and impenetrability of the free surface u(r, t ) • n| S = 0.

Simplifying assumptions

In writing the equations above and in using them below, we make a number of simplifying approximations that we now discuss in detail.

(i) A flat interface. In principle, a free surface can deform in response to the hydrodynamic pressure so that the most general boundary condition for stress there [Eq. ( 5)] would include in the right-hand side a Laplace pressure term -γ Cn. Here C is the local mean curvature of the surface whose shape and position adjust so that the liquid pressure matches the Laplace pressure jump.

Assuming an undeformed flat interface C = 0 implies that the flow pressure cannot overcome surface tension effects. This corresponds to a small capillary number Ca = ηU/γ in the viscous regime and to a small Weber number We = ρU 2 a/γ in the inertial regime. Experimental values are typically U = 7 mms -1 for a = 4 mm so that Ca 10 -3 and We 0.3. The flat interface is thus a reasonable approximation. Consistent with this assumption, the swimmer lower face is chosen to coincide with the interface: the thickness, therefore, plays no role and can formally be taken as zero.

(ii) Fast equilibration of surfactant coverage. To describe the surfactant distribution in the general case, it would be necessary to introduce a surface concentration governed by the advection-diffusion ∂ t + ∇ s • (u ) = D s s + J ex , where s is the surface Laplacian, D s the surface diffusion coefficient, and J ex is the net flux of surfactant coming from the bulk liquid to the surface. In its simplest form, such an exchange flux J ex could be written as

J ex = r a c s ( sat - ) -r d (c sat -c s )
, where sat and c sat are surface and volume concentrations at saturation, and r a and r d are adsorption and desorption rates.

Considering distinct dynamics for the surfactant at the surface and in the bulk would add a level of sophistication to a problem that is already quite complex. For the sake of simplicity, we assume infinitely fast exchange dynamics between the bulk and surface surfactant. In this limit of high sorption rates, equilibration is instantaneous between the volume concentration near the surface c s and the surface concentration ; the former thus fixes the latter. The surface tension may then be written as γ (c s ), i.e., a function of c s only as taken in Eq. [START_REF] Takagi | Surfactant effects on bubbly flows[END_REF].

Three remarks are related to the surfactant dynamics. First, we do no consider evaporation of surfactant into the air. As discussed in Appendix E, measurement of camphor evaporation gives a time constant on the order of hour, and a Biot number much below unity, indicating that its role in self-propulsion is negligible. This often-mentioned peculiarity of camphor is, however, likely to contribute to the long term rejuvenation of the system, and thus to the remarkable longevity of the swimmer activity. Second, our set of assumptions makes our model directly applicable to a swimmer releasing heat rather than a soluble surfactant. Third, another simple limiting case for the surfactant dynamics is to consider a strictly insoluble surfactant [START_REF] Lauga | Viscous Marangoni propulsion[END_REF][START_REF] Mandre | Axisymmetric spreading of a surfactant driven by self-imposed Marangoni stress under simplified transport[END_REF]. Yet, real systems usually require a finite solubility to treat the surfactant release from the core of the swimmer and this limit is therefore difficult to address consistently. This case will not be considered in the following.

(iii) Linear dependence in γ (c s ). Equation [START_REF] Takagi | Surfactant effects on bubbly flows[END_REF] indicates that the surface tension decreases linearly with the concentration c s . Such an assumption is common in the literature and may be justified by starting from a Langmuir isotherm γ = γ 0 + RT sat ln(1 -/ sat ), with R the gas constant and T the temperature. Assuming sat and c s c sat , one finds κ = RT K H , with the Henri length K H = r a sat /r d c sat . This assumption is valid only if the surface and bulk content remain well below their saturation value. To the best of our knowledge, only one experiment reported the fore-aft surface tension asymmetry around a camphor boat, with typical amplitude γ 1 mNm -1 [START_REF] Karasawa | Simultaneous measurement of surface tension and its gradient around moving camphor boat on water surface[END_REF]. This suggests that even outside Henri's regime, a linearized approach around a working point should be justified.

(iv) Imposed flux and punctual release. Equation (3) indicates that the swimmer imposes a fixed flux of surfactant. This corresponds approximately to the experimental situation fulfilled with camphor swimmers of Ref. [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF]. The assumption of punctual release is chosen to allow a direct 104006-5 comparison with a point-source model that we investigated in a previous work [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF], but it will be relaxed in Sec. IV A.

(v) Steady state. Unless mentioned otherwise, the fluid flow is assumed stationary. Unsteady swimming with intermittent motion or time-varying velocity, as sometimes observed experimentally [START_REF] Suematsu | Synchronized intermittent motion induced by the interaction between camphor disks[END_REF][START_REF] Suematsu | Evolution of self-propelled objects: From the viewpoint of nonlinear science[END_REF][START_REF] Akella | Dynamics of a camphoric acid boat at the air-water interface[END_REF], is therefore not explored here. Nonetheless, we will consider briefly in Sec. III B 1 a transient regime, the relaxation toward steady state.

Dimensionless numbers

To make the governing equations dimensionless a classical choice is to set the swimmer radius a as unit length, and the velocity U as unit velocity. This introduces well-known numbers

Re = U a ν , Pe = U a D , ( 6
)
where ν = η/ρ is the kinematic viscosity of the fluid. As usual, the Reynolds number Re indicates the relative magnitude of inertial and viscous effects in the Navier-Stokes equation. The Péclet number Pe compares the surfactant transport by flow advection to transport by diffusion.

In the following, we shall, however, adopt a different nondimensionalizing based on another velocity scale U M . This scale corresponds to the typical Marangoni velocity associated with surface tension gradients U M = κc b /η, with the concentration scale set by c b = J /Da [71]. Finally, with unit stress chosen as ηU M /a, the Eqs. ( 1) and (2) become in steady state

M Sc ũ • ∇ ũ = -∇ p + ˜ ũ, M ũ • ∇ c = ˜ c, (7)
where all quantities are dimensionless, as indicated by the tilde. These equations are complemented by dimensionless BCs for the concentration

c(r)| r→∞ = 0, -∂ z c| D = δs (r), -∂ z c| S = 0, (8)
and for the velocity

ũ(r)| r→∞ = - Pe M e x , ũ(r)| D = 0, τ • n| S = -∇s c, ũ(r) • n| S = 0. (9)
In the above equations, we have introduced two additional dimensionless numbers

Sc ≡ ν D = Pe Re , M = κJ ηD 2 . ( 10
)
The Schmidt number Sc is the ratio of momentum and mass diffusivities in the liquid. The Marangoni number M corresponds to a Péclet number U M a/D based on the Marangoni velocity and compares advective transport by such flows with diffusion. It can also be viewed as a balance between chemical activity and damping effects. The chemical activity of the swimmer combines the solute release rate J and the Marangoni susceptibility κ which quantifies its ability to convert concentration gradients in mechanical stress. It is then moderated by damping effects from transport properties in the fluid, with viscosity reducing the achievable velocities for given applied stresses, and diffusivity smearing out heterogeneities. As a cautionary remark, we note that naming conventions for dimensionless groups are not uniform. The definition used here are also found in Refs. [START_REF] Kang | Forward, reverse, and no motion of Marangoni surfers under confinement[END_REF][START_REF] Vandadi | Reciprocal theorem for convective heat and mass transfer from a particle in Stokes and potential flows[END_REF][START_REF] Yariv | Phoretic self-propulsion at large Péclet numbers[END_REF] but in other works such as Refs. [START_REF] Michelin | Spontaneous autophoretic motion of isotropic particles[END_REF][START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF], our M number was called a Péclet number.

In most real systems, only the Marangoni number M can be controlled. Re and Pe both involve the spontaneous swimming velocity which is an outcome of the experiment, while Sc is set by the nature of the liquid and solute. Considering for concreteness the camphor swimmers of Ref. [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF], typical values of Re, Pe, and M are around 10 2 , 10 5 , and 10 10 , respectively. As can be seen, Re is not small and inertial effects are significant in experimental systems. In the following, we will, however, investigate in detail the limit of Stokes flow reached when Re → 0 (Sc → ∞). This 104006-6 particular case is worth considering on several grounds. First, a quantitative comparison can be made with analytical approaches as developed for such system with a point-source model [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF]. Second, analytical solution for the Stokes flow around a disk may be exploited, as will be shown below, and used to provide physical insight that can help to examine more realistic regimes. Third, this situation is easier to handle with numerical methods as it eliminates one source of nonlinearity. Finite values for the Reynolds (Schmidt) number will also be considered, but less extensively in Sec. IV B. In the numerical computations presented below, the swimmers configuration involve Pe and M values that can reach 10 4 and 10 7 . Though still smaller than experimental value by one or several orders of magnitude, they fall in the correct regime of values much higher than unity.

B. Numerical method

Even within our idealized description, the steady motion of a swimmer is a problem coupling several physical phenomena whose analytical solution appears out of reach. Accordingly, we resort to a numerical approach based on the finite element method (FEM). In practice, the software COMSOL, a multiphysics commercial platform, is used. In this section, we briefly present the details of our implementation and how the reliability of the results was checked.

Geometry and boundary conditions

The configuration and our conventions are shown in Fig. 1. Because the system is symmetric with respect to the y = 0 plane, only one half of the domain is considered. The simulation box is a rectangular cuboid with dimension L x /2 = L y = L z = L. The boundary closest to the swimmer center is thus at a distance L. Remembering that the unit length is the swimmer radius, we set L = 400, a large value so as to minimize finite size effects.

Boundary conditions (BCs) on the box depend on the face considered. For the bottom, front and side faces, the hydrodynamic BC is a fixed velocity -U e x , whereas the BC for surfactant transport is a vanishing concentration. The rear face serves as an outlet for the flow and advected surfactant [START_REF] Besides | the normal diffusive flux of surfactant is set to zero[END_REF]. On the top surface, the no-slip BC below the swimmer imposes a vanishing velocity there. At the flat free surface, the flow BC combines a vanishing normal velocity associated with non penetrability of the interface, with a tangential Marangoni stress given by Eq. ( 7) for τ. Regarding the surfactant, it is released from a point-source at the disk center-except in Sec. IV A where spatially extended source will be considered-and a no-flux "insulating" BC applies everywhere else at the top surface.

Mesh design

An essential step in using the FEM is the choice of a mesh, which always involves a compromise: a finer mesh usually offers better accuracy but the price to pay is a higher computation time. We have adapted our mesh to make it appropriate for the physics of the problem, by locally refining the elements only in those regions where strong gradients exist. More specifically, the mesh is a set of tetrahedrons partitioning space, whose characteristic size ξ and growth rate χ can be adjusted [75]. In practice, we defined three areas where ξ and χ are constrained to ensure sufficient accuracy. In the vicinity of the point source, we set ξ = 5 × 10 -4 . Along the swimmer contour that delineates a sharp border between the no-slip and Marangoni stress BCs, we choose ξ = 2.5 × 10 -3 . Finally, in a large neighborhood around the swimmer-a disk centered at the origin with radius 10-the mesh is constrained to a maximal size ξ m = 5 × 10 -2 and a maximal growth rate χ m = 1.1. The remaining space does not require particular attention and coarser elements with size up to ξ m = L/20 can be used without affecting the calculation. Nevertheless, a rather low growth rate χ m = 1.15 is used to ensure that the mesh close to the swimmer remains sufficiently fine. The resulting mesh usually involves several millions of elements and typically leads to computation times of several hours on a dedicated workstation. 104006-7

Finite-size effects and convergence in mesh

Two requisites must be met by our choice of box size and mesh. First, finite-size effects must be negligible, meaning that our box, though finite, can nevertheless faithfully mimic an infinitely large system. Second, convergence in mesh should be ensured, meaning that upon successive mesh refinement, the results reach a well-defined limit. That both requirements are satisfied is shown in Appendix F. We estimate that forces, the quantities of prime interest in what follows, are obtained with a few percent accuracy.

In contrast to experiments, the swimmer velocity U is imposed at an arbitrary value. For any given choice of parameters, the resulting outcome of a simulation is thus the steady concentration and velocity fields together with the external operator force F op required to satisfy the force balance condition F = 0. The velocity U * that would be adopted by a free swimmer is reached when the external force exerted by the operator vanishes. We now discuss in detail such force balance and the different contributions it involves.

III. FORCE BALANCE AND SWIMMING BEHAVIOR

A. Contributions to force balance

Flow and surface tension contributions

In addition to the external operator force, the force balance in steady state is classically formulated using the capillary force F c generated on the swimmer perimeter and the viscous force F v exerted by the fluid underneath, namely,

F ≡ F c + F v +F op = 0, F c ≡ ∂D γ n ∂D dl, F v ≡ D n D • τdS. ( 11
)
Here F is the total force on the swimmer, n ∂D is the unit vector normal to the contact line along the interface and n D the unit vector normal to the swimmer surface [76]. For better clarity, the triple bar is used for equality that defines a symbol. The influence of surface tension is not limited to the capillary force but also arises indirectly through the fluid: Marangoni stress at the free surface generates a flow that eventually contributes to the force on the disk via viscous interactions. Surface tension effects thus encompass a direct capillary contribution and a flow-mediated Marangoni contribution. To get an in-depth understanding of how swimming spontaneously arises from chemical release, it is meaningful to isolate all contributions from surface tension. This is possible in the Stokes regime by using the Lorentz reciprocal theorem, which has proved very valuable in the context of self-propelled particles in general [START_REF] Masoud | The reciprocal theorem in fluid dynamics and transport phenomena[END_REF][START_REF] Elfring | A note on the reciprocal theorem for the swimming of simple bodies[END_REF] and interfacial swimmers in particular [START_REF] Masoud | A reciprocal theorem for Marangoni propulsion[END_REF][START_REF] Gidituri | Dynamics of a fully wetted Marangoni surfer at the fluid-fluid interface[END_REF].

Applying the Lorentz theorem requires the introduction of a dual problem, which differs from the initial problem only in the choice of the boundary conditions. This is generally a simpler situation where analytical results are available. Here, the dual problem is that of an inert disk moving atop the interface at a velocity U : in the swimmer reference frame considered this corresponds to a fixed disk submitted to a flow -U e x at infinity. It thus differs from the initial problem by the fact that the interface is free of surfactant. The reciprocal theorem then allows rewriting the viscous force as

F v = Fv + F vM , Fv ≡ -CηaU e x , F vM ≡ S û U • ∇ s γ dS e x . ( 12
)
Here any quantity of the dual problem is indicated with a hatted symbol Q. As detailed in Appendix B 2 and Eq. (B9), the associated flow field û for such a disk moving edgewise [81] is known analytically [START_REF] Happel | Low Reynolds Number Hydrodynamics[END_REF], as is the Stokes drag Fv for which C = 16/3. We are now in a position to reformulate the force balance with a term, called the Marangoni force F M , which gathers both capillary and viscous contributions from surface tension. Using Eq. ( 12) together with the divergence theorem, one can write

F M ≡ F c + F vM = - S γ ∇ s • û U dS e x . ( 13
)
Given the integral expression for the Marangoni force F M , the power U e x • F M it would develop for a free swimmer moving at velocity U e x is reminiscent of the rate of change of surface energy:

Ėγ = - S γ ∇ s • udS. ( 14
)
As discussed in details in Appendix A, with ∇ s • u the local rate of area creation, the integrand is the two-dimensional analog to the power of pressure forces in compressible systems. In the free swimming regime (F op = 0), this term quantifies the amount of-chemical-power available to the system. The reason why û appears in Eq. ( 13) in place of the total flow u is that some of the chemical power is lost in extra bulk viscous dissipation : only a fraction of the initial budget is thus available for the swimmer.

Forces evolution and swimming criterion

Magnitude and sign of forces. We first start by examining how surface tension effects evolve for imposed flow of arbitrary velocity, before focusing on equilibrium points where the swimmer is force-free, with force balance achieved without external operator (F op = 0). In the following, all the numerical results for forces are shown as dimensionless numbers. Starting from the motionless state, we can see in Fig. 2 that a finite advection Pe > 0 induces a capillary force F c which is motor, that is with a direction opposite to the flow. The magnitude increases with Pe until reaching a maximum around Pe 10, after which it decreases toward zero. The trend observed for the Marangoni force F M is similar.

Qualitatively, such a behavior can be rationalized as follows. Whereas the system in the rest state is fore-aft symmetric, advection by the fluid will skew the chemical cloud by bringing solute from the front to the rear. Such a convection generates a capillary force that will always be a driving contribution F c = F c e x , F c > 0, as visible from Eq. ( 11). This is the basis for the qualitative understanding of a symmetry breakdown in chemically or thermally active systems [START_REF] Michelin | Spontaneous autophoretic motion of isotropic particles[END_REF]. Regarding the Marangoni force, a similar picture applies. The surface divergence appearing in Eq. ( 13) can be written according to Eq. (B9) as

∇ s • û = -f (r) cos θ with f (r) = 4U/(3π ar √ r2 -1) 104006-9
being positive for all r > 1. Combined with Eq. ( 4) and with the skewed chemical cloud, this implies that in response to convection the Marangoni force is always motor F M > 0. If now we consider specifically the large Pe regime, then we note that while the surfactant is emitted by the swimmer with a constant flux J, it is carried away with an advection velocity U , setting a relevant concentration that scales as J/U . The local concentration is thus a decreasing function of Pe and so is the amplitude of surface tension variations. This directly translates into Marangoni and capillary forces both decreasing toward zero. While a proper treatment would require boundary layer considerations, this simple argument is sufficient to understand the observed trend in the limit of high Pe.

We finally comment on the sign of forces. As may be expected intuitively, the capillary and (total) Marangoni forces are always found to be driving contributions. However, the behavior of the viscous Marangoni contribution is more complex: as can be seen in Fig. 2, the viscous Marangoni force can be either motor or resistive in our system, depending on the strength of imposed convection. A resistive contribution is arguably the natural expectation from previous studies. In particular, for an asymmetric swimmer releasing an insoluble surfactant in the diffusive regime (Pe 1), the viscous Marangoni force was shown to be resistive at the swimming point [START_REF] Lauga | Viscous Marangoni propulsion[END_REF][START_REF] Masoud | A reciprocal theorem for Marangoni propulsion[END_REF]. Though our system is different-symmetric, with soluble surfactant and finite Pe-this observation still holds. As can be seen in Fig. 2, F vM becomes resistive for Pe = 2.4, and remains so at the swimming point (Pe * = 4.8). That flow-mediated Marangoni contribution could be motor was, however, recently identified for partly immersed swimmers [START_REF] Kang | Forward, reverse, and no motion of Marangoni surfers under confinement[END_REF][START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF][START_REF] Vandadi | Reciprocal theorem for convective heat and mass transfer from a particle in Stokes and potential flows[END_REF]. For asymmetric swimmers, the sign reversal occurs in the shallow water regimes and was ascribed to pressure-effects induced by depth-constrained flows. Because we consider only a deep-water configuration and the pressure field does not enter into the force balance [83], the reversal of the flow-mediated contribution must have a different origin. It is actually controlled by the strength of convection, as also observed independently for symmetric half-spheres [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF].

Swimming diagram. The condition for spontaneous motion is the existence of a stable swimming point, i.e., a stable equilibrium point with finite velocity. An equilibrium (or fixed) point is reached when the Marangoni force exactly balances the classical Stokes drag force, as specified in Eq. ( 13), in the absence of operator force. In dimensionless form and with our conventions, this force-free condition reads as

FM (Pe, M) = 16Pe 3M for Pe = Pe * . ( 15
)
As illustrated in Fig. 3, the balance of forces corresponds to the existence of an intersection point between the Marangoni force FM and a linear curve with a slope set by M. Spontaneous swimming also requires that the swimming point is stable, meaning that a perturbation induces a net restoring force:

e x • ∂ (F M + Fv ) ∂Pe Pe=Pe * < 0. ( 16
)
Graphically, this condition is fulfilled when the drag force line intersects the Marangoni force from below.

The number and nature of swimming points depend on the value of M, as illustrated in Fig. 3. For low M, the force-free condition is only satisfied with zero velocity, this is the motionless state. At higher M, both stable and unstable swimming points may be found. In the following, we focus primarily on the resulting swimming diagram. Because the influence of Marangoni flows is very different, we consider separately the regime of low Pe * where the instability sets in (Sec. III B), and the asymptotic regime of high-Pe * (Sec. III C). Before proceeding, we remind that all results within this Sec. III are obtained in the Stokes regime and with a point-source release of surfactant.

104006-10 FIG. 3. Dimensionless force balance condition and swimming points for different Marangoni numbers, as given by Eq. ( 15). Symbols show the numerical results for the dimensionless Marangoni force FM and dashed lines the right-hand side of Eq. ( 15), which is proportional to the Stokes drag. Stable and unstable swimming points are indicated with solid red and empty blue circles, respectively.

B. Bifurcation regime (Pe * small)

A swimming diagram with a bistable region

The swimming diagram in the region of low Pe, shown in Fig. 4, reveals that depending on the Marangoni number M, three distinct behaviors are possible. Below a critical value M c 250, only the motionless state exists, the disk is no swimmer. As M exceeds M c , the number of fixed points actually changes from 1 to 3 (5 if accounting for the two possible directions in swimming but for simplicity only one is considered hereafter, as in Fig. 4). The motionless state remains a stable fixed point [84], but an additional pair of swimming points, one stable and one unstable, have appeared, thus defining a region of bistability where swimming and rest states are both possible. Finally, at a second critical value M c 410, the unstable swimming point merges with the motionless state. Because the latter is now unstable, the only stable state is swimming at a finite velocity.

In line with previous studies on the propulsion of symmetric particles [START_REF] Michelin | Spontaneous autophoretic motion of isotropic particles[END_REF][START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF][START_REF] Yoshinaga | Drift instability in the motion of a fluid droplet with a chemically reactive surface driven by Marangoni flow[END_REF][START_REF] Morozov | Self-propulsion near the onset of Marangoni instability of deformable active droplets[END_REF], spontaneous swimming appears only above a threshold value M c (note, however, the case of a hemispherical swimmer [START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF], for which the threshold M c vanishes as the vessel size L diverges but the approach is only logarithmic M c ∼ 8/ ln L). To date, bifurcations in such interfacial swimming systems were reported as supercritical pitchfork [START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF][START_REF] Yoshinaga | Drift instability in the motion of a fluid droplet with a chemically reactive surface driven by Marangoni flow[END_REF], with a stable swimming state appearing when the motionless state loses its stability. Here we have some evidence of a different bifurcation diagram. The data suggests that the bifurcation is saddle-node at M c and subcritical pitchfork at M c . To the best of our knowledge, the resulting bistability region has not been observed before [87]. Interestingly, it echoes some of the behaviors predicted for active droplets in bulk systems for which a richer panel of transitions have been unveiled [START_REF] Ye Rednikov | On the development of translational subcritical Marangoni instability for a drop with uniform internal heat generation[END_REF][START_REF] Morozov | Self-propulsion near the onset of Marangoni instability of deformable active droplets[END_REF], including bistable regimes and finite amplitude instabilities [START_REF] Ye Rednikov | Drop motion with surfactant transfer in a homogeneous medium[END_REF][START_REF] Ye Rednikov | Active drops and drop motions due to nonequilibrium phenomena[END_REF] The swimming diagram identifies only the fixed point obtained for steady swimming but their stable or unstable nature should indicate what the spontaneous evolution of the system is when initiated from an arbitrary state. To confirm the expected behavior, we extended our numerical calculations to describe time-dependent swimming [90]. Figure 5 the velocity relaxes toward the closest accessible stable state. For M < M c , the motion gradually stops, in line with the absence of spontaneous swimming at small M. When M > M c , whatever the initial value, the velocity is seen to converge to the stable swimming point. The motionless state is now unstable and spontaneous symmetry breaking occurs for infinitesimal perturbations. Finally, in the bistable region of the swimming diagram, the final state depends on the location of Pe in with respect to the unstable fixed point Pe * (u) . If below, then the initial motion is damped until complete stopping; if above, then the stable swimming point is eventually reached.

The bistable region has at least two consequences. First, there is a range of swimming velocity that is forbidden, specifically Pe * ∈ [0, Pe * c ] with Pe * c 3. Since the Péclet number of spontaneous swimming can not be made arbitrarily small, one may not rely on the low-Pe assumption that was used in the past to attack the problem of asymmetric swimmers analytically [START_REF] Lauga | Viscous Marangoni propulsion[END_REF][START_REF] Würger | Thermally driven Marangoni surfers[END_REF]. The second consequence is the existence of hysteresis. Let us assume that M can be varied at will and that perturbations to the swimmer are negligible. Upon increasing M from zero, the disk will switch to motion at M c , while upon decreasing M from a large value, the motion will stop at M c . Both transitions are discontinuous with a finite jump in velocity. If perturbations are present, then the jump might take place anywhere in the [M c , M c ] interval. It must be said that varying the Marangoni number while keeping the same swimmer is a priori difficult in experiments. With fluid and surfactant fixed, one should change the swimmer size or release flux, which is in general impossible for chemical swimmers. We note, however, that our results equally apply to a thermal swimmer, where the heat released would locally lower the surface tension [START_REF] Dietrich | Microscale Marangoni Surfers[END_REF][START_REF] Würger | Thermally driven Marangoni surfers[END_REF][START_REF] Maggi | Micromotors with asymmetric shape that efficiently convert light into work by thermocapillary effects[END_REF]. Changing the power of light illuminating the disk might be a practical way to control the Marangoni number of a given swimmer. Finally, to finish with the implication of bistability, one may wonder whether it plays any role in the intermittent swimming behavior that was observed in several cases [START_REF] Suematsu | Evolution of self-propelled objects: From the viewpoint of nonlinear science[END_REF][START_REF] Akella | Dynamics of a camphoric acid boat at the air-water interface[END_REF].

Discussion: Origin of bistability

Our goal here is to better understand the bistability observed in the swimming diagram of symmetric Marangoni swimmers. To examine when this new feature appears, we consider simpler 104006-12 models, whose common starting point is Eq. ( 13) for the Marangoni force. With Fv given by the Stokes expression, there is in the force balance only one unknown quantity, the surface tension distribution γ involved in the integral expression for F M . We recall that for spontaneous swimming, the system is force-free F op = 0. Our approach is to evaluate this surface tension distribution-or equivalently the surfactant concentration field-by decoupling the transport problem from the Stokes equation for the flow. In the convection-diffusion equation for the surfactant [Eq. ( 7)], the true fluid velocity is thus replaced by an approximate flow field u app . The simplest possibility is to take a uniform flow field u app = -U e x , which, as shown in our recent work, can yield meaningful results [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF], and is thereafter referred to as the uniform convection approximation (UC). A less crude approximation is u app = û, i.e., the flow around a fixed disk and velocity -U e x at infinity in the absence of surfactant, that for convenience is called the translating disk approximation (TD). It can be viewed as the first step of an iterative scheme for solving the coupled flow-transport equations. All calculations involved in these two simplified models are detailed in Appendix B.

Within the uniform convection approximation (u app = -U e x ), all forces can be evaluated in closed form. As in the full problem, one finds that, depending on the imposed Pe, the viscous Marangoni force can be either motor or resistive. Indeed, at low Pe-the regime of interest for the bifurcation-the viscous Marangoni force exceeds the capillary force and is thus the dominant driving contribution. Specifically, we find F M ∼ -Pe ln Pe and F c ∼ Pe at low Pe. Looking at the force balance, this implies that there is no threshold value M c for swimming. Instead, a stable swimming point always exists, which is given by

Pe * UC = 8 exp - 16π M -γ E - 1 2 , for Pe * 1, (17)
where γ E is the Euler γ constant. This prediction is very different from the true swimming diagram, as shown in Fig. 4. It is actually possible to use the UC model to assess the magnitude of Marangoni convection and its relevance to solute transport (see Appendix B 1). It is found that Marangoni convection can not be neglected below a limit Péclet number and would likely alter the bifurcation properties.

104006-13 FIG. 6. Swimming diagram in the high-Pe regime. Points obtained from the full numerical solution are shown with solid red circles when stable and blue circles when unstable. Dash-dot line and grey circles represent the swimming points predicted by the uniform convection and translating disk approximations introduced in Sec. III B 2.

Turning now to the more elaborate translating disk approximation (u app = û), the problem must be solved numerically. Here again, we observe as before a change of sign of the viscous Marangoni force: from driving to resistive when increasing the Pe number. The new feature, however, is that spontaneous swimming is now obtained above a critical Marangoni number M c 8, in closer correspondence to the full problem. The threshold nonetheless remains significantly lower than expected and the associated swimming diagram also misses the bistable regime.

The swimming diagram thus appears to be very sensitive, both quantitatively and qualitatively, to the description of surfactant transport. This is somewhat expected close to a transition where susceptibility to perturbation is maximum. Accounting for Marangoni convection in the surfactant transport appears to be a key ingredient in the low-Pe regime, and its most conspicuous consequence is the bistability region where stable motionless and swimming states coexist. In addition to this qualitative change in the bifurcation, Marangoni flows also make the swimming regime more difficult to reach, with the threshold M c shifted to a higher value. This is in line with qualitative expectations that Marangoni flows will smear out concentration heterogeneities, thus delaying the onset of significant chemical asymmetries. In contrast with such a simple picture, however, we have seen that when considered in detail, Marangoni effects on swimming capabilities are more ambiguous. Besides its role on transport, we evidenced that flow-mediated viscous force systematically promotes symmetry-breaking in the low Péclet regime and this holds true for both approximations considered. As we shall see in the next section, in the high-Pe regime, the role of Marangoni flows is much less significant.

C. Asymptotic swimming regime (large Pe * )

We now focus on swimming characteristics in the large Pe, large M regime, much above the symmetry-breaking threshold. As mentioned above, this is the regime corresponding to most experimental conditions to date with camphor disks and more generally to macroscopic Marangoni swimmers where M M c . Figure 6 shows the stable steady state as a function of the Marangoni number M. In the explored regime, the swimming velocity Pe * is found to increase as a power law Pe * ∼ M α with an exponent α = 0.72 ± 0.02. Here again, such a behavior can be compared with expectations from the simplified models introduced above. The uniform convection approximation leads to an asymptotic swimming regime that, as discussed in Appendix B 1, also obeys a power 104006-14 law,

Pe * UC = (A UC M) 2/3 , with A UC = (3/4) 8π (5/4) . ( 18
)
Despite the drastic approximation involved in the UC model, the asymptotic swimming regime predicted has exponent and magnitude quite close to the full solution, as visible in Fig. 6. We proceed by considering the translating disk approximation which incorporates effects from the noslip BC at the disk surface. Strikingly, the numerical results obtained are very close to those of the full model (Fig. 6), with a difference never exceeding 10 %. Interestingly, within this simplified model it is possible to trace back the origin of the asymptotic power law. In doing so, we focus on the TD approximation on the capillary force, assuming that the total Marangoni force obeys a similar scaling [START_REF]Note that this is true for the UC approximation, where F c , F vM , and F M all follow the same scaling for large Pe[END_REF]. In the large-Péclet regime considered here, the solute emitted in the center of the swimmer will marginally diffuse before reaching the rear edge of the disk. Let us call w and δ, both much smaller than unity, the horizontal and vertical extent of the solute plume at the perimeter. Because of solute conservation, the injected rate J must balance the outgoing leak rate at the rear c wδU , which in dimensionless form yields a typical edge concentration,

c = Pe wδ -1 . ( 19
)
With w < 1, the capillary force is immediately obtained as Fc = c w so that, assuming that the total Marangoni force obeys a similar scaling as the capillary force, the force balance yields

16Pe * 3M = 1 Pe * δ(Pe * ) . ( 20
)
The swimming velocity is therefore fully determined once the thickness of the solute boundary layer is known as a function of Pe. We now consider two possibilities for the boundary layer thickness δ(Pe), the details of which are gathered in Appendix B 2. In the classical treatment [START_REF] Levich | Physicochemical Hydrodynamics[END_REF], one obtains δ(Pe) ∼ Pe -1/3 and therefore Pe * ∼ M 3/5 . However, such an exponent does not match the one obtained numerically in Fig. 6 and suggests an underestimation of the driving force. In fact, the boundary layer considered so far is relevant for a plate that is semi-infinite. What we need instead is the boundary layer just below the disk edge. As shown in Appendix B 2, the flow field has then a different structure, because vertical convection becomes dominant and squeezes the boundary layer, leading to

δ(Pe) ∼ Pe -2/3 , Pe * TD ∼ M 3/4 . ( 21
)
This regime is very close to the one observed in numerical simulation, both for the translating disk model and for the full problem. We note that Eq. ( 21) differs from the prediction obtained by Ender et al. [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF] for hemispherical swimmers. This shows that the swimmer shape, by altering the boundary layer, may also influence the swimming properties.

To conclude this section, the asymptotic swimming regime is very different from the low-Pe regime in that it depends only weakly on the level of description retained. While the oversimplified but fully tractable uniform convection model predicts a slightly lower exponent than the numerical value, a decoupled model approximating the advective transport by the surfactant-free flow field yields nearly quantitative agreement with the exact solution and allows rationalizing the scaling of swimming velocity from boundary layer considerations. In stark contrast with the bifurcation regime, this demonstrates that in the asymptotic swimming regime, forces acting on the swimmer, and therefore the propulsion velocity, are barely influenced by transport from Marangoni flows. This is consistent with the previous evidence that force balance in the asymptotic regime is ruled by solute boundary layer properties, and can be further undertood in the following way. In the limit of Pe 1 and M 1, where fore-aft asymmetry is fully established, one expects-in line with recent velocimetry measurement on asymmetric swimmers [START_REF] Sur | Translational and rotational motion of disk-shaped Marangoni surfers[END_REF]-that Marangoni flows uû have for typical magnitude the propulsion velocity U , because both have the same physical origin. Therefore, the Péclet number will remain the relevant parameter for the boundary layer, whose properties remain essentially unchanged by Marangoni transport. The same holds true for the conservation arguments and the associated force balance developed above, thus explaining the minor role of Marangoni flows. As a last remark, estimation of Marangoni flows in the UC model (see Appendix B 1) also suggests that they would marginally affect swimming properties in the asymptotic regime.

IV. EXTENSIONS

A. Swimmer with spatially extended source

Our swimmer so far involves a point-source of surfactant, with emission only at the disk center. This assumption was convenient to build a toy model amenable to analytical treatment [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF]. If localized release is found in asymmetric swimmers such as camphor boats [START_REF] Shimokawa | Power law observed in the motion of an asymmetric camphor boat under viscous conditions[END_REF][START_REF] Suematsu | Quantitative estimation of the parameters for self-motion driven by difference in surface tension[END_REF], then the most common situation for symmetric swimmers is a uniform release all over the disk surface [START_REF] Nagayama | A theoretical and experimental study on the unidirectional motion of a camphor disk[END_REF][START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF][START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Nakata | Self-motion of a camphanic acid disk on water with different types of surfactants[END_REF][START_REF] Suematsu | Self-motion of a Benzoquinone disk coupled with a redox reaction[END_REF][START_REF] Nakata | Mode change in the self-motion of a benzoquinone disk coupled with a NADPH system[END_REF][START_REF] Iida | Experimental and theoretical studies on the self-motion of a phenanthroline disk coupled with complex formation[END_REF]. We now examine the effect of such an extended chemical source. The surfactant is emitted over a disk with radius a s with a surfacic flux J = J /π a s 2 so that the total flux J is conserved. For convenience, we introduce the notation ε = 1a s /a and refer to the limit ε → 0, where the entire disk acts as a source, as "full release."

The swimming diagram for swimmers with fixed total flux but an increasingly larger source is shown in Fig. 7. Compared to punctual source case, the qualitative behavior is the same but there are some quantitative changes as the source grows in size. The most conspicuous is the shift in threshold M c and M c to much higher values, discussed below, and which results in a growing extent of the bistability region. Besides, the branch of unstable swimming tends to approach the stable branch above, thus widening the basin of attraction of the motionless state. By contrast, one observes at large Péclet numbers only a small shift toward lower swimming velocity but the asymptotic behavior shows little change. In fact, all the swimming points approach a master curve, which is again a power law Pe * ∼ M α , with an exponent slightly weaker than in the point source regime α = 0.67 ± 0.02 instead of 0.72. It could a priori be expected that the source size influences the portion of perimeter where capillary forces apply, thus possibly altering the scaling law derived above. Yet, well above the threshold, swimming is overall very similar with a punctual or an extended source.

The major effect of an extended source is to make the swimming more difficult, and we now focus on this shift in threshold. When the source radius approaches the swimmer radius (ε → 0), the threshold actually diverges (see inset in Fig. 7), implying that a fully releasing disk would not swim. To understand the origin of this phenomenon, we used again the UC and TD approximations which allow to test the relevant contributions. It turns out that the divergence appears only when including the Marangoni flow contribution to the surfactant transport. As discussed in detail in Appendix D, we ascribe the divergence to our specific choice of boundary conditions and geometry.

The assumption of a fixed normal flux on the source and zero outside imposes a discontinuity at the source edge and results in a divergence of the concentration gradient. As the source edge approaches the disk perimeter, the surfactant concentration becomes more and more homogenous and less asymmetric, decreasing the driving forces that eventually becomes insufficient to maintain motion. This phenomenon is peculiar to our configuration because the discontinuity in flux normal to the swimmer (and surface) can approach the triple line. It would disappear for swimmers that are spherical, elliptical or cylindrical in shape. The ε divergence is thus a combined effect of our fixed flux boundary condition and zero-thickness disk geometry. As such, it is a consequence of our idealized description and is probably not relevant to experimental system [START_REF] Speculatively | however, it could broaden the bistability region, for very thin swimmers[END_REF].

If the -divergence can be regarded as a specificity of the geometry, then Fig. 7 suggests that the key quantity controlling the swimmer velocity is the total release rate J and not the source size. We have tested experimentally this prediction by using core-crown swimmers where camphor is emitted only in the central region. Practical details and data are provided in Appendix G. In the asymptotic regime where Pe * ∼ M α , and for swimmers of similar size, the ratio U/J α should remain constant. Choosing α = 0.68 in agreement with experiments (see below discussion of Fig. 8), we indeed found a good rescaling for the velocity of core-crown swimmers, even though error bars preclude a stringent test of the exponent value. Overall, this suggests that as a first approximation, the pointsource assumption is a reasonable one.

B. Inertial effects: Finite Schmidt and Reynolds numbers

Results presented so far all assumed a Stokes flow, i.e., a vanishing Reynolds number corresponding to a negligible inertial term in the Navier-Stokes equation. However, a typical value for camphor disks [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF] is Re 3 × 10 2 , implying that inertial effects are significant and that the fluid 104006-17 flow should be described with the Navier-Stokes equation. As can be seen in our Eq. ( 7) for ũ, it is natural within the chosen dimensionless form to characterize inertial effects with the Schmidt number Sc = Pe/Re and we do so thereafter. Note that it is also consistent with the fact that we focus on a velocity which is selected by the swimmer and not imposed externally: Reynolds number will vary while maintaining a fixed Schmidt number set by the fluid properties. On the practical side, viscous boundary layers develop below the swimmer and at the liquid-gas interface behind. A specific mesh was designed to accommodate those effects.

The influence of inertial effects on the swimming curve is shown in Fig. 8. Overall, it is quite modest although quantitatively significant especially in the asymptotic regime. In the low-Pe regime close to the bifurcation, there is only a 25% reduction of the swimming threshold M c and a small increase of swimming velocity. In both cases, occurrence of spontaneous swimming is facilitated by the intrinsic nonlinearity of the flow equation. An opposite effect is seen in the high-Pe regime where compared to the Stokes limit, swimming velocities are reduced. The asymptotic regimes nonetheless preserves a power law Pe * ∼ M α with an apparent exponent decreasing slightly as inertial effects grow: α 0.72, 0.65, 0.62, 0.60 for Sc = ∞, 10 3 , 10 2 , 10, respectively. As we detail in Appendix C, finite Re numbers induce a viscous boundary layer and an even thinner concentration boundary layer. A scaling argument then predicts for the swimming velocity Pe * ∼ M 1/2 Sc 1/6 , as found in Ref. [START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF]. Such an exponent should be considered as a lower bound. Indeed, we saw in the Stokes regime that proper scaling requires accounting for edge effects which make the boundary layer thinner. If a similar phenomenon occurs here, then one would expect a higher exponent, as observed numerically.

We are finally in a position to compare our numerical prediction to experimental data. Figure 8 includes results from experiments we carried on with camphor disks [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF]. The M values could be computed from measurements of the dissolution rate for immersed swimmers. Our numerical findings indicate that the swimming curve in the high-Pe regime is a power law with exponent α in the range 0.6-0.75. In particular, for the experimentally relevant Schmidt numbers around Sc = 10 3 , we find an exponent α = 0.65. This is certainly close to the value α = 0.68 found in experiments and to the prediction α = 2/3 of the toy model put forward in Ref. [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF], which is essentially equivalent to the present UC approximation. Thus, the predictions are very similar in an oversimplified picture and in a full numerical treatment accounting for all complexities, including Marangoni flows and inertial effects, and both are compatible with experiments. In light of these results, we conclude that the success of the toy model originates in a series of compensations between neglected physical effects (including no-slip boundary condition, Marangoni transport, finite Reynolds) and as such should be viewed as somewhat incidental. A final comment is in order regarding the prefactor in the power law. Extrapolation of the Stokes results to the much higher M values typical in experiments leads to an overestimation of swimming velocities by almost an order of magnitude. Accounting for inertial terms significantly reduces the discrepancy but the predictions nonetheless appear to remain above experimental data. This suggests that other modeling assumptions-from fixed-flux BC to the exchange kinetics or surface distortion-may need to be revisited to reach a perfect quantitative match.

C. Concentration swimmer

We have focused so far on swimmers emitting surfactant at a prescribed flux. As a final extension of our study, we now investigate a different type of swimmer where the boundary condition for surfactant release is instead a fixed concentration. For brevity, we refer to those two classes as "flux swimmers" and "concentration swimmers." They correspond to the two simplest choices for boundary conditions, Neumann and Dirichlet, respectively. The assumption of fixed concentration is indeed expected to be relevant in the numerous instances where the swimmer is entirely made of the surfactant itself. These include droplets [START_REF]The droplets considered here float are the air-liquid interface and are not moving through a fluid[END_REF] such as aniline oil and pentanol [START_REF] Chen | Self-propelled motion of a droplet induced by Marangoni-driven spreading[END_REF][START_REF] Nagai | Mode selection in the spontaneous motion of an alcohol droplet[END_REF] or swimmers made of solid, such as benzoquinone, aspirin, or pure camphor [START_REF] Satterwhite-Warden | Co-operative motion of multiple benzoquinone disks at the air-water interface[END_REF][START_REF] Bánsági | Motion and interaction of aspirin crystals at aqueous-air interfaces[END_REF][START_REF] Nakata | Self-rotation of a camphor scraping on water: New insight into the old problem[END_REF]. Applying the same method as above and going back to the Stokes regime, we find that concentration swimmers with partial release do swim and their swimming velocity Pe * was obtained for different source radii (a s = 0.2, 0.4, 0.6). Because concentration swimmers involve a Marangoni number different than flux swimmers, comparing their swimming curves is not straightforward. Nevertheless, one can measure the total flux J * released by the concentration swimmer in its swimming state and associate the effective Marangoni number M * = κJ * /ηD 2 , as defined previously. Doing so, one can compare the swimming curves for flux and concentration swimmers, as shown in Fig. 9(a). The curve Pe * (M) for flux swimmers and the curve Pe * (M * ) for concentration swimmers are in fact identical [START_REF]A different notation was introduced to emphasize that the former and the latter are different functions[END_REF]. Said otherwise, the velocity is the same for a swimmer with prescribed flux J and a concentration swimmer releasing a flux J in its swimming state. This suggests that at least in the conditions considered-sufficiently high M and release only partial-the release mechanism influences the velocity only through the total flux emitted, which is in agreement with finding of Refs. [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF].

Exploiting this observation, we can now predict the velocity of a concentration swimmer. Let us assume that the velocity of flux swimmer is Pe * = F (M), with F a known function and that the total flux emitted by a concentration swimmer moving at an arbitrary velocity Pe is a s Dc o G(Pe), where a s Dc o is the total flux in a purely diffusion situation and G is a function. Given those assumptions, the velocity Pe * of the concentration swimmer is solution of

Pe * = F[ M G(Pe * )], M ≡ κa s c o ηD , ( 22
)
where M defines a Marangoni number for a concentration swimmer with fixed value c o within the source radius a s . Let us apply Eq. ( 22) to the asymptotic swimming regime, where both F and G are power laws. In this case, Pe * = F (M) ∼ M α as demonstrated previously. Regarding the flux released by a concentration swimmer, there are two limiting cases to consider.

(i) Swimmer with nearly full release (a s a). One expects a local flux ∼Dc o /δ(Pe), with δ(Pe) = a s Pe -μ the thickness of the diffusive boundary layer and μ = 1/3 [106] leading to a total flux ∼ a s 2 × Dc o Pe μ /a s , such that G(Pe) ∼ Pe μ . Equation ( 22) then gives the swimming velocity Pe (ii) Swimmer with very small source (a s a). As the source size vanishes, at a given Pe one can expect the convection to become less significant. We have observed numerically that in the range of Pe considered, the flux released actually shows limited variations with Pe, it changes by 20% for the source radius a s = 0.6 and even less for smaller radius. A simple approximation is then to consider the flux as independent of Pe, which leads to

Pe * ∼ M α . ( 24
)
The numerical data for Pe * (M ) is shown in Fig. 9(b). It lies in between the limiting cases (i)

and (ii) and as the source size decreases, it comes closer to the later. We conclude that at least in an approximate manner, the swimming curve of a concentration swimmer with partial release may be deduced from the swimming curve of a point-source swimmer. As a final case, we consider the fully releasing swimmer, where concentration is imposed all over the surface, including at the perimeter. We find a very different behavior: there is no spontaneous motion. Because the concentration along the swimmer edge is fixed everywhere to the same value, the capillary force vanishes. The viscous Marangoni force alone remains to yield a motor contribution but is insufficient to maintain steady motion. Accordingly there is no swimming point over the range of M and Pe explored. We surmise that this absence of swimming is not specific to the disk geometry considered here and that it would also apply to systems with other shapes, such as hemispherical or ellipsoidal, but this would have to be proven. In any case, our observation raises some questions: do genuine concentration swimmers actually exist? For all experimental system mentioned above where the swimmer is made of the surfactant itself and that do swim, is it legitimate to assume that fixed concentration is imposed? This remains to be clarified.

V. CONCLUSIONS AND PERSPECTIVES

To summarize, we investigated the steady motion of interfacial disks releasing isotropically a soluble surfactant. Using a full numerical description to account for all couplings at work, we characterized the spontaneous swimming velocity and assessed the contribution of various forces by comparison with simplified models. In this respect, it turns out that surface tension effects are more complex than usually expected and are not subsumed in the capillary force; the contribution induced by viscous Marangoni flows may be motor as well. Our main findings can be summarized as follows.

(i) The transition to swimming is discontinuous, occurs above a critical Marangoni number, and involves a forbidden range of velocity and a bistability region, where steady motion and rest state are both possible. A full account of Marangoni effects is crucial to capture these features, which are absent in simplified models. (ii) In the regime of large Pe and M relevant to most experiments, the swimmer velocity follows the power law Pe ∼ M α with an exponent α close to 3/4. Here Marangoni flows have little effect on the force balance for the swimmer and accordingly simplified descriptions may capture the swimming behavior.

Throughout this study, we have focused exclusively on the various contributions that enter the force balance and on the most accessible observable: the swimming velocity. Yet, a deeper understanding may be possible by examining local properties such as the flow velocity and the concentration fields. Preliminary inspection suggests that characterizing those fields is not obvious especially when it comes to evaluate the accuracy of simplified models. Indeed, relevant characteristics of flow and concentration fields may depend subtly on the observable of interest, be it the force on the particle, or the solute dispersion, etc. We defer such an endeavor to a subsequent work.

Even for individual self-propulsion, there are several questions left open. First, is the discontinuous swimming transition observable? Swimmers of colloidal size [START_REF] Dietrich | Microscale Marangoni Surfers[END_REF] might be plausible candidates to reach low values of Marangoni number at threshold, but only if the symmetry can be maintained accurately as it is precisely expected to matter at the bifurcation [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF]. Second, we addressed here the limit case of soluble surfactants with fast equilibration kinetics but exploring more general situations is also interesting. Indeed, it is known that in the case of a fixed source, the nature of the surfactant-insoluble or soluble-modifies the scaling of the velocity field [START_REF] Bandi | Hydrodynamic Signatures of Stationary Marangoni-Driven Surfactant Transport[END_REF][START_REF] Mandre | Axisymmetric spreading of a surfactant driven by self-imposed Marangoni stress under simplified transport[END_REF], and that for hemispherical swimmers, evaporation effect may be relevant [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF]. Third, we have focused essentially on steady swimming. Experimentally, the swimming modes are much richer with the observation of intermittent motions [START_REF] Suematsu | Synchronized intermittent motion induced by the interaction between camphor disks[END_REF][START_REF] Suematsu | Evolution of self-propelled objects: From the viewpoint of nonlinear science[END_REF][START_REF] Akella | Dynamics of a camphoric acid boat at the air-water interface[END_REF] whose description would deserve dedicated investigations.

Finally, beyond the individual behavior, a broader perspective is to address the role of Marangoni forces in the other situations that have been explored experimentally, including interaction with a wall or with another swimmer, and collective effects that emerge in assemblies. Given the numerical demand and difficulty of a complete description, the only route to address such phenomena is to design simplified models for swimmers. We believe point source models discussed here remain useful approaches for guiding the exploration of more complex situations [START_REF] Lippera | Alignment and scattering of colliding active droplets[END_REF] with the present work delimiting their capabilities and helping for further improvements.

where γ E is the Euler γ constant. In the opposite limit of large Péclet, we have

Pe 1, Fc 1 √ π Pe + O(Pe -3/2 ), FM A √ π Pe + O(Pe -3/2 ), A≡ 2 3 √ π (3/4)
(5/4) 0.509. (B4) Examination of limit cases reveals that the viscous Marangoni force F vM = F M -F c changes behavior with Pe. At low Pe, it is a driving contribution, being actually dominant over the capillary force because of the Pe ln Pe scaling. However, it changes sign to become a resistive, and at higher Péclet cancels approximately half of the driving capillary force. Finally, combining Eq. (B3) for FM and Eq. ( 15) leads to Eq. ( 17) for the swimming curve at low Pe * , whereas in the asymptotic regime, one recovers Pe * ∼ M 2/3 , as stated by Eq. ( 18).

b. Relevance of Marangoni convection

The uniform convection model reduces the convective contribution to solute transport to that of an homogeneous velocity field, thus neglecting transport by Marangoni flows. Here we briefly discuss the range of validity of this assumption.

One simple criterion to assess when Marangoni convection necessarily becomes relevant is to compare the drift velocity U and the expected Marangoni velocities at the front u Mar. of the swimmer. Indeed, when u Mar. > U , Marangoni flows are strong enough to revert the advection transport, thus inducing solute transport toward the front that will significantly dampen the fore-aft asymmetry required for swimming. The BC at the free surface is expressed as ∂ z ũ = -∇ s c. Assuming that the scale of variation in the vertical direction is a-no other one exists at low Pe-gives the typical surface velocity ũMar. ∼ -∇ s c. Using the solute concentration Eq. (B1) to estimate the Marangoni velocity close to the front (r 1, θ = 0), flow reversal occurs for

e -Pe 2π [1 + Pe] - Pe M > 0. (B5)
This equation provides a minimum Péclet number below which Marangoni convection cannot be neglected. It can be solved in asymptotic regimes of low and high Pe:

Pe > M 2π
for Pe 1, Pe > ln M 2π for Pe 1.

(B6)

Figure 10 shows where this condition is satisfied in the swimming diagram of the UC model. Clearly this is not the case in low Pe region. However, in the asymptotic regime Pe * 1, the predicted swimming velocity becomes much larger than the limit where Marangoni effects come into play thus suggesting that this simple approach becomes valid, in agreement with the discussion given in the main text (Secs. III B 2 and III C). Note finally that when the estimated Marangoni convection exceeds the moving velocity, the distribution of surfactant becomes more isotropic and bifurcation to a swimming state is expected to be delayed. In this regard, it is interesting to remark that the actual (discontinuous) bifurcation occurs at a critical Pe * c 3 which is close to the value Pe 3.4 where the swimming curve becomes invalid according to the criterion of Eq. (B5).

Translating disk approximation

We now examine a more elaborate treatment of the advective transport, where u app = û, i.e., the Stokes flow generated by an inert moving disk. Working in the disk reference frame, the disk is fixed at the frame origin and submitted to a flow with velocity -U e x at infinity. This approximation thus retains the no-slip boundary condition on the disk while discarding the effect of Marangoni flows on the solute transport.

The interest of the approach is once again to obtain a transport problem which can be solved independently, with solution parametrized by the velocity U , using the known flow field û. Instead of the surface flow given in the main text by Eq. (B7), we now use the full tridimensional solution. Using Cartesian coordinates, it reads as [START_REF] Happel | Low Reynolds Number Hydrodynamics[END_REF] û

U = 2 3π 2x 2 λ 3/2 (λ + 1) 2 - √ λ λ + 1 + 3 arcsin 1 √ λ + 1 -1, v U = 2 3π 2x ỹλ 3/2 (λ + 1) 2 , ŵ U = 2 3π 2xz √ λ (λ + 1) , ( B7)
with the notations

= -1 + x2 + ỹ2 + z2 , = 4z 2 + 2 , λ = 1 2 ( + ). ( B8
)
Let us note that at the free surface z = 0, the flow field can be conveniently expressed in polar coordinates (see notations in Fig. 1) as

û U = 2 3π 3 arcsin 1 r - π 2 (cos θ e r -sin θ e θ ) + √ r2 -1 r2 (cos θ e r + sin θ e θ ) . (B9)
In the regime of large Péclet, one expects a concentration boundary layer to appear whose thickness results from a balance between advective and diffusive terms. We now discuss the possible scalings for this thickness, which were used in Sec. III C to obtain prediction on the asymptotic swimming regime. Consider first the central region, with = 1 -(x 2 + y 2 ) = O(1), where solute is released. Equation (B8) can be expanded in the vicinity of the surface to give

û U = - 8 3π |z| √ + O(z 2 ), ŵ U = 4 3π xz|z| 3/2 + O(z 3 ). (B10)
Following classical treatment of boundary layers [START_REF] Levich | Physicochemical Hydrodynamics[END_REF], one rescales the z direction by the boundary layer thickness δ according to Z = z/ δ and look for the scaling δ ∼ Pe -m such that advective diffusive terms appear at the same order Note that in this standard configuration, vertical and horizontal advective terms are involved at the same order. As explained in the main text, Eq. (B12) yields for the swimming velocity predictions that are not consistent with numerical results. The reason is that the discussion so far has focused on the central region of the swimmer, whereas we need to evaluate the concentration near the edge of the disk. The velocity BC switches at the disk edge from no-slip to stress-free and as a consequence, the velocity scalings are modified, leading to a different thickness for the solute boundary layer. Specifically, when x2 + ỹ2 = 1, Eq. (B8) gives

Pe(ZPe -m ∂ x c + xZ 2 Pe -m ∂ Z c) ∼ Pe 2m ∂ 2 ZZ c, ( B11
û U = - 2 3π (3 + ỹ2 ) |z| + O(z), ŵ U = 2 3π |z| + O(z). (B13)
Performing a balance between advection and diffusion in this region now leads to

Pe( √ ZPe -m /2 ∂ x c + √ ZPe m /2 ∂ Z c) ∼ Pe 2m ∂ 2 ZZ c. (B14)
As a consequence, convection is now dominated by advection in the vertical direction and the exponent for which it is balanced by diffusion is twice higher than previously: m = 2/3. The vertical flow pinches the boundary layer, making it thinner δ(Pe) ∼ Pe -2/3 . (B15)

As discussed in the main text, this modified scaling provides a swimming behavior in the large-Péclet regime that is consistent with the numerical result. As explained in Sec. IV A, the swimming threshold M c diverges when the source radius approaches the swimmer radius (a s → a or ε → 0). Here we discuss a possible origin for this effect.

We consider a motionless disk releasing a chemical solute all over its surface with prescribed flux and ask what the concentration field is when diffusion alone is at play. Choosing units so that disk radius, diffusion coefficient, and flux are all unity, the answer is

c(r, z) = ∞ 0 dkJ 1 (k)J 0 (rk)e -kz dk k , z > 0, (D1)
with J n the Bessel function of order n. In the limit z → 0, the integral can be performed explicitly and the result for the concentration at the surface c(r, 0

) ≡ c(r) is then c(r) = 2 π E (r 2 ), r < 1, = 2 π r [r 2 E (r -2 ) -(r 2 -1)K (r -2 )], r > 1, (D2)
where K and E are the elliptic integrals of first and second kind, respectively. A noticeable feature is the logarithmic divergence of c (r) at the disk edge, namely for r → 1,

-c (r) = - 1 π 2 + ln |r -1| 8 + O(r -1). (D3)
We show in Fig. 11 the concentration along the x axis, for a fully releasing swimmer (a s = 1). A divergence is indeed seen at the disk edge and the concentration below the disk (|x| < 1) is not far from the pure diffusion prediction approximation given by Eq. (D2) [112].

What consequences are to be expected when the divergence in concentration gradient approaches the disk edge? On the one hand, the surfactant is released closer and closer to the edge. On the other 104006-26 hand, the strong Marangoni stresses acting on the fluid in the immediate vicinity of the wall presumably induce some fluid motion outward. Both effects combine to make the surfactant distribution along the edge more and more homogeneous as ε → 0. This trend is visible in Fig. 11(b), which displays the concentration profile at a fixed M. As the source radius a s increases, the concentration profile along the x axis gradually loses its asymmetry and come closer to the pure diffusion profile. As a consequence, the surfactant concentration tends to be less heterogeneous along the edge, and the Marangoni flows also becomes more symmetric, with increasingly canceling contributions. With the capillary and viscous Marangoni forces both subsiding, steady motion becomes impossible. For instance, this happens for a s = 0.75 in Fig. 11. Because of the coupling between diffusion, Marangoni flows and the influence of the edge, it appears difficult to put this qualitative picture on a quantitative basis. A justification of the scaling observed for M c (ε) in Fig. 7 thus remains elusive.

We note that the behavior observed here is particular to our configuration. Would the releasing surface be a hemisphere as in Refs. [START_REF] Ender | Surfactant-loaded capsules as Marangoni microswimmers at the air-water interface: Symmetry breaking and spontaneous propulsion by surfactant diffusion and advection[END_REF][START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF], the concentration field in the pure diffusion problem would be c(r) ∼ a 2 J /(Dr), which is regular everywhere for r a. The singularity would also disappear for an ellipsoid and for a cylindrical swimmer of finite thickness, because the flux at the triple line is parallel to the interface. We conclude that the ε divergence is associated to a mathematical divergence induced by the combination of our fixed flux boundary condition and zerothickness disk geometry. It must be added that due to a number of possible factors (slightly deformed interface, swimmer partly immersed, finite exchange kinetics, etc) regularization is expected in a real swimmer.

APPENDIX E: ON CAMPHOR EVAPORATION

We did not account for surfactant evaporation in our modeling. Even though evaporation is often mentioned for explaining surface self-cleaning and the camphor boats ability to cross its own wake, experiments with surfactant injection [START_REF] Roché | Marangoni Flow of Soluble Amphiphiles[END_REF] show that bulk dissolution alone can provide an efficient surface-cleaning pathway at short time. We thus chose to keep the parameters minimal by not introducing this additional mechanism. We present here experimental measurement of the evaporation rate and a discussion that further rationalize this assumption.

We measured the change in concentration in a camphor solution due to evaporation. In practice, we used the same conditions as those used for our experimental study of camphor swimmers, as described in Ref. [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF]. A rectangular pool (17 × 27 cm) with height h = 1 cm is filled with a solution of camphor (1 g L -1 ) and the mean concentration C m (t ) in the solution is measured by spectroscopy. The time evolution of C m (t ) is shown in Fig. 12.

These data can be rationalized with a simple model. We assume that the flux of camphor leaving the liquid through evaporation can be written as

J = k[C(0) -C gaz ] = -D ∂ z C| z=0 , (E1)
Here, C(z) is the concentration at position z, the interface is at z = 0, the concentration in the gaz C gaz is taken as negligible, and the evaporation constant k quantifies the interfacial transport from bulk liquid to the gaz. Equation and a fit to the data yields k = 6 × 10 -7 m s -1 . Two remarks are in order. First, we note that such a value leads to λ = 10 -3 m and Bi = 10. A Biot number much above unity indicates that the concentration can not be approximated as uniform, if diffusion alone is at play. However, in our experiments, periodic sampling of the solution and spontaneous convection due to external perturbation ensure some mixing of the solution, so the uniform concentration remains a reasonable approximation. Second, our value for k is well below the estimate made in Ref. [START_REF] Soh | Dynamic self-assembly in ensembles of camphor boats[END_REF], which suggests an apparent evaporation rate k = 10 -4 m s -1 . Whereas the method there was indirect, our measurement is direct and involves only evaporation.

We can now discuss the relevance of evaporation effects for the swimmer. Taking for the typical size the radius a = 4 mm yields a Biot number Bi 3.4. However, this value is only an upper estimate. Indeed, the size which is relevant in this case is the thickness of the concentration boundary layer. As discussed in the main text and Appendix B, it scales as δ = aPe -1/3 in the high-Pe regime. For the typical camphor swimmer (a = 4 mm, v = 70 mm s -1 , Pe = 3 × 10 5 ), we find Bi = 4 × 10 -2 .

(E4)

While this result was derived in the Stokes regime, on can expect even smaller values for finite Reynolds numbers. In any case, a small Biot number indicates that diffusion acts on smaller timescale than evaporation: the concentration in the solute layer released by the swimmer will thus decrease primarily because the layer spreads by diffusion and not because some solute is lost by evaporation. Finally, as a last argument, we also stress that the scaling arguments developed in the main text describes what happens below the swimmer, where evaporation does not occur. Only behind the swimmer could evaporation becomes significant, but the force balance is dominated by the concentration in the vicinity of the edge. This again suggests that evaporation should not play a key role in our system. Note finally that Ender et al. [START_REF] Ender | From diffusive mass transfer in Stokes flow to low Reynolds number Marangoni boats[END_REF] showed that for a hemispherical swimmer in the asymptotic regime, evaporation can be neglected if Bi Pe 1/5 . Should the same criterion be applicable to our disk-shaped swimmer, it would clearly be fulfilled.

In conclusion, the small value of Biot number indicates that evaporation has no significant influence on the self-propulsion and accordingly it was discarded in our modeling. Nevertheless, we believe evaporation does play a role for experiments where the water pool is finite. By eliminating some of the camphor released by the swimmer, it reduces accumulation within the pool and 104006-28 delays the point at which the interface is so charged in camphor that self-propulsion stops. This phenomenon is relevant on the timescale of hours.

APPENDIX F: FINITE-SIZE EFFECTS AND CONVERGENCE IN MESH

Our default configuration for the FEM calculations relies on the box and mesh defined in Sec. II B. To evaluate the finite size effects and ascertain the convergence in mesh, we have changed systematically the box size L by a factor S f , and the mesh size by a factor M f . The refinement in mesh applies everywhere, so the relative resolution between the various areas defined is maintained. For a given quantity, we compute the relative difference with the result obtained in a reference configuration with the largest box and finest mesh (S f = 2, M f = 1.5). The result for the various forces considered in this work are shown in Fig. 13, for fixed physical parameters Pe = 3, M = 230. We note that though the ranges of S f and M f are quite broad (0.1-2 and 0.05-1.5, respectively), most forces deviate form their reference value by less than 3%. The only exception is the viscous force F v when computed with Eq. [START_REF] Dai | The mechanism of surfactant effects on drop coalescence[END_REF]. However, using the alternative expression given by Eq. ( 12), the convergence is much faster. Therefore, this is the route that was used throughout this work. Overall, we see that the forces in the default configuration (S f = M f = 1) have reached values that depends very little on the box size and mesh refinement, with variations smaller than a few percent.

APPENDIX G: VELOCITY OF CORE-CROWN SWIMMER

We investigated in Sec. IV A how the swimmer velocity depends on the spatial extent of the source, i.e., the chemically active area. Here we address this question experimentally.

To modulate in a controlled manner the extent of the chemical source, we designed "core-crown" swimmers, made of two assembled parts. The core is loaded with camphor while the crown is not. We considered swimmers with radius a = 5 mm and an increasing source radius a s . In contrast to the numerical model, the total flux J is not constant but grows with the area of the source A s . Accordingly, the velocity increases as well. We can, however, account for this variable total flux. Assuming Pe * ∼ M α one would expect for v/A α s a constant value. Using for α the exponent found from the experimental velocity of full-disk swimmer with various radii (α = 0.68 ± 0.02, see Fig. 8) we indeed find a good rescaling. As shown in Fig. 14, data points normalized in this way remain close to unity (dashed line), even though the error bars prevent a definitive assessment of the exponent value. The main factor governing the velocity seems to be the total flux, the influence of source size appears to be much weaker.

104006-29

FIG. 1 .

 1 FIG. 1. (a), (b) Schematic of the symmetric Marangoni swimmer studied (a) and of the simulation box (b) together with notations used. The swimmer reference frame is considered with origin at the disk center O, and with a flow velocity -U e x imposed far away from it. Note that the picture is not to scale, since L/a = 400 in practice. (c) Snapshot from simulation illustrating the self-generated surfactant asymmetry around the swimmer and the associated flow streamlines. The red wake represents the accumulation of surfactant at the back.

104006- 8 FIG. 2 .

 82 FIG. 2. Evolution of dimensionless forces: the capillary F c , viscous Marangoni F vM , Marangoni F M , and Stokes drag Fv forces as a function of Pe for M = 300. Positive forces (solid symbols) correspond to motor contributions, and negative (empty symbols) to resistive ones: the viscous Marangoni force F vM changes sign around Pe = 2.4. Stable spontaneous swimming occurs at Pe * = 4.8.

11 FIG. 4 .

 114 FIG. 4. Swimming diagram in the low-Pe regime, near the onset of spontaneous motion. Stable and unstable swimming points obtained from full numerical solutions are shown, respectively, with solid red and empty blue circles. A region of bistability (purple colored area) occurs for Marangoni number in the range [M c 255, M c 410]. The dash-dot line and grey circles represent the swimming diagram predicted by simpler models, respectively, uniform convection and translating disk approximations (see Sec. III B 2). Note that the gray hatched low-Pe region corresponds to a parameter range for which finite size effects may affect the data [84].

FIG. 5 .

 5 FIG. 5. Relaxation toward stable state, as probed by time-dependent numerical calculations. (a) Swimming diagram showing steady-state solutions Pe * (M) and initial conditions Pe in spanning the different regimes. Relaxation from a point initially in the gray hatched area is expected toward the lower value of Pe * . (b) Relaxation of the velocity Pe(t ) for the different initial conditions. Top: M > M c ; middle: M c < M < M c ; bottom: M < M c .

104006- 15 FIG. 7 .

 157 FIG. 7. Swimming diagram for a swimmer with spatially extended source. The source radius is a s = 0 (•), 0.5 ( ), 0.6 ( ), 0.7 ( ), 0.75 (♦); as in previous figures red solid (respectively blue empty) symbols stand for stable (respectively unstable) swimming points. (Left inset) Threshold M c as function of ε = 1a s /a. Source radius symbols are identical to main figure, with two additional values a s = 0.85 ( ) and 0.91 (8). (Right inset) Swimming curve in the high-Pe regime. Gray circles represent point source results for the translating disk approximation.

104006- 16 FIG. 8 .

 168 FIG. 8. Swimming diagram at finite Schmidt number Sc (or nonzero Reynolds number Re). (a) Low-Pe: evolution of the near-bifurcation regime. Stable and unstable points are indicated with filled and void symbols, respectively. (b) High-Pe: asymptotic swimming regime. Experimental data is obtained from Ref. [47]. The dashed line is a power-law fit, giving an exponent 0.68 ± 0.02.

104006- 18 FIG. 9 .

 189 FIG. 9. (a) Swimming velocity Pe * (M * ) for concentration swimmers of various source radius a s . The red square symbols show Pe * (M) for flux swimmer. (b) Swimming curve Pe * (M ) for concentration swimmer. The two dashed lines show the power law expected in the asymptotic regime, with exponent 1 and 3/4 (see text).

104006- 23 FIG. 10 .

 2310 FIG. 10. Relevance of Marangoni convection in the uniform convection model. The solid line is the swimming curve and the shaded area is the domain where Marangoni convection can not be neglected.

  ) thus implying m = 1/3 and δ(Pe) ∼ Pe -1/3 .

FIG. 11 .

 11 FIG. 11. Surfactant concentration c(x) along the x axis (y = z = 0) as found by numerical calculation, for different source radius a s . M = 2000 is fixed and Pe is close to the swimming value Pe * (between 22.7 and 27) when defined, and Pe = 23 otherwise. (a) Full releasing swimmer with a s = 1. The dashed line is the prediction of Eq. (D2). (b) Increasing source radius. The hatched areas correspond to the swimmer extent.

λ 27 FIG. 12 .

 2712 FIG.12. Time evolution of mean camphor concentration in a solution subject to evaporation. The line is a fit to Eq. (E3).

FIG. 13 .

 13 FIG. 13. Finite size effects (left) and convergence in mesh (right) on the various forces considered in this work. S f and M f are the size and mesh factors, respectively, and F is the difference with the value in the reference configuration F ref . The physical parameters are Re = 0, M = 230, and Pe = 3.

ACKNOWLEDGMENTS

We thank F. Nadal for fruitful discussion upon initiating FEM simulations. This work was supported by the French research programs Grant No. ANR-16-CE30-0028, and IDEX-LYON of the University of Lyon in the framework of the French program "Programme Investissements d'Avenir" (Grant No. ANR-16-IDEX-0005). FIG. 14. Rescaled swimming velocity of core-crown swimmers made of a chemically active core (orange) and passive crown (grey). A s /A is the ratio between the source and disk areas. Inset: raw data.

APPENDIX A: MARANGONI FORCE AND POWER BALANCE

We examine here the power balance of mechanical forces acting in the system. As we will see, this involves the classical viscous bulk dissipation and the forces acting on the swimmer, together with the rate of change of interfacial energy, that we introduce below to quantify the chemical power available for the system.

The total interfacial energy in the system, and its time derivative-the rate of change in interfacial energy-are defined as [START_REF] Dangla | 2D droplet microfluidics driven by confinement gradients[END_REF][START_REF] Titta | Level-set simulations of a 2D topological rearrangement in a bubble assembly: Effects of surfactant properties[END_REF]

Following Ref. [START_REF] Dangla | 2D droplet microfluidics driven by confinement gradients[END_REF], the rate of change in surface energy Ėγ can be expressed as

The first equality is generic, with the the last two terms corresponding, respectively, to the power of Laplace forces, and to the one developed by Marangoni stresses (or equivalently by viscous stresses at the interface). The second equality is peculiar to our configuration, as we used the assumption of flat interface, that is, of vanishing curvature C. Because ∇ s • u gives the local rate of surface creation, the last expression is a two-dimensional analog to the work of pressure in compressible three-dimensional systems.

Let us now consider the situation of interest. To do so, in the remainder of this Appendix, we switch to the laboratory reference frame where the disk moves at velocity U = U e x and the fluid is at rest at infinity. As already explained in the main text, maintaining steady-state motion for arbitrary values of U generally requires an external force F op to be applied. Because the associated Stokes flow verifies ∇ • σ = 0, the total power of mechanical forces

introducing the classical total rate of viscous dissipation

lead eventually to

The power delivered by the external operator to the system thus decomposes into a bulk viscous term-always dissipative-and a surface term which can act as a power source. Indeed, for spontaneous swimming, F op = 0, so that it is clear that the viscous bulk dissipation is sustained by the chemical power injected in the system in the form of the rate of change of interfacial energy. Of course, the same calculation can be carried for the dual problem, for which the rate of change of interfacial energy vanishes so that

Combined with Eq. (A4) and the force balance condition on the swimmer at steady-state F op + F M + Fv = 0, this yields for the Marangoni force power

Overall, we find that the available power for the Marangoni force acting on the disk amounts to the total available interfacial power diminished by the extra bulk dissipation between the actual and the dual problems.

In addition, comparing with the Lorentz theorem expression for F M [Eq. ( 13)] and the definition of Ėγ [Eq. (A2)], one can obtain an explicit expression for δD v in terms of a rate of change of interfacial energy,

As is seen, it turns out that the fraction of total interfacial power lost as extra bulk dissipation δD v , is exactly the rate of change of interfacial energy due to the extra flow δu.

APPENDIX B: SWIMMING PROPERTIES OF SIMPLIFIED MODELS

As introduced in Sec. III B 2, we consider here models whose basic idea is to decouple the surfactant transport from other processes. In the solute transport equation, the true-and coupledflow field u is replaced by an approximation u app .

Uniform convection approximation a. Concentration field and forces

We start with a toy model where the approximate flow field is uniform u app = -U e x . In this simple framework [START_REF] Boniface | Self-propulsion of symmetric chemically active particles: Point-source model and experiments on camphor disks[END_REF], the solute concentration from a point-source and the capillary force given by Eq. ( 11) are readily obtained as

with I 1 the modified Bessel function of order 1. Using Eq. ( 13) and denoting as G the Meijer's function [START_REF] Gradshteyn | Table of Integrals, Series, and Products[END_REF], the total Marangoni force is

In the regime of small Péclet, the limiting behaviors are In the regime of large Péclet number, a concentration boundary appears, whose thickness results from a balance between advection and diffusion. The case of Stokes flow was treated above. Here we consider the regime of finite Reynolds number. A viscous boundary layer appears, whose thickness is classically given by [START_REF] Guyon | Physical Hydrodynamics[END_REF] δv ∼ Re -1/2 .

(C1)

For real swimmers, the Schmidt number Sc is typically of order 10 3 . The concentration boundary layer is then much thinner than the viscous boundary layer and at the scale of the former, one can approximate the velocity field as varying linearly in z: