M Gauding

F Thiesset

E Varea

L Danaila

Structure of iso-scalar sets

An analytical framework is proposed to explore the structure and kinematics of iso-scalar fields. It is based on a two-point statistical analysis of the phase indicator field which is used to track a given iso-scalar volume. The displacement speed of the iso-surface, i.e. the interface velocity relative to the fluid velocity, is explicitly accounted for, thereby generalizing previous two-point equations dedicated to the phase indicator in two-phase flows. Although this framework applies to many transported quantities, we here focus on passive scalar mixing.

Particular attention is paid on the effect of Reynolds (the Taylor based Reynolds number is varied from 88 to 530) and Schmidt numbers (in the range 0.1 to 1), together with the influence of flow and scalar forcing. It is first found that diffusion in the iso-surface tangential direction is predominant, emphasizing the primordial influence of curvature on the displacement speed. Second, the appropriate normalizing scales for the two-point statistics at either large, intermediate and small scales are revealed and appear to be related to the radius of gyration, the surface density and the standard deviation of mean curvature, respectively. Third, the onset of an intermediate 'scaling range' for the two-point statistics of the phase indicator at sufficiently large Reynolds numbers is observed. The scaling exponent complies with a fractal dimension of 8/3. A scaling range is also observed for the transfer of iso-scalar fields in scale-space whose exponent can be estimated by simple scaling arguments and a recent closure of the Corrsin equation. Fourth, the effects of Reynolds and Schmidt numbers together with flow or scalar forcing on the different terms of the two-point budget are highlighted.

Introduction

There exist a large variety of physical situations in which a description in terms of curved surfaces or interfaces instinctively emerges. Leaving aside some fields of physics such as soft-matter physics [START_REF] De Gennes | Capillarity and wetting phenomena: drops, bubbles, pearls, waves[END_REF], heterogeneous materials [START_REF] Torquato | Random Heterogeneous Materials. Microstructure and Macroscopic Properties[END_REF] or biological/chemical-physics [START_REF] Garcia-Ruiz | Growth patterns in physical sciences and biology[END_REF], : Email address for correspondence: thiesset@coria.fr Abstract must not spill onto p.2 this concept applies naturally to fluid flows. A two-phase flow is the first case that comes to mind since the surface formed by the liquid-gas interface can even be observed with the naked eye [START_REF] Dumouchel | On the experimental investigation on primary atomization of liquid streams[END_REF]. A description in terms of curved surfaces is also widely encountered in reacting flows (diffusion and premixed flames), where the chemical reactions occur in thin layers. This has inspired the flamelet model [START_REF] Peters | Laminar flamelet concepts in turbulent combustion[END_REF]) which considers reacting zones as a collection of thin layers, whose inner structure is identical to a laminar flame, propagating normal to themselves in the direction of the unburned turbulent mixture. In single-phase non-reacting flows, there are situations where a thin interface separates some zones of irrotational motion to some zones of strong vortical intensity [START_REF] Da Silva | Interfacial layers between regions of different turbulence intensity[END_REF]. This layer, which can be observed in many archetypal flow configurations such as wakes, jets or boundary layers, is referred to as the turbulent/non-turbulent interface, abbreviated TNTI. Mixing can also be treated using some geometric measures of iso-scalar surfaces [START_REF] Catrakis | Mixing in turbulent jets: scalar measures and isosurface geometry[END_REF][START_REF] Dimotakis | Turbulence, fractals, and mixing[END_REF]. There are also a variety of natural situations, related to e.g. clouds and precipitations, dunes, coasts erosion, ocean mixing, ice melting, aquifers which can properly be described through a morphological analysis of moving interfaces.

For all such situations, the macroscale features of the interface are of great interest. In two-phase flows, the surface area or surface density (surface area per unit volume) of the liquid-gas interface is generally the parameter one seeks to optimize by resorting to the creation of a spray [START_REF] Ashgriz | Handbook of atomization and sprays: theory and applications[END_REF]. This parameter also controls the evaporation rate in flows with phase change [START_REF] Lebas | Numerical simulation of primary break-up and atomization: Dns and modelling study[END_REF][START_REF] Jay | Combined surface density concepts for dense spray combustion[END_REF]. It is also a key parameter in climate change studies for which the processes taking place at the air-sea interface are primordial [START_REF] Liss | Ocean-atmosphere interactions of gases and particles[END_REF]. In premixed flames, the flame surface area is an important parameter as it appears in the expression of the volume integrated burning rate and heat release (e.g. [START_REF] Trouvé | The evolution equation for the flame surface density in turbulent premixed combustion[END_REF]. For the TNTI, the surface area allows the rate of entrainment of irrotational zones into the turbulent flow to be estimated [START_REF] Sreenivasan | Mixing, entrainment and fractal dimensions of surfaces in turbulent flows[END_REF][START_REF] Krug | The turbulent/non-turbulent interface in an inclined dense gravity current[END_REF].

The versatility of the notion of curved surface finds its foundation on some mathematical grounds. Given any field variable ξ (e.g. temperature, concentration, enstrophy, etc) that varies in space, one can take any iso-value ξ 0 to define an interface which separates the regions where ξpxq ą ξ 0 from the regions where ξpxq ă ξ 0 . The kinematic equations for both the interface position and its geometrical features (surface density, curvatures) are known [START_REF] Pope | The evolution of surfaces in turbulence[END_REF][START_REF] Drew | Evolution of geometric statistics[END_REF][START_REF] Vassilicos | Moving surfaces in turbulent flows[END_REF] thereby embedding in a single mathematical framework single-or two-phase, reacting or non-reacting flows, in presence or absence of phase change.

The wrinkling of the interface is related to intrinsic instabilities and to inhomogeneities (specifically the turbulence) of the carrier environment which itself reveals some multi-scale fluctuations. This means, not only the macroscale features (i.e. measured at scales larger than a typical integral correlation lengthscale) are important, but also the microstructural characteristics (measured at a scale r) are worth being explored. In this respect, it is now well known that interfaces that one may find in turbulence and turbulent mixing [START_REF] Sreenivasan | The fractal facets of turbulence[END_REF][START_REF] Sreenivasan | Mixing, entrainment and fractal dimensions of surfaces in turbulent flows[END_REF][START_REF] Catrakis | Mixing in turbulent jets: scalar measures and isosurface geometry[END_REF][START_REF] De Silva | Multiscale geometry and scaling of the turbulent-nonturbulent interface in high Reynolds number boundary layers[END_REF], turbulent premixed flames [START_REF] Gouldin | An application of fractals to modeling premixed turbulent flames[END_REF][START_REF] Gouldin | Experimental evaluation of the fractal geometry of flamelets[END_REF], two-phase flows [START_REF] Le Moyne | Fractal dimension and scale entropy applications in a spray[END_REF][START_REF] Dumouchel | Application of the scale entropy diffusion model to describe a liquid atomization process[END_REF][START_REF] Grout | Fractal analysis of atomizing liquid flows[END_REF], and material lines evolving in turbulent flows [START_REF] Villermaux | Line dispersion in homogeneous turbulence: stretching, fractal dimensions, and micromixing[END_REF]) reveal some

"fractal facets". The fractal dimension of interfaces was predicted analytically [START_REF] Mandelbrot | On the geometry of homogeneous turbulence, with stress on the fractal dimension of the iso-surfaces of scalars[END_REF][START_REF] Constantin | Fractal geometry of isoscalar surfaces in turbulence: theory and experiments[END_REF][START_REF] Grossmann | Fractal-dimension crossovers in turbulent passive scalar signals[END_REF][START_REF] Iyer | Fractal iso-level sets in high-Reynolds-number scalar turbulence[END_REF]. For instance, [START_REF] Mandelbrot | On the geometry of homogeneous turbulence, with stress on the fractal dimension of the iso-surfaces of scalars[END_REF] proved that the fractal dimension of isoscalars is 2.5 for Burgers turbulence and 8{3 for Kolmogorov turbulence. Using tools from geometric measure theory, [START_REF] Constantin | Fractal geometry of isoscalar surfaces in turbulence: theory and experiments[END_REF] showed that the fractal dimension of iso-scalars might evolve between 7{3 near the TNTI to 8{3 in fully turbulent regions. Later, Constantin (1994b,a) showed that a value of 8{3 holds for iso-scalars in the limit of small molecular diffusivities (high Schmidt numbers) while flame fronts exhibit a fractal dimension of 7{3. The dimensional analysis of [START_REF] Hawkes | A petascale direct numerical simulation study of the modelling of flame wrinkling for largeeddy simulations in intense turbulence[END_REF]; Thiesset et al. (2016a) showed that premixed flame fronts have a fractal dimension of 7{3 (8{3) in low (high) Karlovitz number combustion regimes. The Prandtl (or Schmidt) number dependence of the fractal dimension of iso-scalars is predicted by [START_REF] Grossmann | Fractal-dimension crossovers in turbulent passive scalar signals[END_REF]. The fractal dimension of clouds were also investigated [START_REF] Liss | Ocean-atmosphere interactions of gases and particles[END_REF][START_REF] Hentschel | Relative diffusion in turbulent media: the fractal dimension of clouds[END_REF]. Using numerical data of scalar mixing with an imposed mean gradient at relatively high Reynolds numbers, [START_REF] Iyer | Fractal iso-level sets in high-Reynolds-number scalar turbulence[END_REF] recently showed that the fractal dimension is 2 (8{3) for scalar iso-values far away from (close to) the mean.

Note that we omitted here the possibility that the fractal dimension might be scale-dependent as argued in the review by [START_REF] Dimotakis | Turbulence, fractals, and mixing[END_REF]. Note also that what is here simply referred to as a fractal dimension may recover different mathematical definitions (Hausdorff dimension, Kolmogorov capacity [START_REF] Vassilicos | Fractal dimensions and spectra of interfaces with application to turbulence[END_REF][START_REF] Vassilicos | The multispiral model of turbulence and intermittency[END_REF]).

Characterizing and predicting the microscopic scale-dependent features (= the microstructure) of interfaces requires the coupling between the interface and the surrounding medium to be well understood. In this goal, one needs to identify the range of scales over which some characteristic physical parameters (e.g. surface tension, fluid viscosity, scalar diffusivity, etc) or some physical processes (turbulent straining, production by mean scalar gradient or interface reactivity, etc), have an influence. It is also worth drawing the connections between the typical length-scales of the dynamical or scalar field (integral, Taylor, Corrsin, Kolmogorov, Batchelor, Obukhov length-scales) to those of the interface (inner and outer cutoff scales, radius of curvature, surface density length-scale). All these questions necessitate a scale-by-scale description of the processes at play in the kinematic evolution of contorted iso-surfaces or iso-volumes. To the best of our knowledge, there does not exist such a theoretical framework that may be valid at all scales, irrespectively of the flow configuration and flow regime.

The present study is an attempt to fill this gap. We propose an analytical description that relies on a two-point statistical analysis (correlation and/or structure functions) of the phase indicator function. The latter field variable is used as a 'marker' or 'localizer' of the fluid iso-volume formed by a given iso-scalar value. Such two-point statistics are employed in different branches of physics, generally to gain information about the morphological content (the microstructure) of heterogeneous materials [START_REF] Adler | Flow in simulated porous media[END_REF][START_REF] Torquato | Random Heterogeneous Materials. Microstructure and Macroscopic Properties[END_REF][START_REF] Teubner | Scattering from two-phase random media[END_REF][START_REF] Kirste | Röntgenkleinwinkelstreuung an kolloiden systemen asymptotisches verhalten der streukurven[END_REF][START_REF] Frisch | Contribution to the statistical geometric basis of radiation scattering[END_REF][START_REF] Berryman | Relationship between specific surface area and spatial correlation functions for anisotropic porous media[END_REF] or fractal aggregates [START_REF] Sorensen | Light scattering by fractal aggregates: a review[END_REF][START_REF] Morán | FracVAL: An improved tunable algorithm of cluster-cluster aggregation for generation of fractal structures formed by polydisperse primary particles[END_REF]. In fluid mechanics, there are only few papers dealing with these aspects [START_REF] Hentschel | Relative diffusion in turbulent media: the fractal dimension of clouds[END_REF][START_REF] Vassilicos | Fractal dimensions and spectra of interfaces with application to turbulence[END_REF][START_REF] Vassilicos | The multispiral model of turbulence and intermittency[END_REF][START_REF] Vassilicos | Moving surfaces in turbulent flows[END_REF][START_REF] Lu | Direct numerical simulations of multifluid flows in a vertical channel undergoing topology changes[END_REF], 2019;[START_REF] Elsas | Geometry and scaling laws of excursion and iso-sets of enstrophy and dissipation in isotropic turbulence[END_REF].

Here, the main originality of the present work is that this morphological descriptor is supplemented by an exact transport equation which allows the different physical process acting on the iso-scalar volumes to be characterized. It generalizes the equation proposed by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF] firstly dedicated to two-phase flows, to cases where the interface possesses an intrinsic displacement speed (as for premixed flames or diffusive scalars). The machinery for obtaining two-point statistical equations is the same as the one used to derive the scale-byscale budgets of the dynamical or scalar field (see [START_REF] Hill | Exact second-order structure-function relationships[END_REF]Danaila et al. 2004, among others). We will also resort to some analytical studies emanating from the fields of heterogeneous materials [START_REF] Adler | Flow in simulated porous media[END_REF][START_REF] Torquato | Random Heterogeneous Materials. Microstructure and Macroscopic Properties[END_REF][START_REF] Teubner | Scattering from two-phase random media[END_REF]; [START_REF] Kirste | Röntgenkleinwinkelstreuung an kolloiden systemen asymptotisches verhalten der streukurven[END_REF][START_REF] Frisch | Contribution to the statistical geometric basis of radiation scattering[END_REF][START_REF] Berryman | Relationship between specific surface area and spatial correlation functions for anisotropic porous media[END_REF]) and aggregates [START_REF] Sorensen | Light scattering by fractal aggregates: a review[END_REF][START_REF] Morán | FracVAL: An improved tunable algorithm of cluster-cluster aggregation for generation of fractal structures formed by polydisperse primary particles[END_REF] allowing the phase indicator structure function to be related to some integral geometric measures of the interface (surface density, mean and Gaussian curvatures) and some fractal characteristics. Although the proposed theory may apply to very different situations, we focus here on passive scalar mixing which is explored using Direct Numerical Simulation data covering a wide range of Reynolds and Schmidt numbers. By doing so, we expect emphasizing the key physics that ought to be accounted for e.g. a geometrical closure to the turbulent scalar flux in the equation for the mean scalar.

The present study has four objectives. Firstly, it aims at generalizing the equations firstly derived by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF] to the case of diffusive scalars. The new set of equations reveals the importance of the interface displacement speed which embeds different physics depending on the flow configuration.

The second objective of our work is to characterize the influence of some nondimensional numbers (Reynolds and Schmidt numbers) and some geometrical features (e.g. the mean curvature) on the different components of the interface displacement speed. Thirdly, it aims at identifying the characteristic length-scales, asymptotic scaling and normalizing quantities of the phase indicator structure functions. Fourthly, it intends to explore the effect of Reynolds and Schmidt numbers together with other effects (decay, mean scalar/velocity gradient) on the different processes revealed by the scale-by-scale budgets of the phase indicator field.

The paper is organized as follows. The equations for the transport of the iso-scalar surface and the corresponding phase indicator structure functions are derived in section 2. We also derive the asymptotic limits at either large or small scales revealing the link between two-point statistics of the phase indicator and some integral geometric measures of the iso-surface (volumes, surface area, mean and Gaussian curvature). The numerical database and post-processing procedures are portrayed in section 3. Results are presented in 4. Technicalities are gathered in the Appendix. Conclusions are drawn in a last section.

Analytical considerations

Kinematics of iso-scalar excursion sets

Consider the scalar field ξpx, tq whose transport equation is

B t ξ `∇x ' uξ " ∇ x ' D∇ x ξ `Ξ `9 ω ξ .
(2.1)

Here, u is the fluid velocity, D the scalar diffusivity and 9 ω ξ is the scalar reaction rate. Ξ represents any other source term such as production by a mean gradient.

The iso-scalar ξ 0 forms an interface which separates the zones where ξpx, tq ą ξ 0 and ξpx, tq ă ξ 0 . It is then worth defining the phase indicator function φ "

Focus on Fluids articles must not exceed this page length Hpξ ´ξ0 q, where H denotes the Heaviside function. We then simply have:

φpx, tq " " 1 when ξ ą ξ 0 0 elsewhere.

(2.2) φpx, tq is sometimes referred to as the excursion set of ξpx, tq ą ξ 0 , i.e. the probability that ξpx, tq ą ξ 0 [START_REF] Elsas | Geometry and scaling laws of excursion and iso-sets of enstrophy and dissipation in isotropic turbulence[END_REF]. The transport equation for φpx, tq writes [START_REF] Hirt | Volume of fluid (VOF) method for the dynamics of free boundaries[END_REF][START_REF] Drew | Evolution of geometric statistics[END_REF][START_REF] Vassilicos | Moving surfaces in turbulent flows[END_REF]:

B t φ `u ' ∇ x φ " S d |∇ x φ|.
(2.3)

|∇ x ' | denotes the norm of the gradient of any quantity ' and S d is known as intrinsic displacement speed of the interface. Note that Eq. ( 2.3) is valid only in the sense of distributions [START_REF] Drew | Evolution of geometric statistics[END_REF], i.e. outside from the interface both B t φ and ∇ x φ are zero, while they are equal to the Dirac delta function at the interface.

Similarly, S d is defined only at the surface ξpxq " ξ 0 . Eq. ( 2.3) shows that in the laboratory coordinate system, the observer sees the interface moving at a speed w " u `Sd n where n " ´∇x ξ{|∇ x ξ| is the unit vector normal to the iso-scalar surface. When S d " 0 (as in two-phase flows in absence of phase change), the velocity of the interface w is equal to the fluid velocity at the interface u.

The displacement speed S d is defined by [START_REF] Gibson | Fine structure of scalar fields mixed by turbulence. i. zero-gradient points and minimal gradient[END_REF][START_REF] Pope | The evolution of surfaces in turbulence[END_REF][START_REF] Gran | Negative flame speed in an unsteady 2-D premixed flame: A computational study[END_REF][START_REF] Peters | Statistics of flame displacement speeds from computations of 2-D unsteady methane-air flames[END_REF]:

S d " ∇ x ' D∇ x ξ |∇ x ξ| loooooomoooooon S d d `Ξ |∇ x ξ| loomoon S s d `9 ω ξ |∇ x ξ| loomoon S r d . ( 2
S d d " n ' ∇ x pn ' D∇ξq |∇ x ξ| loooooooooomoooooooooon S n d `2DH lo omo on S c d , (2.5)
where S n d and S c d are the normal and tangential diffusion contributions to the displacement speed. H " ∇ x ' n{2 is the mean curvature of the iso-surface.

It is negative when the surface is concave in the direction of ξpxq ą ξ 0 and convex in the opposite case. The rightmost term in Eq. (2.5) reveals that S c d depends linearly on the mean curvature of the scalar iso-surface. In presence of heat release, it may be more convenient to define S d in a density weighted formulation [START_REF] Giannakopoulos | Consumption and displacement speeds of stretched premixed flames-theory and simulations[END_REF][START_REF] Gran | Negative flame speed in an unsteady 2-D premixed flame: A computational study[END_REF][START_REF] Peters | Statistics of flame displacement speeds from computations of 2-D unsteady methane-air flames[END_REF]. [START_REF] Yu | Evolution equations for the decomposed components of displacement speed in a reactive scalar field[END_REF] derived the transport equations for the different components of the displacement speed.

General two-point equations

The machinery for obtaining the two-point equations of the phase indicator field when S d " 0 is described in details by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF]. Here, we aim at generalizing such equations for cases where the displacement speed is not zero. In this goal, we start by writing the transport equations for φpx, tq at a point x `and x ´, arbitrarily separated in space (Fig. 1). Hereafter, the and ´superscripts denote the quantity at the point x `and x ´, respectively.

Multiplying the equation at x `by φ ´and the one at x ´by φ `, yields:

φ ´Bt φ ``u `' φ ´∇x `φ`" φ ´Sd |∇ x φ| `, (2.6a)
φ `Bt φ ´`u ´' φ `∇x ´φ´" φ `Sd |∇ x φ| ´.
(2.6b)

Since for any quantity r's, we have ∇ x `r's ´" ∇ x ´r's `" 0, one obtains:

φ ´Bt φ ``u `' ∇ x `φ`φ´" φ ´Sd |∇ x φ| `, (2.7a) φ `Bt φ ´`u ´' ∇ x ´φ`φ´" φ `Sd |∇ x φ| ´.
(2.7b)

Summing up these two equations gives,

B t φ `φ´`u`' ∇ x `φ`φ´`u´' ∇ x ´φ`φ2 " φ ´Sd |∇ x φ| ``φ `Sd |∇ x φ| ´.
(2.8)

We now define the mid-point X " px ``x ´q{2 and separation vector r " x `´x (see Fig. 1). Using the relations ∇ x `" 1 2 ∇ X `∇r and ∇ x

´" 1 2 ∇ X ´∇r [START_REF] Hill | Exact second-order structure-function relationships[END_REF][START_REF] Danaila | Progress in studying small-scale turbulence using exact two-point equations[END_REF], we obtain:

B t φ `φ´" ´∇X ' pσuqφ `φ´´∇ r ' pδuqφ `φφ ´Sd |∇ x φ| ``φ `Sd |∇ x φ| 2φ `φ´p σt∇ x ¨uuq
(2.9) Eq. (2.9) is the transport equation for the correlation function of the phase indicator field where S d can take any values. pσ'q " p' ``' ´q{2 and pδ'q " p' `´' ´q. Note that Eq. (2.9) also considers flows in which the velocity divergence may not be zero. This is accounted for in the rightmost term on right-hand side of Eq. (2.9) which reads as the product of φ `φ´a nd the average of ∇ x ¨u between the two-points x `and x ´.

The squared increment of φ, i.e. pδφq 2 " pφ `´φ ´q2 is related to φ `φ´b y φ `φ´" (2.12)

Eq. (2.12) is the general expression for the unaveraged squared increments pδφq 2 . It generalizes the equation derived by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF] to the case where S d ‰ 0 and/or ∇ x ¨u ‰ 0. Compared to the equations detailed by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF], where only the unsteady and the two transfer terms (in r and X space) were present, Eq. (2.12) reveals an additional source term which notably depends on the correlation between a quantity related to the bulk phase φ and a surface quantity, namely S d |∇ x φ|. When estimated numerically, this type of correlation requires a specific treatment which will be described later.

The right-hand side of Eq. (2.12) also contains an additional non-linear forcing term which depends on the velocity divergence. In incompressible flows, this term vanishes.

Eq. (2.12) can be applied to very different types of scalar, either passive, active or reacting, in either decaying or forced turbulence, representative of either single or two-phase flows:

' In two-phase flows with no phase change, S d " 0 and one recovers the equation first derived by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF]. Material surfaces share also the property S d " 0 [START_REF] Pope | The curvature of material surfaces in isotropic turbulence[END_REF]). across the interface that originates from the evaporation rate and the density jump between the two-phase should also be accounted for.

'

' Eq. (2.12) also applies to the enstrophy field. In this situation, S d contains a contribution due to diffusion effects and an additional forcing term due to vortex stretching [START_REF] Krug | The turbulent/non-turbulent interface in an inclined dense gravity current[END_REF]. The present framework is thus likely to help scrutinizing the structure and kinematics of the turbulent/non-turbulent interface which is often defined through a given iso-enstrophy value.

To summarize, the new framework proposed here is very general and enables to treat a variety of different scalars (passive or reacting scalars, with or without forcing) in different flow situations (single or two-phase flows, in forced or decaying turbulence, in presence of phase change).

Because the flows under consideration can be turbulent, it is worth supplementing Eqs. (2.9) and (2.12) by some averaging operators. The choice of a specific average generally depends on the flow situations [START_REF] Hill | Exact second-order structure-function relationships[END_REF]Thiesset et al. 2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF]. One can simply apply an ensemble average operator, noted x'y E , which has the advantage of commuting with time t, spatial X and scale r derivatives.

Hence, the ensemble average of Eq. (2.12) is:

B t xpδφq 2 y E " ´∇X ' xpσuqpδφq 2 y E ´∇r ' xpδuqpδφq 2 y E `2xσtS d |∇φ|uy E ´2pxφ ´Sd |∇ x φ| `yE `xφ `Sd |∇ x φ| ´yE q `2xpδφq 2 pσt∇ x ¨uuqy E .
(2.13)

In the present study, we further exploit the statistical symmetry of the flow (see section 3.2) and we will consider a spatial average over a periodic domain of volume V box :

x'y R " 1 V box ¡ X ' dX.
(2.14) Spatial averages commute with time t and r derivatives, but not with the X divergence operator. However, by periodicity, the fluxes pσuqpδφq 2 normal to the domain boundaries vanish [START_REF] Hill | Exact second-order structure-function relationships[END_REF]Thiesset et al. 2020) where the set of solid angles Ω " tϕ, θ | 0 ď ϕ ď π, 0 ď θ ď 2πu with ϕ " arctanpr y {r x q and θ " arccospr z {|r|q (r x , r y , r z denotes the components of the r vector in x, y, z directions, respectively (2002) who discusses in great details these aspects.

Among this wide corpus of literature, it is worth mentioning the work by [START_REF] Kirste | Röntgenkleinwinkelstreuung an kolloiden systemen asymptotisches verhalten der streukurven[END_REF]; [START_REF] Frisch | Contribution to the statistical geometric basis of radiation scattering[END_REF] who proved that for isotropichomogeneous media, and by further assuming that the interface separating the two phases is of class C 2 , the limit of xpδφq 2 y R at small scales is given by:

lim rÑ0 xpδφq 2 y R " Σ r 2 " 1 ´r2 8 ˆxH 2 y s ´xGy s 3 ˙ .
(2.17)

Here, H and G are the mean and Gaussian curvatures, respectively. x'y s is used to denote the surface area weighted average. [START_REF] Berryman | Relationship between specific surface area and spatial correlation functions for anisotropic porous media[END_REF]; Thiesset et al.

(2021) proved that Eq. (2.17) remains valid in anisotropic media by applying an additional angular average to xpδφq 2 y R . When |r| Ñ 8, Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF] showed that:

lim rÑ8 xpδφq 2 y R " 2xφy R p1 ´xφy R q.
(2.18)

The limit of xφ `|∇ x φ| ´yR as |r| tends to zero can be expressed as follows [START_REF] Teubner | Scattering from two-phase random media[END_REF]:

lim rÑ0 xφ `|∇ x φ| ´yR " Σ 2 " 1 `r 2 xHy s ı . (2.19)
As far as we are aware, the next terms of the small scale expansion of xφ `|∇ x φ| ´yR are not known. In the limit of large separations, we have [START_REF] Teubner | Scattering from two-phase random media[END_REF])

lim rÑ8 xφ `|∇ x φ| ´yR " xφy R Σ.
(2.20) (2.21) apply up to a separation r, which is twice the 'reach' of surface. The 'reach' is a notion that pertains to non-convex bodies. It is defined as the minimal normal distance between the surface and its medial axis (see e.g. [START_REF] Federer | Curvature measures[END_REF]).

The medial axis of a given body is the set of all points having more than one closest point on the object's boundary. It can also be seen as the location of centers of all bi-tangent spheres, i.e. the spheres that are tangent to the surface in at least two points on the surface. The reach is thus given by the minimal radius of these bi-tangent spheres. In some special situations (in absence of narrow throats or necks), it can be related to the minimal radius of curvature.

Otherwise [START_REF] Morán | FracVAL: An improved tunable algorithm of cluster-cluster aggregation for generation of fractal structures formed by polydisperse primary particles[END_REF] showed that when increasing the ratio between the largest scales (the aggregate radius of gyration) and the smallest scales (the radius of the primary particle), the correlation function reveals an increasing range of scales complying with a fractal scaling (a power law). [START_REF] Vassilicos | Fractal dimensions and spectra of interfaces with application to turbulence[END_REF]; [START_REF] Vassilicos | The multispiral model of turbulence and intermittency[END_REF]; [START_REF] Vassilicos | Moving surfaces in turbulent flows[END_REF] showed that xpδφq 2 y E might reveal a power law behavior whose exponent is related to the fractal dimension (more precisely the Kolmogorov capacity) of iso-scalar surfaces. This was investigated in great details by [START_REF] Elsas | Geometry and scaling laws of excursion and iso-sets of enstrophy and dissipation in isotropic turbulence[END_REF] for the enstrophy, dissipation and velocity gradient invariants. When several structures are present, xpδφq 2 y E,R also depends on the way the different fluid structures are organized in space. The reach of the surface plays an important role here since it is the scale which separates the zones in scale space where xpδφq 2 y R depends only on integral geometric measures (Σ, xH 2 y s , xGy s ) and the range of scales for which two-point statistics become a morphological descriptor [START_REF] Torquato | Random Heterogeneous Materials. Microstructure and Macroscopic Properties[END_REF] for which both the geometry and the additional information about the medial axis is required for the structure to be characterized. For scales larger than the reach, the separation r cannot be interpreted as the size of the structure under consideration (as it will be seen later the correlation xφ `φ´y R tends to 0 when the scale r is similar to the size of the structure), but should rather be referred to as the morphological parameter as it is generally done in morphological analysis using e.g. integral geometrical measures (the Minkowski functional) of parallel sets [START_REF] Arns | Characterisation of irregular spatial structures by parallel sets and integral geometric measures[END_REF][START_REF] Dumouchel | Morphology of contorted liquid structures[END_REF].

Given the asymptotic limits detailed in previous section, it seems natural to examine the limit at small and large-scales of Eq. (2.12) (or Eq. (2.9)). In [START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF], it was demonstrated that Eq. (2.12) naturally converge to the transport equation for the surface density when |r| Ñ 0. The latter can be written in the form [START_REF] Pope | The evolution of surfaces in turbulence[END_REF][START_REF] Candel | Flame stretch and the balance equation for the flame area[END_REF][START_REF] Drew | Evolution of geometric statistics[END_REF][START_REF] Blakeley | On the kinematics of scalar iso-surfaces in decaying homogeneous, isotropic turbulence[END_REF] Rapids articles must not exceed this page length as in two-phase flows) and showed that, in the limit of small separations, the X-transport term in Eq. (2.12) asymptotes the convection process of the surface density Σ (the rightmost term on LHS Eq. (2.23)) while the unsteady term in Eq.

(2.12) obviously tends to the unsteady term in Eq. (2.23). [START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF] argued that by difference, the r-transfer term is proportional to the strain rate

K T , viz. lim rÑ0 " ´∇r ' xpδuqpδφq 2 y R ‰ " K T Σ r 2 .
(2.24)

Note that here again, this holds true in anisotropic configuration by using an additional angular average. When the interface displacement speed is not zero, the RHS in Eq. ( 2.12) has the following asymptotic limit:

lim rÑ0 " 2 pσxS d |∇φ|y R q ´2 `xφ ´Sd |∇ x φ| `yR `xφ `Sd |∇ x φ| ´yR ˘‰ " 2xS d y s Σ ´2Σ ´xS d y s `xS d Hy s r 2 " ´2xS d Hy s Σ r 2 " K C Σ r 2 .
(2.25)

Consequently, the additional term in the two-point budget due to the presence of an interface displacement speed asymptotes, in the limit of small separations, to the curvature component of the stretch rate.

Similarly, given that at large scales, the r-transfer term tends to zero, Eq. (2.12) should provide insights into the excursion set volume equation. We indeed obtain that at large scales, the budget simplifies to

lim rÑ8 B t xpδφq 2 y R " 2p1 ´2xφy R qB t xφy R " 2p1 ´2xφy R qxS d y s Σ, (2.26)
where use was made of the equation for the volume [START_REF] Drew | Evolution of geometric statistics[END_REF]):

B t xφy R " xS d y s Σ. (2.27)
Note that the volume of the excursion set xφy R reads as the probability that ξpxq ą ξ 0 . It is thus related to the cumulative distribution of ξpxq. Assuming a Gaussian distribution for the scalar field ξpx, tq, xφy R can be written analytically as

xφy R " 1 2 ˆ1 ´erf ˆξ0 ? 2ξ rms ˙˙(2.28)
where the subscript "rms" stands for the standard deviation of the considered quantity. To dig a little deeper in the interpretation of xpδφq 2 y R , we can resort to some tools from mathematical morphology. This analysis is carried in Appendix A, where we provide a handy demonstration that at intermediate scales, xpδφq 2 y R measures the morphological content of the sets under consideration. The structure function is thus here aptly named since it is a function that depends on the structure (actually the micro-structure) of the fluid elements. It allows some geometrical features of iso-scalar sets to be inferred such as its volume, the surface area, mean and Gaussian curvatures and its transport equation naturally approaches the transport equations of surface density and volume at respectively small and large scales. Note that by virtue of the Gauss-Bonnet theorem, the presence of xGy s in Eq. (2.17) indicates that xpδφq 2 y R depends also on the topology (the Euler characteristic) of the field under consideration.

Consequently, the present set of equations allows not only the morphology of scalar excursion sets to be described, it also accounts for its kinematic evolution through Eqs. (2.9) or (2.12). The interface retains through its geometry and kinematics the signature of the flow dynamics and, in some instances (e.g. twophase flows, TNTI), may even influence the whole flow dynamics. Therefore, we like referring this framework to as a morphodynamical theory since it is likely to provide insights into the morphological evolution of fluid elements.

Numerical database and post-processing

Direct Numerical Simulation of scalar mixing in decaying turbulence

The present analytical framework is appraised using data from Direct Numerical Simulations (DNS). We studied two flow configurations, i.e. forced turbulence (denoted by the "F" letter in Table 1) and decaying turbulence (denoted by the letter "D" in Table 1). We have also one case, noted "T", which was used for tests and validation purposes. Table 1 gathers all important simulation parameters and related statistical quantities, where N denotes the number of grid points along one coordinate axis, ν is the kinematic viscosity, and

R λ " u rms λ ν (3.1)
is the Reynolds number based on the Taylor micro-scale λ " a 15νu 2 rms {x y, where u rms is the root-mean-square velocity, xky " xu i u i y{2 is the mean kinetic energy, and x y " 2νxS ij S ij y is the mean energy dissipation rate, with the strain rate tensor given by S ij " pBu i {Bx j `Bu j {Bx i q{2. Further, xξ 2 y denotes the mean scalar variance, x ξ y " 2DxpBξ{Bx i q 2 y is the mean scalar dissipation rate and λ ξ " a 6Dxξ 2 y{x ξ y denotes the Corrsin length-scale. L t and L ξ denote the integral length scales of the velocity and scalar fields, respectively. The Kolmogorov and the Batchelor length scales are defined as

η " ˆν3 x y ˙1{4 and η B " ˆνD 2 x y ˙1{4 " η Sc 1{2 , (3.2)
respectively. The turbulent Péclet number Pe λ ξ is defined in Section 4.5.

' The forced turbulence database encompasses 6 different values of Taylor based Reynolds numbers R λ (cases F0 to F5 in Table 1) from 88 to 530. For the lowest Reynolds number (F0, R λ " 88), we carried out another simulation with four different scalar fields with different diffusion coefficients D, which correspond to Schmidt number variations from Sc " 0.1 to 1.0. The numerical database for case F0 (Sc " 1.0) to F5 is the same as the one used in [START_REF] Gauding | Line segments in homogeneous scalar turbulence[END_REF][START_REF] Gauding | High-order structure functions for passive scalar fed by a mean gradient[END_REF]. To maintain a statistically steady state, an external stochastic forcing is applied to the velocity field [START_REF] Eswaran | An examination of forcing in direct numerical simulations of turbulence[END_REF]). The forcing is statistically isotropic and limited to low wave-numbers to avoid the small scales to be affected by the forcing scheme. The passive scalar field is fed by a uniform mean scalar gradient G ξ which is applied on the y-direction. Hence, the scalar field ξ can be decomposed into a mean field G ξ y and a fluctuating field ξ, i.e., ξ " G ξ y `ξ .

(3.

3)

The value of the mean scalar gradient G ξ is set to unity without loss of generality.

The indicator function is defined on the fluctuating field ξ, which is statistically homogeneous but not isotropic. The statistical anisotropy that is induced by the mean scalar gradient is further discussed in appendix F. A resolution condition of κ max η ą 2.5 (where κ max is the maximum wavenumber achievable on the numerical grid and η the Kolmogorov length-scale) is maintained for all cases. As a consequence, the number of grid points has been increased to as high as N " 4096 3 for case F5. The statistics presented in Table 1 and throughout the paper correspond to spatial and ensemble averages over M statistically independent snapshots. M varies between 6 for case F5 to 106 for case F0. We have checked that the number of indepedent snapshots M was sufficient for two-point statistics to be converged. Some tests are reported in Appendix G ' For the decaying turbulence case, we explored two distinct situations, the first where the uniform imposed mean scalar gradient is maintained (case D0)

and another where both the velocity and scalar field are decaying (case D1).

For both D0 and D1, we have carried out DNS for two values for the Schmidt number (0.2 and 1.0). The initial velocity field is generated in spectral space to be random and statistically isotropic. It satisfies incompressibility and obeys a prescribed energy spectrum, i.e.,

Epκ, t " 0q9κ 4 exp ˜´2 ˆκ κ p ˙2¸, (3.4)
where κ p is the wave-number at which the initial energy spectrum has its peak.

We chose κ p " 15 as a compromise between limiting confinement effects and the goal of reaching a high Reynolds number. The initial mean kinetic energy xky equals 10 leading to an initial Reynolds number, defined as u rms pt " 0q{pνκ p q, as large as 689. For case D0, the scalar field is initialized to zero, allowing scalar structures to develop naturally from the injection of energy through the imposed mean scalar gradient. For case D1, the scalar field decays freely from a prescribed spectrum, which is identical to the energy spectrum of the velocity field given by Eq. (3.4) with the same initial peak wave-number κ p . Values reported in Table 1 and throughout the paper were obtained at time t " 10, which is about one decade after the onset of the exponential decay of the kinetic energy. At this time, turbulence is highly resolved with a resolution condition κ max η " 12.5. For more details on the setup of the simulations, see [START_REF] Gauding | On the self-similarity of line segments in decaying homogeneous isotropic turbulence[END_REF].

' The test case T0 corresponds to decaying turbulence, but with a mesh size of 512 3 . This allowed us to test the appropriateness of the post-processing procedures. These validations are presented in the Appendix B and C and discussed hereafter in the paper.

The present DNS data were obtained by solving the Navier-Stokes equations and a scalar advection-diffusion equation using a dealiased pseudo-spectral approach. For dealiasing, a filter procedure proposed by [START_REF] Hou | Computing nearly singular solutions using pseudo-spectral methods[END_REF] is used, which ensures stability and inhibits spurious oscillations in real space. For cases F0-F5, a second order semi-implicit Adams-Bashforth/Crank-Nicolson method is used for temporal integration. For the decaying turbulence simulations D0-D1, a low-storage, stability preserving, third-order Runge-Kutta scheme is employed, where for stability, the viscous and diffusive terms are treated by an integrating factor technique. For all cases, the numerical domain is a triply periodic box with length L box " 2π. The simulations have been carried out with an in-house hybrid MPI/OpenMP parallelized simulation code on the supercomputer JUQUEEN at research center Jülich, Germany.

We show some typical snapshots of the ξpxq " 0 iso-surface for different values for the Schmidt and Reynolds numbers in Fig. 2. We show only a 2π ˆ2π ˆπ{2 subset of the simulated domain. The interface is colored by the displacement speed S d which is normalized by the velocity standard deviation u rms " 2xky{3.

The color scale covers the range ´1 ď S d {u rms ď 1. Note that although Fig. 2 gives another impression (remember that only a subset of the domain is presented here), the volume fraction formed by the iso-scalar ξpxq " 0 is the same for all cases and is equal to 0.5. We note that while keeping R λ constant (the three leftmost figures in Fig. 2), a Schmidt number variation from 1 to 0.1 yields a substantial decrease of surface density. The interface is less wrinkled and the corrugation covers a narrower range of scales. This highlights the role of diffusion on the iso-surface geometrical quantities. On the other hand, an increase of the Reynolds number from (the four rightmost figures in Fig. 2) is followed by the creation of smaller and smaller wrinkles and an increase of the morphological content of the iso-scalar volume. xpδφq 2 y R,E is thus believed to widen with R λ .

Fig. 2 also reveals that the displacement speed S d varies mostly in zones of high curvature. This suggests that the curvature of the interface might play a crucial role for understanding the variations of the displacement speed and its different components.

Post-processing procedure

The computation of two-point statistics is challenging as it involves the execution of a convolution operation. We compute two-point statistics accurately in real space by a hybrid MPI/OpenMP parallelization employing the two-dimensional pencil domain decomposition of the DNS code. The partial angular average is approximated by averaging over the r x -, r y -, and r z -directions. Special attention is required for the transfer term, which involves the divergence of a two-point quantity. To avoid the assumption of isotropy, the transfer term is approximated by a second-order finite difference scheme. For instance, in r x -direction, the transfer term reads

∇ r ' xpδuqpδφq 2 y R pr x , 0, 0q « 1 2∆r x " xpδu x qpδφq 2 y R pr x `∆r x , 0, 0q ´xpδu 1 qpδφq 2 y R pr x ´∆r x , 0, 0q ‰ 1 2∆r y " xpδu y qpδφq 2 y R pr x , ∆r y , 0q ´xpδu 2 qpδφq 2 y R pr x , ´∆r y , 0q ‰ 1 2∆r z " xpδu z qpδφq 2 y R pr x , 0, ∆r z q ´xpδu 3 qpδφq 2 y R pr x , 0, ´∆r z q ‰ , (3.5)
where ∆r i is the grid spacing and u x , u y , u z are the velocity component in x, y, and z directions, respectively. The transfer terms in r y -and r z -directions are obtained by a similar procedure.

Some two-point statistics were also computed over the whole r space using the routines available through the increments library of the project pyarcher

(Thiesset & Poux 2020).
Because six nested loops are needed to cover the whole pX, rq-space, we make use of an openMP parallelization for enhancing the calculation speed. Full 3D two-point distributions were estimated only for case F0 and T0 and limited to the range of scales ´96dx ď pr x , r y , r z q ď 96dx where most of the processes take place. By doing so, we are able to check that the partial angular average operated only over r x -, r y -, and r z -directions was leading to similar results than those obtained from the full angular average over whole set of solid angles. These tests are presented in Appendix C and show that the partial angular average yields similar results as the full angular average. In what follows, only the results for the partial angular average will be presented.

Compared to the equation derived by Thiesset et al. (2020), the influence of the interface displacement results in an additional source term in Eq. (2.12).

The latter highlights a correlation of the bulk phase φ with the surface quantity

S d |∇ x φ|.
Hence, this term requires a special treatment. Here, we adapt and develop a procedure inspired by the method of [START_REF] Seaton | Spatial correlation functions from computer simulations[END_REF]. The reader is advised to refer to Appendix B for a description and a validation of the method.

The geometrical properties of the iso-surface (local surface area, mean and Gaussian curvatures, surface conditional statistics) are extracted using the surface operators routines of pyarcher.

Results

Surface conditional statistics

The analytical section presented above reveals that the displacement speed S d contains many of the key physics in the behavior of iso-scalar surfaces. Depending on the situation, it may incorporate different processes such as diffusion, chemical reactions or any other source terms acting on the scalar field. It is thus important to understand how S d evolves along the iso-surface when R λ and/or Sc are varied.

In our situation, the displacement speed contains up to three components, one due to diffusion in the normal direction noted S n d , another due to diffusion in the tangential direction S c d " 2DH and a last contribution associated with the imposed mean scalar gradient noted S s d . Their respective expressions are summarized in Table 2.

As it is generally done in reacting flows [START_REF] Gran | Negative flame speed in an unsteady 2-D premixed flame: A computational study[END_REF][START_REF] Peters | Statistics of flame displacement speeds from computations of 2-D unsteady methane-air flames[END_REF], the variations of S d along the surface ξpx, tq " ξ 0 are analyzed through surface weighted average of its different components conditioned on the mean curvature - -

H.
2 -1 0 1 2 2DH -2 -1 0 1 2 (a) S d |H s S d d |H s S c d |H s -2 -1 0 1 2 2DH -0.2 -0.1 0.0 0.1 0.2 (b) S s d |H s S n d |H s
2 -1 0 1 2 2DH -2 -1 0 1 2 3 (a) S d |H s S d d |H s S c d |H s -2 -1 0 1 2 2DH -0.2 0.0 0.2 0.4 0.6 0.8 (b) S s d |H s S n d |H s

Second-order structure functions

The present framework is first invoked to explore the dynamics of iso-scalars with an imposed mean scalar gradient in stochastically forced turbulence. In this situation, both the flow and scalar characteristics are statistically stationary and the displacement speed has two components, one arising from the diffusive term S d d , and one due to the imposed mean gradient S s d . We focus on the particular effect of the Reynolds and Schmidt numbers.

The Reynolds number dependence of xpδφq 2 y R,E,Ω for different values of the iso-scalar is first considered. Results are presented in Fig. 7 for R λ ranging from 88 to 530 and for two values for the iso-scalar, i.e. ξ 0 " 0 (Fig. 7(a))

and ξ 0 " ξ rms (Fig. 7(b)). xpδφq 2 y R,E,Ω is normalized by its asymptotic value at large scales, i.e. 2xφy E,R p1 ´xφy E,R q while the separation r is normalized by L Σ " 4xφy E,R p1 ´xφy E,R q{Σ. Using the normalization with L Σ " 4xφy E,R p1 ´xφy E,R q{Σ, the asymptotic limits at large and small scales intersect at r{L Σ " 1. Our definition for L Σ finds its inspiration in [START_REF] Lebas | Numerical simulation of primary break-up and atomization: Dns and modelling study[END_REF]; Thiesset et al. (2020). A somehow similar definition for L Σ was conjectured by [START_REF] Peters | A spectral closure for premixed turbulent combustion in the flamelet regime[END_REF] for premixed flames. In the context of the Bray-Moss-Libby model [START_REF] Bray | A unified statistical model of the premixed turbulent flame[END_REF][START_REF] Libby | Implications of the laminar flamelet model in premixed turbulent combustion[END_REF], L Σ is known as the wrinkling scale [START_REF] Kulkarni | Reynolds number scaling of burning rates in spherical turbulent premixed flames[END_REF], although the surface density is defined differently by the latter authors.

10 -2 10 -1 10 0 10 1 r/L Σ 10 -2 10 -1 10 0 R λ (a) F0 F1 F2 F3 F4 F5 10 -2 10 -1 10 0 10 1 r/L Σ
One observes that all curves collapse at small scales where xpδφq 2 y R,E,Ω Ñ Σr{2 and at large scales where xpδφq 2 y R,E,Ω Ñ 2xφy E,R p1 ´xφy E,R q. The most expert readers will probably notice that this behavior is also observed when two-point statistics of the velocity (or scalar) field are normalized by the Taylor (or Corrsin) microscale [START_REF] Thiesset | On self-preservation and log-similarity in a slightly heated axisymmetric mixing layer[END_REF]. Speculatively, this indicates that L Σ plays for [START_REF] Sorensen | Light scattering by fractal aggregates: a review[END_REF][START_REF] Morán | FracVAL: An improved tunable algorithm of cluster-cluster aggregation for generation of fractal structures formed by polydisperse primary particles[END_REF]. Note that what we call here a fractal dimension should rather be identified to a Kolmogorov capacity [START_REF] Vassilicos | Fractal dimensions and spectra of interfaces with application to turbulence[END_REF][START_REF] Vassilicos | The multispiral model of turbulence and intermittency[END_REF]). We also exclude the possibility that the fractal dimension of the intersection of the iso-scalar volume with a line (what we actually measure using two-point statistics of the phase indicator) may be different from the fractal dimension of the iso-volume itself [START_REF] Vassilicos | The multispiral model of turbulence and intermittency[END_REF]).

The numerical value for D f is found to be in the range t2.62 ´2.64u, in quite good agreement with the DNS value reported by Iyer et al. (2020) in the exact same numerical configuration (they find D f " 2.67 at R λ " 650 for ξ 0 " 0 using the box-counting method) and the theoretical analysis of [START_REF] Mandelbrot | On the geometry of homogeneous turbulence, with stress on the fractal dimension of the iso-surfaces of scalars[END_REF] or [START_REF] Grossmann | Fractal-dimension crossovers in turbulent passive scalar signals[END_REF] (also reproduced by [START_REF] Iyer | Fractal iso-level sets in high-Reynolds-number scalar turbulence[END_REF]) providing a value of 8{3. When ξ 0 " ξ rms , the value for the fractal dimension is roughly the same although [START_REF] Iyer | Fractal iso-level sets in high-Reynolds-number scalar turbulence[END_REF] showed that the fractal dimension decreases when the threshold ξ 0 is moved away from the mean value ξ 0 " 0. It is also worth noting that the scale dependence of the local slope in the range L Σ ď r ď 10L Σ does not exceed 5%. This is in contrast with the results presented by [START_REF] Iyer | Fractal iso-level sets in high-Reynolds-number scalar turbulence[END_REF] (see their Fig. 2) where there is no distinct scaling range for ξ 0 " 0. This suggests that, in agreement with the observations of [START_REF] Elsas | Geometry and scaling laws of excursion and iso-sets of enstrophy and dissipation in isotropic turbulence[END_REF], measuring the fractal dimension using the second-order structure function of the iso-scalar field is probably more robust that the one inferred from the box counting method.

The flow under consideration is anisotropic due to the presence of the mean scalar gradient. We here coped with this by employing the partial angular average along three coordinates of r. It is thus worth evaluating if the above features for xpδφq 2 y E,R are retrieved along the different directions of the separation vector r. This point is addressed in Appendix F where similar trends are observed, irrespectively of the direction. Only some small differences between the directions parallel and perpendicular to the mean scalar gradient appear at the large scales.

In appendix D, we also consider normalizing the separation r by the radius of gyration R g . The latter can be computed directly from xpδφq 2 y R,E,Ω by [START_REF] Sorensen | Light scattering by fractal aggregates: a review[END_REF][START_REF] Yon | From monomers to agglomerates: A generalized model for characterizing the morphology of fractal-like clusters[END_REF]:

R 2 g " 1 2 ş 8 0 r 4 Aprqdr ş 8 0 r 2 Aprqdr (4.1)
where Aprq is the correlation function normalized in such a way that Aprq " 1 at r " 0 and Aprq " 0 at large scales, viz.

Aprq " 1 ´xpδφq 2 y R,E,Ω 2xφy E,R p1 ´xφy E,R q (4.2)
The results presented in Fig. 20 of appendix D show that the radius of gyration is a characteristic scale of the distribution xpδφq 2 y R,E,Ω at large scales. Hence, R g plays for xpδφq 2 y R,E the same role as the integral length-scale for normalizing the two-point statistics of the velocity (or scalar) field.

In Appendix E, we also test the appropriateness of using the standard deviation of mean curvature H rms as a similarity variable. This type of normalization is expected to hold at small scales. Indeed, going back to Eq. (2.17), and further assuming xGy s ! xH 2 y s , we have

xpδφq 2 y E,R,Ω ΣH ´1 rms " 1 2 rH rms ˜1 ´prH rms q 2 8 ¸, (4.3)
which is thus expected to be independent of Reynolds and Schmidt numbers when plotted in terms of rH rms . The evolution of xpδφq 2 y E,R,Ω for different R λ when the separation r is normalized by rH rms is presented in Fig. 24 of Appendix E. It appears that rH rms plays for the phase indicator field the same role as the Kolmogorov (or Batchelor) length-scale for normalizing the two-point statistics of the velocity (scalar) field.

The effect of Schmidt number on xpδφq 2 y E,R,Ω is plotted in Fig. 8 for R λ "

10 -2 10 -1 10 0 10 1 r/L Σ 10 -2 10 -1 10 0 Sc (a) Sc = 1.0 Sc = 0.4 Sc = 0.2 Sc = 0.1 10 -2 10 -1 10 0 10 1 r/L Σ Sc (b)
Figure 8: Evolution of xpδφq 2 y E,R,Ω {2xφy E,R p1 ´xφy E,R q with increasing Sc. The separation r is normalized by LΣ " 4xφy E,R p1 ´xφy E,R q{Σ. The dotted gray lines represent the asymptotic theoretical limits at large and small scales. (a) ξ0 " 0, (b) ξ0 " ξrms.

10 -2 10 -1 10 0 10 1 r/L Σ 10 -3 10 -2 10 -1 R λ (a) 10 -2 10 -1 10 0 10 1 r/L Σ R λ (b) 10 -1 10 0 -1 0 1 -0.21 10 -1 10 0 -1 0 1 -0.23
Figure 9: Evolution of the transfer term ´x∇r ¨xpδuqpδφq 2 y E,R yΩ normalized by 2KT xφy E,R p1 ´xφy E,R q with increasing R λ . The separation r is normalized by LΣ. The local scaling exponent is also plotted in the inset. The dotted gray lines represent the asymptotic theoretical limits at small scales. (a) ξ0 " 0, (b) ξ0 " ξrms. The color legend is the same as in Fig. 7.

88 and Sc ranging from 0.1 to 1 and ξ 0 " 0 (Fig. 8(a)) and ξ rms (Fig. 8(b)).

We observe again that normalizing the separation r by L Σ and xpδφq 2 y E,R,Ω by 2xφy E,R p1´xφy E,R q yields a convincing collapse of the different curves at both small and large scales. The Schmidt number effects are perceptible only at intermediate scales where the scale distribution widens with increasing Sc. This means that increasing the diffusivity of the scalar tends to decrease the morphological content of the iso-scalar fields. This is the first evidence that diffusion acts as a restoration effect which counteracts the influence of turbulent straining. More details on this aspect will be given later when examining the budgets of xpδφq 2 y E,R,Ω .

Transfer term

We now address the influence of Reynolds and Schmidt numbers on the transfer term ´x∇ r ¨xpδuqpδφq 2 y E,R y Ω . In Fig. 9, we consider the case where R λ is ranging from 88 to 530 and Sc " 1.0. The transfer term is normalized by 2K T xφy E,R p1 xφy E,R q while the separation is divided by L Σ . This normalization is found to yield a good collapse of all curves at small scales. The specific evolution of the strain rate K T with respect to R λ and Sc will be discussed later. Although visible, the influence of the iso-value ξ 0 is rather limited, at least when the latter is moved from ξ 0 " 0 to ξ 0 " ξ rms .

In Fig. 21 of appendix D, we also report that the transfer term is independent of R λ in the large-scales limit when the separation is normalized by R g while x∇ r xpδuqpδφq 2 y E,R y Ω is divided by a sort of turbulent strain felt at a scale R g which can be written as a xky{R g (xky is the turbulent kinetic energy). Consequently, R g and a xky{R g are appropriate for normalizing the transfer term in the large scale limit. Fig. 25 of appendix E proves that the small-scale similarity variables for the transfer term are H ´1 rms and K T ΣH ´1 rms which plays for φ the same role as the Kolmogorov (Batchelor) scales for the velocity (scalar) field.

In Fig. 9, we also show the local scaling exponent of the transfer term. Although the scaling range appears more restricted than the one observed for secondorder moments, there seems to be a plateau forming around a value of about t´0.21; ´0.23u at the larger R λ . Let us naively assume that, at intermediate scales, the flux can be written as: If we account for internal intermittency, i.e. ζ u ą 2{3, the predicted exponent of the transfer term is closer to the numerical value. Note that this reasoning holds also in the small scale limit, where xpδu || q 2 y 1{2 " xpδφq 2 y " r 1 and hence the transfer term should scale as r 1 , which is numerically observed. To give this scaling analysis a bit more strength, we use the closure proposed by de Divitiis (2014,2016,2020) which has the advantage of not relying on a parametrized turbulent diffusion hypothesis. When adapted to xpδφq 2 y E,R , the latter writes,

xpδuqpδφq 2 y E,R " xpδu || q 2 y 1{2 R,E xpδφq 2 y R,E " r ζu{2`ζ . ( 4
´x∇ r ¨xpδuqpδφq 2 y E,R y Ω " 1 2 xpδu || q 2 y 1{2 E,R,Ω B r xpδφq 2 y E,R,Ω . (4.5)
Assuming again that xpδu || q 2 y E,R y Ω scales as r 2 at small scales and r 2{3 at intermediate scales, we obtain that the transfer term should scale as r 1 and r ζ´2{3 " r 7{3´D f at small and intermediate scales, respectively. This reasoning is in agreement with the numerical data.

While the Kolmogorov four-fifth law and Yaglom four-third law are known to provide a r 0 scaling for the transfer term of either velocity or scalar in the inertial range, the one pertaining to the phase indicator is substantially different and is proved to relate to the fractal dimension of the iso-surface. According to our elaborations, a fractal dimension of 8{3 translates into a r ´1{3 scaling for

10 -2 10 -1 10 0 10 1 r/L Σ 10 -3 10 -2 10 -1 Sc (a) 10 -2 10 -1 10 0 10 1 r/L Σ Sc (b)
Figure 10: Evolution of the transfer term ´x∇r ¨xpδuqpδφq 2 y E,R yΩ normalized by 2KT xφy E,R p1 ´xφy E,R q with increasing Sc. The separation r is normalized by LΣ. The dotted gray lines represent the asymptotic theoretical limits at small scales. (a) ξ0 " 0, (b) ξ0 " ξrms. The color legend is the same as in Fig. 8.

the transfer term of iso-volumes while a fractal dimension of 7{3 results in a r 0 scaling.

The effect of Schmidt number on the transfer term is displayed in Fig. 10(a) for ξ 0 " 0 and Fig. 10(b) for ξ 0 " ξ rms . In both cases, R λ " 88. When Sc increases from 0.1 to 1.0, all curves collapse well at small scales thereby complying with the tL Σ , 2K T xφy E,R p1 ´xφy E,R qu scaling. It is further observed that decreasing the diffusivity of the scalar field (i.e. increasing the Schmidt number) acts in widening the range of scales over which the transfer term operates. The same trend was observed for increasing Reynolds numbers. This suggests that the appropriate non-dimensional number for characterizing the phase indicator scale distribution and its transfer is likely to be the Péclet number. This assertion will be discussed in more details later in this paper.

Two-point budget

The different terms of the angularly averaged budget Eq. (2.15) for Sc " 1.0, 88 ď R λ ď 530 and for two values for the iso-scalar ξ 0 " 0 and ξ 0 " ξ rms are presented in Fig. 11. Here again, the different terms are normalized by 2K T xφy E,R p1 ´xφy E,R q while the separation is divided by L Σ . The normalization by the large-scale quantities R g and a xky{R g is reported in Fig. 22 of Appendix D while the one based on H rms and K T Σ{H rms is plotted in Fig. 26 of Appendix E.

Fig. 11 reveals that the transfer term is positive which means that, as expected, the action of turbulence stirs, stretches and folds the scalar field thereby increasing its tortuousness and its morphological content. The diffusive component of the interface propagation term, i.e. the term due to S d d , is negative and thus acts in smoothing the interface. [START_REF] Peters | A spectral closure for premixed turbulent combustion in the flamelet regime[END_REF] used to refer to the process associated with S d as a kinematic restoration effect which appears indeed aptly named as it tends to counteract the influence of turbulent strain by smoothing the interface. The term in the two-point budget associated with the imposed mean gradient, i.e. the term due to S s d , is negative for an iso-scalar value ξ 0 " 0 and positive for ξ 0 " ξ rms . This indicates that the imposed mean gradient decreases the morphological content of the iso-scalar close to the ξ 0 " 0 iso-value and redistributes it to the scalar iso-values away from the mean.

10 -1 10 0 10 1 r/L Σ -0.2 -0.1 0.0 0.1 0.2 R λ (a) Transfer S d d -term S s d -term 10 -1 10 0 10 1 r/L Σ R λ (b)
The influence of R λ is also visible in Fig. 11. When normalized by 2K T xφy E,R p1829 xφy E,R q, the different terms collapse at small scales and decrease in amplitude in the intermediate range of scales. When R λ increases, the range of scales over which the different terms of the budget are contributing, increases. Using the large-scale normalization (see Fig. 22 in appendix D), the terms collapse at large scales and their respective amplitude increases with R λ . The evolution of the different terms of the budget using the small-scale similarity variables is presented in Fig. 26 of Appendix E. It is also worth noting that, at small, up to intermediate scales, the relative influence of the S s d -term compared to e.g. the transfer term decreases when R λ increases. This suggests that in the limit of very large R λ , the strain rate and the diffusion terms balance while the source term due to the imposed mean gradient remains concentrated only at the large scales. In other words, at large The effect of Schmidt number on the different terms of the budget at constant R λ " 88 is displayed in Fig. 12(a) for ξ 0 " 0 and Fig. 12(b) for ξ 0 " ξ rms . Here again, the normalization using L Σ and K T yields a good collapse of all curves in the limit of small separations. For ξ 0 " 0, decreasing the Schmidt number from 1.0 to 0.1, i.e. increasing the scalar diffusivity by a factor of 10, leads to a smaller amplitude of the diffusion term in the intermediate range of scales. In the same range of scales (i.e. up to r « L Σ ), the transfer term normalized by K T appears rather insensitive to Schmidt number variations. The influence of Sc on the transfer term is perceptible only at large scales where it is observed that the scale at which the transfer term approaches zero decreases with increasing scalar diffusivity. The third term in the budget due to the imposed mean gradient is plotted as dash-dotted lines. When ξ 0 " 0, the latter is negative, progressively tends to zero when Sc decreases, and becomes even slightly positive for Sc " 0.1. This means that when Sc " 0.1, the scalar is so diffusive that the imposed mean gradient becomes a gain in the budget for this particular iso-value. For ξ 0 " ξ rms , the term due to the imposed mean gradient acts at rather large scales, is positive and increases in amplitude when Sc decreases. At large scales the budget is composed of only the diffusion and imposed mean gradient terms, while the transfer term is zero.

Strain and curvature components of the stretch rate

We now investigate the evolution of the different normalizing quantities with respect to R λ and Sc. When both R λ and Sc vary, it may be more appropriate to define the turbulent Péclet number:

P e λ ξ " a 2xky{3λ ξ D (4.6)
With this definition, the Péclet number is related to the Schmidt and Reynolds number by

P e λ ξ " ˆ6R 10 ˙1{2 R λ Sc 1{2 (4.7)
where R is the scalar to mechanical time-scale ratio, i.e. Results are portrayed in Fig. 13 where all quantities are made non-dimensional by multiplying by the Kolmogorov time scale τ η " pν{x yq 1{2 .

For P e λ ξ ą 50, the normalized strain rate K T τ η is nearly constant around a value of about 0.25. Since this observation holds for both ξ 0 " 0 and ξ 0 " ξ rms , this means that the different iso-scalars experience nearly the same turbulent straining. The value of 0.25τ η is in agreement with the finding of Yeung et al.

(1990) who report a value of 0.28 for material surfaces. It is also consistent with the phenomenological model of Thiesset et al. (2016b) which gives roughly the same value of 0.28. It is argued by [START_REF] Girimaji | Propagating surfaces in isotropic turbulence[END_REF] that the strain rate experienced by propagating surfaces is smaller than that acting on material surfaces since there will be less time for the iso-surface to align with strain. Our value of 0.25 for K T τ η instead of 0.28 is thus consistent with this argument. In addition, when the scalar diffusivity is increased so that P e λ ξ ă 50, we note a substantial decrease of K T τ η which drops down to 0.2. This indicates that the higher the diffusivity, the larger is the displacement speed, and the smaller the time for the iso-surface to align with strain. The constancy of K T τ η at large Péclet numbers is also predicted by the closure of de Divitiis (2014,2016,2020) given by Eq. (4.5). Indeed, Kolmogorov's first similarity hypothesis implies for r Ñ 0:

xpδu || q 2 y E,R,Ω " x y ν r 2 , (4.9)
which by virtue of Eqs. (2.17) and (2.24) gives K T τ η " const..

A careful analysis of Fig. 13 further shows that K d C τ η is always negative and thus counteracts the effect of turbulent straining. It also approaches a constant value when Péclet increases, but at a smaller rate than K T . The influence of the imposed mean gradient on the source component of the stretch rate K s C τ η is perceptible at finite Péclet number and it is observed to be negative for ξ 0 " 0 and positive for ξ 0 " ξ rms . Here again, this suggests that the imposed mean gradient acts in redistributing the surface density from the mean iso-value ξ 0 " 0 to isovalues away from the mean. K s C τ η approaches zero at the highest Péclet number. This indicates that in the limit of very high P e λ ξ , the imposed mean gradient does not influence the evolution of the iso-scalar surface density, the latter being driven only by diffusion and straining effects.

Characteristic length-scales

The previous analysis of two-point statistics highlighted the existence of three characteristic length-scales for the excursion set φ. The first one, L Σ , is relevant for normalizing the two-point statistics at both small and large scales. The second R g applies in the large scale limit while the third H ´1 rms is relevant at small up to intermediate scales. Some speculations about the connection between L Σ (H ´1 rms )

and the Corrsin (Batchelor) microscales have already been stated earlier in this paper. We now provide more rigorous evidence for this.

First, it is worth recalling that, for isotropic media, the surface density Σ is related to the number of zero-crossings of the field ξpxq ´ξ0 [START_REF] Torquato | Random Heterogeneous Materials. Microstructure and Macroscopic Properties[END_REF]. On the other hand, there exist numerous studies that highlight a close multiplied by roughly 20, the ratio R g {L ξ and λ ξ {L Σ remain close to 0.5 with some small departures which do not exceed ˘20%. Hence, it appears that L Σ " λ ξ and R g " L ξ are good approximations. The small departures in the scaling between R g and L ξ can be due to several effects. First, some confinement due to the finite ratio between L ξ and the simulation box size can be at play. Second, the scaling between R g and L ξ can also be altered by some finite Péclet numbers effects. The latter are likely to reveal themselves in the evolution of R g {L ξ which is first decreasing before reaching a plateau for Péclet numbers above 50. Such finite Péclet number effects were also noticed by [START_REF] Shete | Area of scalar isosurfaces in homogeneous isotropic turbulence as a function of reynolds and schmidt numbers[END_REF]. As far as the scaling between L Σ and λ ξ is concerned, it is first worth recalling that the Rice's theorem is valid only if both ξ and ∇ξ have Gaussian probability density functions and are statistically uncorrelated. While ξ is normally distributed, internal intermittency leads to significant departure from Gaussian distributions for the scalar gradient.

K T τ η K d C τ η K s C τ η
At finite Reynolds or Péclet numbers, the assumption of statistical independence between ξ and its derivatives is not likely to hold. As a consequence, all data from the literature (e.g. [START_REF] Sreenivasan | Zero-crossings in turbulent signals[END_REF] indicate that the Rice theorem is a good approximation for turbulent signals verified within 20%. Departures are thus of same magnitude here. Finally, we explored here only two scalar iso-values.

It is not excluded that the proportionality between L Σ and λ ξ ceases to apply for some higher iso-values of ξ 0 . Let us now focus on the scaling of H rms . It is first worth recalling that the relation K T « K d C " τ η holds relatively well, except maybe at the smallest Péclet numbers (Fig. 13). On the other hand, the analysis performed in section 4.1

indicates that S d « 2DH is a rather safe approximation. With these relations, one can easily conclude that the standard deviation of the mean curvature is proportional to the inverse of the Batchelor length-scale η B . This result is tested with success in Fig. 14(a) which proves that the standard deviation of mean curvature is indeed proportional to the Batchelor length-scale (Due to computational limitations, we were unable to estimate H rms for case F5.

Hereafter, we will assume H rms " 0.22{η B for this case). It is worth noting that DNS data give η B H rms " const. even at low Péclet numbers where the approximation K T « K d C " τ η does not hold anymore. We thus believe that the relation between the Batchelor length-scale and the standard deviation of mean curvature is likely to be a general result which does not necessarily requires

K T « K d C " τ η .
We have also tested the scaling of H rms with the Obukhov-Corrsin length-scale η OC " pD 3 {x yq 1{4 " η{Sc 3{4 (not shown), but the latter was found to be inappropriate. This conclusion is consistent with the findings of e.g. [START_REF] Antonia | Effect of Schmidt number on small-scale passive scalar turbulence[END_REF]; [START_REF] Donzis | Scalar dissipation rate and dissipative anomaly in isotropic turbulence[END_REF] which prove that the Batchelor lengthscale is the appropriate normalizing the scalar fluctuations scale distributions irrespective of the Schmidt number.

To conclude, for the values of ξ 0 investigated here, the geometrical characteristics scales R g , L Σ and H ´1 rms of the iso-scalar fields can be related to the somehow more 'usual' characteristic scales of the scalar field: the integral, Corrsin and Batchelor length-scales, respectively. There thus exists an intimate connection between the geometrical and hydrodynamic characteristic scales.

In the appendix, we plot the local scaling exponent of xpδφq 2 y R,E,Ω when the separation is made non-dimensional using either R g (Fig. 20) or H ´1 rms (Fig. 24).

We observe in Fig. 24 that the onset of the fractal scaling range starts at a scale proportional to H ´1 rms . On the other side (Fig. 20), the end of the scaling range appears at a scale proportional to R g . Therefore, the inner cutoff of the fractal range is related to H ´1 rms (and hence η B ) while the outer cutoff is given by R g (and hence L ξ ). Note that [START_REF] Sreenivasan | Mixing, entrainment and fractal dimensions of surfaces in turbulent flows[END_REF] conjectured that for small Schmidt number, the inner cutoff should be related to the Obukhov-Corrsin length-scale η OC . Our findings indicate that the inner cutoff should better be scaled with the Batchelor length-scale.

A fractal scaling should result in [START_REF] Sreenivasan | Mixing, entrainment and fractal dimensions of surfaces in turbulent flows[END_REF],

L box Σ " k f ˆRg H ´1 rms ˙Df ´2 , (4.11)
where k f is the fractal prefactor and L box " 2π is used for normalization as in [START_REF] Shete | Area of scalar isosurfaces in homogeneous isotropic turbulence as a function of reynolds and schmidt numbers[END_REF]. Figure 14(b) portrays the evolution of surface density Σ with respect to R g H rms for all Péclet numbers. The log-log representation clearly indicates that a power-law is at play with an exponent in very close agreement with D f " 2.62 and D f " 2.64 for ξ 0 " 0 and ξ 0 " ξ rms respectively. Surprisingly, Eq. (4.11) holds even at the lowest Péclet numbers, although there is no clear scaling range for xpδφq 2 y R,E,Ω . We also note that the fractal prefactor k f in Eq. (4.11) depends on particular the choice of the iso-scalar but not on the Péclet number.

Expressing R g and H rms in terms of L ξ and η B in Eq. (4.11) yields:

Σ " ˆC L ξ L t P e 3 2 λ Sc ´1˙D f ´2 (4.12)
where C is the kinetic energy dissipation constant and P e λ is the Taylor based Péclet number P e λ " R λ Sc. When expressed in terms of P e λ ξ , Eq. (4.12) writes:

Σ " ˆC L ξ L t R ´3{4 P e 3 2 λ ξ Sc ´1{4 ˙Df ´2 (4.13)
Eq. (4.12) indicates that for Sc " 1, and further omitting the dependence of C and L ξ {L t to R λ , the surface density of iso-scalars Σ should grow as P e . If D f " 8{3, we have Σ " P e λ " P e λ ξ " R λ while

D f " 7{3 leads to Σ " P e 1{2 λ " P e 1{2 λ ξ " R 1{2 λ .
The P e λ scaling derived in Eq. (4.12) is different from the one observed by [START_REF] Shete | Area of scalar isosurfaces in homogeneous isotropic turbulence as a function of reynolds and schmidt numbers[END_REF] in a configuration similar to the present one (although the velocity forcing and dealiasing procedures were different). They obtained that Σ " P e 1{2 λ over an impressive range of Péclet numbers. This result was obtained by averaging the surface areas over 20 different iso-levels covering the range of scalar fluctuations.

Since the fractal dimension D f is known to vary for different iso-levels [START_REF] Iyer | Fractal iso-level sets in high-Reynolds-number scalar turbulence[END_REF], lumping together the values for 20 iso-surface areas could have resulted in a different scaling with respect to the Péclet number. Further investigations are required to confirm whether or not Eq. (4.12) holds true for the data of Shete & de Bruyn Kops (2020).

Decaying turbulence

We now proceed to the analysis of the scale-by-scale budgets in decaying turbulence. In this situation, the time derivative term in Eq. (2.15) contributes to the budget. We consider two cases where the production term associated with the imposed mean scalar gradient is either retained or deactivated. For each situation, a Schmidt number of 1.0 and 0.2 is analyzed.

We start with the case where the kinetic energy is freely decaying but the imposed mean scalar gradient is maintained. Therefore, the time derivative term from 0 to ξ rms . We observe that when normalized by 2K T xφy E,R p1 ´xφy E,R q and L Σ , the transfer term is almost independent of the choice of the iso-scalar value.

10 -2 10 -1 10 0 10 1 r/L Σ -0.3 -0.2 -0.1 0.0 0.1 0.2 0.3 ξ 0 (a) Transfer S d d -term S s d -term dt-term 10 -2 10 -1 10 0 10 1 r/L Σ ξ 0 (b)
This means that the different iso-surfaces experience the same scale dependent turbulent straining. The Schmidt number variations are similar to the one already documented in the previous section. The contribution due to scalar diffusion is always negative which means that the restoration effect acts in counteracting the turbulent straining. At large scales, increasing ξ 0 and decreasing Sc is followed by a increasing amplitude of the S d d -term which is consistent with an increase of S d and D. The term associated with the imposed mean scalar gradient, is positive, acts at larger scales and its amplitude increases with ξ 0 . The balance between all these terms yields the time variations of xpδφq 2 y E,RΩ . This term has a rather small amplitude compared to the three others. We note however that at large scales and for ξ 0 ą 0, the time variations of xpδφq 2 y E,RΩ are positive meaning that the volume xφy E,R tends to increase with time. This is consistent with the observation that in decaying turbulence in presence of an imposed mean scalar gradient, the variance of ξ grows in time. We also note that at large scales, all terms of the budget tend to zero when ξ 0 " 0 meaning that the volume xφy E,R for ξ 0 " 0 is conserved.

Deactivating the source term leads to S s d " 0. In this situation, the time variations of xpδφq 2 y E,R,Ω is due to the unbalance between the transfer term and diffusion term. Results for Sc " 1.0 and 0.2 and for ξ 0 " 0, 0.5ξ rms and ξ rms are presented in Fig. 16. Here again, we observe that the particular choice of ξ 0 does not change drastically the transfer term when the latter is scaled in terms of from 0 to ξ rms . We also observe that the slope of the time derivative term at small scales is negative meaning that the surface density is decreasing in time.

At large scales, the unsteady term is negative, asymptotes the diffusion term and increases in amplitude with ξ 0 and Sc.

In summary, forced and decaying turbulence share some common behaviors:

' The transfer term is positive in both case, meaning that turbulent straining acts in increasing the morphological content of iso-scalar sets.

' The contribution due to diffusion (the S d d -term) is always negative. It thus plays a restoration effect that counteracts the effect of turbulent straining.

' Increasing scalar diffusivity (decreasing Sc) always leads to a decrease of the iso-scalar morphological content. This is mainly due to the above-mentioned restoration effect, where the latter increases with D, which is consistent with the assumption S d d « 2DH.

In contrast, there are key differences between forced and isotropic turbulence which are summarized below:

' In forced stationary turbulence, the time derivative term in Eq. (2.15) is by definition zero for all scales.

' In decaying situations, depending on the iso-scalar value and the presence or absence of mean scalar gradient, the time derivative term can be either positive or negative as detailed below:

˝In decaying turbulence, in absence of mean scalar gradient, the time derivative term is systematically negative, meaning that both the surface density and iso-scalar volume are decreasing.

˝In presence of mean scalar gradient, the time derivative term can be either positive or negative depending on the iso-scalar value.

A new theory is proposed to characterize the time evolution of iso-scalar volumes.

It is based on the two-point transport equation of the phase indicator field.

The main analytical tools that are convoked emanate from two, apparently disconnected, fields of physics. On the one hand, using known analytical results from the field of heterogeneous media and fractal aggregates, we have shown that two-point statistics of the phase indicator allow some integral geometric quantities (volume, surface area and curvatures) and some morphological characteristics (reach, inner/outer cutoff, fractal dimension) to be measured. On the other hand, we invoked the machinery of the scale-by-scale budgets which is adapted to the kinematic equation of iso-volumes. Combining such two approaches allows not only the geometry, morphology and topology of the fluid structures to be assessed; it also embeds their scale/space/time evolution. As a consequence, we like referring to this framework as a morphodynamical theory. It also naturally degenerates to the transport equations for the volume and surface density in the limit or large and small scales, respectively, thereby offering promising perspectives for modelling either scalar mixing [START_REF] Catrakis | Mixing in turbulent jets: scalar measures and isosurface geometry[END_REF], two-phase flows [START_REF] Lebas | Numerical simulation of primary break-up and atomization: Dns and modelling study[END_REF] or combustive flows [START_REF] Trouvé | The evolution equation for the flame surface density in turbulent premixed combustion[END_REF]) using a

volume-surface density approach. The new set of equations derived in the present work generalizes some previous analysis by Thiesset et al. (2020[START_REF] Thiesset | Space-scale-time dynamics of liquidgas shear flow[END_REF]. In the present work, light is shed on scalar mixing in either forced or decaying turbulence using state-of-the-art DNS data covering a large range of Reynolds and Schmidt numbers. We paid attention to the correlation between the different components of the displacement speed and the mean curvature of the interface.

It is shown that the tangential diffusion contribution dominates, meaning that, as a first approximation, S d is proportional to the scalar diffusivity D and the mean curvature H. This is a result of major importance which proves that there exists an intimate relation between the geometry of interface and some dynamical processes such as diffusion. Further, the geometry of the interface at a microscale (i.e. at a scale r) has an influence on some macroscopic (i.e. when r Ñ 8)

processes such as the conservation of scalar iso-volumes.

The search for the appropriate similarity variables have shown that there exists three important characteristic length-scales for the second-order structure function of the phase indicator field.

' The first one corresponds to the inverse of the mean curvature standard deviation H ´1 rms . This scale together with ΣK T are the similarity variables at small scales up to intermediate scales. The existence of this normalizing scale is justified by the small-scale expansion of xpδφq 2 y E,R,Ω (Eq. (2.17)). We also discovered that the scale beyond which a fractal scaling range starts to appear is proportional to H ´1 rms suggesting that the inner cutoff is related H ´1 rms . It was also observed that H rms η B " const. over the range of Péclet numbers investigated here. An explanation for this observation when K T " K C " τ ´1 η is provided. However, the proportionality between H ´1 rms and η B is a likely more general result which holds even at low Péclet numbers. This means that the assumption K T " K C " τ ´1 η is not likely to be a necessary condition.

' The second set of normalizing scales can be expressed in terms of surface density Σ and volume xφy R,E . They arise very naturally from the small-scale expansion (Eq. (2.17)) and the large scale limit of xpδφq 2 y E,R,Ω , respectively. It was previously identified by Thiesset et al. (2020). These quantities are appropriate for normalizing two-point statistics in the limit of either small and large scale.

It was found that, for the range of ξ 0 investigated here, L Σ is proportional to the Corrsin length-scale, in agreement with previous analysis based on the Rice's theorem.

' The last set of normalizing scales is provided by R g , the radius of gyration and the strain felt at scale R g , viz. xky 1{2 {R g . These are the characteristic quantities of the large scale processes. It is found that R g is proportional to the integral length-scale of the scalar field L ξ . The local scaling exponent of xφy R,E

was further shown to depart from a constant scaling exponent at a scale similar to R g . Hence, the outer cutoff is related to R g and thus L ξ .

At sufficiently large R λ , the distribution xpδφq 2 y E,R,Ω was shown to behave according to a power law in the intermediate range of scales. The corresponding scaling exponent is related to the fractal dimension that is found to be close to 8{3 in agreement with the theoretical analysis of [START_REF] Constantin | Fractal geometry of isoscalar surfaces in turbulence: theory and experiments[END_REF]; [START_REF] Mandelbrot | On the geometry of homogeneous turbulence, with stress on the fractal dimension of the iso-surfaces of scalars[END_REF]; [START_REF] Grossmann | Fractal-dimension crossovers in turbulent passive scalar signals[END_REF]. The fractal dimension together with the inner and outer cutoff allows the surface area of the iso-scalars to be estimated. Surprisingly, it was observed that this applies even at low Péclet where a fractal scaling range is not likely to hold. The transfer term was also shown to possess a scaling range, and we provided the value for the scaling exponent by resorting to the closure proposed by de Divitiis (2014,2016,2020).

The effect of Reynolds and Schmidt number on the different contributions of the stretch rate, viz. the strain rate K T , the curvature term associated with diffusive effect K d C and the curvature term associated with the forcing K s C is explored. It is shown that in forced turbulence with an imposed mean scalar gradient, K T is positive and compensated by both K d C , which is systematically negative, and K s C whose sign depends on the iso-scalar value ξ 0 . In the limit of large R λ , K s C Ñ 0 which proves that the geometry of the interface at the smallest scales tends to be independent of the type of forcing. The closure of de [START_REF] De Divitiis | Finite scale Lyapunov analysis of temperature fluctuations in homogeneous isotropic turbulence[END_REF][START_REF] De Divitiis | von Kármán-Howarth and Corrsin equations closure based on Lagrangian description of the fluid motion[END_REF][START_REF] De Divitiis | von Kármán-Howarth and Corrsin equations closures through Liouville theorem[END_REF] was also invoked to prove that K T is proportional to τ ´1 η .

Finally, we examined the scale distribution of the different terms of the scaleby-scale budget, see Eq. (2.15), for different values of R λ and Sc, in either forced or decaying turbulence. It was shown that the transfer term, which measures the interaction between the interface and velocity field at scale r, is systematically positive. This means that turbulence acts in increasing the morphological content of the interface. On the other hand, the term associated with the diffusive component of the displacement speed is always negative, meaning that diffusion acts in counteracting turbulent straining through a so-called kinematic restoration effect. Although this conclusion appears rather intuitive, it is here significantly strengthened by a quantitative and analytical framework based on two-point statistical equations. The last term is due to the forcing imposed by a mean scalar gradient. The latter can be either positive or negative depending on the iso-scalar value. This term tends to zero in the small and intermediate range of scales when R λ increases and its contribution is progressively pushed towards the largest scales. This proves again that the geometry of the interface tends to be independent of the type of forcing at sufficiently large R λ . We also explored the case of decaying turbulence, with and without the imposed scalar gradient. We showed that the transfer term remains roughly independent to the iso-scalar value meaning that the different iso-scalars experience the same turbulent straining.

The time-evolution of the phase indicator structure function is thus given by the balance between the S d d and S s d terms. Without the imposed mean gradient, the unsteady term was found negative meaning that both the surface density, morphological content and volume are decreasing during the decay.

The present work opens up attractive perspectives. First, given the amount of computational time needed to post-process this database, we have restricted ourselves to a limited range of iso-scalar values. A more systematic study of the evolution of xpδφq 2 y E,R,Ω and the different terms of the budget with respect to ξ 0 is now needed. Second, we have explored the effect of Reynolds and Schmidt numbers independently. However, our database does not address the case of varying Sc and R λ while keeping the Péclet number P e λ ξ constant. This would allow one to conclude about the similarity of two-point statistics with respect to the Péclet number. Third, the flow configuration explored here is statistically homogeneous and hence, the scale-by-scale budgets are independent on the flow position. The next step for addressing some complex flow configurations will thus consist in better characterizing the effect of inhomogeneities and anisotropy.

We have focused here only on Schmidt numbers Sc ď 1. For very high Schmidt numbers, we expect results to be quite different. Indeed, for Sc " 1 and R λ " 1, scalar spectra are expected to reveal two distinct scaling ranges: the inertialconvective scaling range with an exponent close to ´5{3 which extends up to the Kolmogorov scale followed by a viscous-convective scaling range with an exponent of ´1 which ends at the Batchelor scale. Consequently, for high Schmidt numbers, we expect the structure function of φ to reveal two distinct scaling ranges: the first with an exponent of 3 ´Df,1 which corresponds to the inertial-convective range, preceded by a scaling range with an exponent 3 ´Df,2 corresponding to the viscous-convective range. At very small-scales, the local scaling exponent of xpδφq 2 y R,Ω should approach 1 as shown by Eq. (2.17). Therefore, it would be of great interest in this context to measure the values for D f,1 and D f,2 together with their respective inner and outer cutoff length-scales using the present framework.

However, numerical data at high Schmidt numbers are particularly challenging to obtain since one faces numerical issues to achieve both high Reynolds (for the inertial-convective range to establish) and high Schmidt numbers (for the diffusive-convective range to be sufficiently large). We hope that in the future, such a dataset will be available so that to test the ability of the present framework to infer the scaling and the terms of the transport equation of xpδφq 2 y R,E .

Finally, besides turbulent mixing, the present mathematical framework should now be harnessed for giving insights into the physics of other type of interfaces such as reacting fronts, turbulent/non-turbulent layers, or two-phase flows in presence of evaporation.

green sets in Fig. 17. Given Eq. (2.10), the spatially averaged structure function xpδφq 2 y R is thus given by the volume of φpxq, plus the volume of φpx `rq, minus two times the intersection. It thus reads as the disjunctive union (or symmetric difference) of φ ´and φ `which is graphically represented by the orange sets in Fig. 17.

For small values of the separation vector, one sees that the orange set in Fig. 17 delineates the contours of the excursion set while for large scales, the disjunctive union writes as twice its volume. For intermediate scales, the correlation and structure function might depend on the morphology of the structures under hand.

Therefore, the graphical representation of xpδφq 2 y R in Fig. 17 The additional source term in Eq. (2.12) highlights the correlation of φ with S d |∇ x φ|. Before describing the procedure we employed for computing this term, we first recall available methods to infer numerically the surface-bulk correlation xφ ´|∇ x φ| `yR . One method is presented in [START_REF] Seaton | Spatial correlation functions from computer simulations[END_REF] which is also briefly described in [START_REF] Ma | Precise algorithms to compute surface correlation functions of two-phase heterogeneous media and their applications[END_REF]. It consists in computing the correlation function xφ ` φ ´ y R pr, q of the fields φ pxq which denotes either the dilated (when the scale ą 0) or eroded (when the scale ă 0) version of φpxq.

Here the eroded/dilated objects can simply be defined from the excursion set, i.e.

φpx, q " HpΥ pxq ´ q, where Υ pxq is the level-set field of the iso-scalar under consideration. Then, the surface-bulk correlation function can be proven to be equal to

xφ ´|∇ x φ| `yR " lim Ñ0 1 2 xφ ` φ ´ y R pr, q (B 1)
Numerically, this translates into

xφ ´|∇ x φ| `yR " xφ ` φ ´ y R pr, q ´xφ ` φ ´ y R pr, ´ q 4 , ( B 2)
where should be chosen sufficiently small for the eroded/dilated sets remain topologically equivalent to the actual set. In practice, we found that when remains in the range 2dx ď ď 4dx, results are very similar and the surface-bulk correlation function is in very good agreement with the theoretical limits at large and small scales.

The method for computing xφ ´Sd |∇ x φ| `yR is somehow similar. It relies on the idea that one can define a local dilatation/erosion scale that depends on local values of S d . At a fictive time t`τ , the level-set field will be given by Υ px, τ q`S d τ while the one obtained at time t ´τ is Υ px, ´τ q ´Sd τ . S d τ plays here the same role as which thus corresponds to the special case where S d " const. Then, one can compute xφ ´Sd |∇ x φ| `yR as d -term as estimated from the level-set or its approximation using Eq. (B 5). Case T0 Sc " 1.0, ξ0 " ξrms. In (a) and (b), the S i d -terms are normalized by KT and the grey dotted lines represent the theoretical asymptotic limits at large and small scales.

xφ ´Sd |∇ x φ| `yR " lim τ Ñ0 1 2 xφ τ φ τ y R pr, τ q τ (B 3)
or numerically

xφ ´Sd |∇ x φ| `yR " xφ τ φ τ y R pr, τ q ´xφ τ φ τ y R pr, ´τ q 4τ (B 4)
where φ τ pxq " HpΥ pxq ´Sd τ q. Here again, τ should be chosen sufficiently small for the dilated/eroded system to remain topologically equivalent to the original set. Elaborating on the same reasoning as for , we have chosen τ in such a way that maxp|S d |τ q remains of few dx. Instead of the level-set function Υ , whose computation from our data might be particularly expensive, we have used the following approximation of the level-set in the vicinity of ξ 0 :

Υ pxq " ξpxq ´ξ0 |∇ x ξ|pxq (B 5)
We have compared the results obtained for the S d -term by using a level-set field or the ansatz given by Eq. (B 5). The level-set field is computed from the scalar field ξ by use of the reinitialization procedure of [START_REF] Sussman | A level set approach for computing solutions to incompressible two-phase flow[END_REF]. We start with Eq. (B 5) as an initial condition and the algorithm was run over a sufficiently large number of iterations for the level-set to be a signed distance over the whole domain. The reinitialization procedure was solved using the two-phase flow solver archer [START_REF] Ménard | Coupling level set/VOF/ghost fluid methods: Validation and application to 3D simulation of the primary break-up of a liquid jet[END_REF].

Results are presented in Fig. 18(b) where the S d -term estimated either from the level-set or its approximation (Eq. (B 5)) are compared. Here we consider the T0 dataset where there is no forcing for either the dynamical or scalar field. corresponding to the different R λ collapse well for scales above r " 0.1R g . When R λ increases, the scale distributions widen when R λ increases and moves towards the small scales. It is further worth noting that the fractal scaling ends at a scale r « R g which means that R g plays the role of the outer cutoff of the fractal scaling.

We also carried out the same analysis for the transfer term ´x∇ r xpδuqpδφq 2 y E,R y Ω which is normalized using 2xφy E,R p1 ´xφy E,R qxky 1{2 E,R {R g . This quantity can be understood as the strain acting at a scale R g . In Fig. 21, we observe that this normalization leads to a good collapse of the different curves in the large scales r ą 0.1R g . For smaller scales, the different curves depart from each other and as R λ increases, the transfer rate term acts over a wider range of scales while its peak value moves towards smaller scales.

The different terms of the budget Eq. (2.15), normalized by the same quantities obey the same trend, i.e. the different curves collapse relatively well for scales larger than 0.1R g and move towards smaller scales when R λ is increased.

We also plot the budget for different Schmidt numbers in Fig. 23 using the large-scale similarity variables. Noticeable is the shift of the transfer term towards smaller r{R g when the Schmidt number increases from 0.1 to 1.0. In contrast, the diffusive term moves towards larger scales and increases in amplitude. A final observation of Fig. 23 is that, for ξ 0 " ξ rms , the forcing term due to the mean scalar gradient increases when Sc decreases from 1.0 to 0.1. For ξ 0 " 0, the budget for Sc " 0.1 is composed of the transfer term and diffusive term while the forcing term is almost negligible. Appendix E. Small-scale similarity

We here report the evolution of xpδφq 2 y E,R,Ω {ΣH ´1 rms with increasing R λ while the separation r is normalized by H ´1 rms . This scenario is tested in Fig. 24, where one observes a remarkable degree of similarity for both ξ 0 " 0 and ξ 0 " ξ rms . The range of scales over which this small scale similarity applies tends to increase with increasing R λ . It is further worth noting that the fractal scaling starts at a scale r « 20H ´1 rms (see the inset in Fig. 24). Hence, H ´1 rms appears to be proportional the inner cutoff of the fractal scaling.

We have observed that xpδφq 2 y E,R,Ω can be well represented by the following parametric expression

xpδφq 2 y E,R,Ω ΣH ´1 rms " r 2 " 1 `´r ηi ¯2 pζ´1q{2 " 1 `´r ηo ¯2 ζ{2 (E 1)
where ' " '{H ´1 rms . This parametrization is inspired by the one proposed by [START_REF] Batchelor | Pressure fluctuations in isotropic turbulence[END_REF] for representing the two-point statistics of the velocity field.

The main interest in deriving Eq. (E 1) is that it allows the fractal exponent ζ, the outer and inner cutoff η o and η i to be estimated using non-linear least-square curve fitting. By doing so, these parameters are gathered in an unambiguous way which does not imply any degree of arbitrariness, notably in the estimation of the best range of scaling. This is even more relevant for low to moderate Reynolds numbers. A similar approach was employed by Thiesset et al. (2016a); [START_REF] Krug | Fractal scaling of the turbulence interface in gravity currents[END_REF]. The appropriateness of this parametric expression is demonstrated in Fig. 24 at the largest and smallest Reynolds numbers. We obtain η i « 10η B which is in close agreement with the prediction of Thiesset et al. (2016b).

The small-scale similarity is now tested for the transfer term x∇ r xpδuqpδφq 2 y E,R y Ω which is normalized by K T ΣH ´1 rms while r is divided by H ´1 rms .

Results are presented in Fig. 25 confirming this small-scale similarity for both ξ 0 " 0 and ξ 0 " ξ rms . In particular, we note that it applies, at least up to the scale where the transfer term is maximum, and the range of scales over which the small scale similarity applies widens with increasing R λ . The maximum of the transfer term given in terms of K T ΣH ´1 rms is roughly constant and equals to 0.38.

The different terms of the budget Eq. (2.15) normalized using the smallsimilarity variables are presented in Fig. 26. Here again, we see that the transfer term complies well with a small-scale similarity when plotted using K T ΣH ´1 rms and H ´1 rms as similarity variables. However, the similarity holds only at the smallest scales for the S d d -term, where the term associated with S s d is negligible. When R λ increases, this term progressively tends to zero thereby leading to a closer degree of similarity for the S d d -term.

In conclusion, H ´1 rms and K T Σ{H rms play the same role as the Kolmogorov (or and perpendicular (dashed lines) to the mean scalar gradient. The insets represent the local scaling exponent. For the sake of clarity only case F1, F3 and F5 are presented. (a) ξ0 " 0, (b) ξ0 " ξrms. The small scale similarity variables have been used for normalization.

Batchelor) scales for normalizing the two-point statistics of the velocity (scalar) field. The statistical error on the calculation of the mean x'y R,E is then computed from ' rms { ? M where ' rms is the standard deviation computed from M snapshots.

We repeated this for xpδφ 2 qy R , its scaling exponent, and the different terms of its transport equation.

Case F5 is the most prone to numerical errors since it is the one for which we dispose of the smallest number of snapshots (M " 6, see Table 1). Hence, we infer errors for this particular case only.

Statistical errors are displayed in Fig. 28 as the blue error-bars while the mean value x'y E,R is plotted as the green curves. We observe that these errors are particularly small. The typical statistical errors on each of the plotted quantities are 1-2% for xpδφq 2 y R,E and its scaling exponent, and 3-4% for the terms of its transport equation. For the other forced cases (F0-F4), statistical errors are expected to be smaller. It is thus rather safe to conclude that statistical errors are marginal in our study.

  .4) Recall that S d is defined only at ξpxq " ξ 0 . Eq. (2.4) shows that the displacement speed can be decomposed into a scalar diffusion component S d d , a scalar reaction rate contribution S r d , and a scalar source term part S s d . Gran et al. (1996); Peters et al. (1998) further showed that the diffusion component of the displacement speed S d d can be further decomposed by projecting the diffusion term along the iso-scalar surface normal and tangential directions:

Figure 1 :

 1 Figure 1: Schematic representation of two points x `and x ´, the midpoint X " pX, Y, Zq and the separation vector r " prx, ry, rzq.

  In case of a passive scalar with no forcing, only S d d contributes to S d . When e.g. a scalar mean gradient G ξ in a given direction α is superimposed, then another contribution emerges from S s d " u α G ξ {|∇ x ξ|. ' The present framework also applies to premixed flames. Then, ξ can be associated to the fuel or oxidizer mass fraction or to the temperature field. In this situation, S d incorporates both S d d and S r d . Diffusion flames can also be analyzed using the present framework. In either premixed or diffusion flames, one generally define S d in the density weighted manner (Giannakopoulos et al. 2019). Note also that due to heat release, one should also account for the additional term due to ∇ x ¨u. ' In two-phase flows with phase change, ξ can for instance represent the liquid volume fraction and S r d relates to the evaporation rate. Here again, one should account for the term due to non-zero velocity divergence. The velocity jump

Figure 2 :

 2 Figure 2: Iso-surface colored by the displacement speed S d for increasing Schmidt numbers or Reynolds numbers. (a) F0, Sc " 0.1, (b) F0, 0.4, (c) F0, 1.0, (d) F2, (e) F4, (f) F5.

Figure 3 :Figure 4 :

 34 Figure 3: Surface averaged displacement speed and its components conditioned on the mean curvature H for ξ0 " 0. From dark to light, case F0 (Sc " 1.0) to F4. In (a) S d (full lines), S d d (dashed lines) and S c d (dotted lines) are portrayed while S n d (dashed lines) and S s d (full lines) are displayed in (b).

Figure 5 :Figure 6 :

 56 Figure 5: Surface averaged displacement speed and its components conditioned on the mean curvature H for ξ0 " ξrms. From dark to light, case F0 (Sc " 1.0) to F4. In (a) S d (full lines), S d d (dashed lines) and S c d (dotted lines) are portrayed while S n d (dashed lines) and S s d (full lines) are displayed in (b).

Figure 7 :

 7 Figure7: Evolution of xpδφq 2 y E,R,Ω {2xφy E,R p1 ´xφy E,R q with increasing R λ . The separation r is normalized by LΣ " 4xφy E,R p1 ´xφy E,R q{Σ. The local scaling exponent is also plotted in the inset. The dotted gray lines represent the asymptotic theoretical limits at large and small scales. (a) ξ0 " 0, (b) ξ0 " ξrms.

  xpδφq 2 y R,E,Ω the same role as the Taylor (Corrsin) microscale for normalizing the two-point statistics of the velocity (scalar) field.At intermediate scales, the influence of the Reynolds number is perceptible. It is observed that a pseudo 'inertial range' is forming whose extent increases with R λ . The local scaling exponent B logpxpδφq 2 y R,E,Ω q{B logprq is plotted in the inset of Fig.7and reveals that a power law with an exponent of about ζ « 0.36 ´0.38 applies over about one decade at R λ " 530. As shown by[START_REF] Vassilicos | Moving surfaces in turbulent flows[END_REF];[START_REF] Elsas | Geometry and scaling laws of excursion and iso-sets of enstrophy and dissipation in isotropic turbulence[END_REF], the distribution of xpδφq 2 y R,E at intermediate scales contains information about the fractal characteristics of iso-scalar surfaces. The scaling exponents relate to the fractal dimension of the iso-surface by D f " 3 ´ζ. The same relation is used in the community of fractal aggregates (see e.g.

Figure 11 :

 11 Figure 11: Budget of xpδφq 2 y E,R,Ω with increasing R λ . Full lines: transfer term, dashed lines: S d d -term, dash-dotted lines: S s d -term. All contributions are normalized by 2KT xφy E,R p1 ´xφy E,R q while the separation r is normalized by LΣ. (a) ξ0 " 0, (b) ξ0 " ξrms. The color legend is the same as in Fig. 7.

R

  λ , the strain rate K T and the curvature component of the stretch rate due to S d d , i.e. K d C " ´2xS d d Hy s balance, whilst the curvature component of the stretch rate due to S s d , i.e. K s C " ´2xS s d Hy s tends to zero. For ξ 0 " 0, the area weighted averaged displacement speed is zero. Fig. 11(a) confirms that all terms tend to zero at large scales in agreement with Eq. (2.22). For ξ 0 " ξ rms , only the transfer term approaches zero when r Ñ 8 while the S d d -and S s d -terms balance. The limit of the S d d -term at large scales is also in agreement with Eq. (2.22) which is displayed by the horizontal gray dotted lines in Fig. 11(b). The balance between the S d d -and S s d -terms at large scales suggests that the volume of the excursion set which naturally decreases due to diffusion effect is exactly compensated by the imposed mean gradient.

Figure 12 :

 12 Figure 12: Budget of xpδφq 2 y E,R,Ω with increasing Sc. Full lines: transfer term, dashed lines: S d d -term, dash-dotted lines: S s d -term, dotted lines: ´Sd d -term. All contributions are normalized by 2KT xφy E,R p1 ´xφy E,R q while the separation r is normalized by LΣ. (a) ξ0 " 0, (b) ξ0 " ξrms. The color legend is the same as in Fig. 8.

Figure 13 :

 13 Figure 13: Scaling of the different components of the stretch rate K " KT `KC with respect to R λ .(a) ξ0 " 0, (b) ξ0 " ξrms.

Figure 14 :

 14 Figure 14: (a) Ratio of the length-scales L ξ {Rg, λ ξ {LΣ and HrmsηB, with respect to P e λ ξ for both ξ0 " 0 and ξ0 " ξrms. The dotted lines represent a ratio of 0.5 and 0.22. (b) Surface density Σ versus RgHrms. The lines represent the expectations using a fractal dimension D f " 2.62.

Figure 15 :

 15 Figure 15: Budget of xpδφq 2 y E,RΩ with increasing ξ0 in decaying turbulence with imposed mean scalar gradient. Full lines: transfer term, dashed lines: S d d -term, dash-dotted lines: S s d -term, dotted lines: dt-term. All contributions are normalized by 2KT xφy E,R p1 ´xφy E,R q while the separation r is normalized by LΣ. (a) Sc " 1, (b) Sc " 0.2. Three values of ξ0 " 0, 0.5ξrms, ξrms are displayed from dark to light.

2KFigure 16 :

 16 Figure 16: Budget of xpδφq 2 y E,R,Ω with increasing ξ0 in decaying turbulence without imposed mean scalar gradient. Full lines: transfer term, dashed lines: S d d -term, dash-dotted lines: S s d -term, dotted lines: dt-term. All contributions are normalized by 2KT xφy E,R p1 ´xφy E,R q while the separation r is normalized by LΣ. (a) Sc " 1, (b) Sc " 0.2. Three values of ξ0 " 0, 0.5ξrms, ξrms are displayed from dark to light.

  It is now possible to cope with diffusive and/or reactive scalars, in presence or absence of source terms. All these processes are embedded in the interface displacement speed that may possess different contributions and different physical origins depending on the flow situation. It is an exact framework and has the potential of being applied to different flow variables in different flow situations. Hence, it is believed to offer promising perspective to probe the physics of interfaces in a broad sense.

Figure 18 :

 18 Figure 18: (a) Comparison of the spatially and angularly averaged S i d -terms for different values of τ , from max |S d |τ " 1 to max |S d |τ " 4dx. Case F0 Sc " 1.0, ξ0 " 0. (b) Spatially and angularly averaged S id -term as estimated from the level-set or its approximation using Eq. (B 5). Case T0 Sc " 1.0, ξ0 " ξrms. In (a) and (b), the S i d -terms are normalized by KT and the grey dotted lines represent the theoretical asymptotic limits at large and small scales.

  Fig. 18(a) displays the S d -term for different values of maxp|S d |τ q ranging from 1 to 4dx. We considered the forced turbulence, forced scalar case (F0-Sc " 1.0) for which S d is constituted of both a diffusion component S d d and a source component S s d . Results show that the computed S d -term is identical irrespective of the chosen value for τ . They also follow the expected asymptotic limit at small scales. Throughout the present study, we have chosen maxp|S d |τ q " 3dx.

Figure 20 :Figure 21 :

 2021 Figure 20: Evolution of xpδφq 2 y E,R,Ω {2xφy E,R p1 ´xφy E,R q with increasing R λ . The scale r is normalized by Rg. The local scaling exponent is also plotted in the inset. (a) ξ0 " 0, (b) ξ0 " ξrms.

Figure 22 :Figure 23 :

 2223 Figure 22: Budget of xpδφq 2 y E,R,Ω with increasing R λ . Full lines: transfer term, dashed lines: S d d -term, dash-dotted lines: S s d -term. All contributions are normalized by 2xφy E,R p1 ´xφy E,R qxky 1{2 E,R {Rg. (a) ξ0 " 0, (b) ξ0 " ξrms.

Figure 24 :

 24 Figure 24: Evolution of xpδφq 2 y E,R,Ω {ΣH ´1 rms with increasing R λ . The scale r is normalized by H ´1 rms . The local scaling exponent is also plotted in the inset. (a) ξ0 " 0, (b) ξ0 " ξrms. The gray dash-dotted lines denote the Batchelor type parametrization given by Eq. (E 1)

Figure 25 :Figure 26 :

 2526 Figure 25: Evolution of the transfer term ´x∇r ¨xpδuqpδφq 2 y E,R yΩ normalized by KT ΣH ´1 rms while r is divided by H ´1 rms . Dark to light represent F0 to F5. The local scaling exponent is also plotted in the inset. (a) ξ0 " 0, (b) ξ0 " ξrms.

Figure 27 :

 27 Figure 27: Anisotropy of xpδφq 2 y R,E between the direction parallel (full lines)and perpendicular (dashed lines) to the mean scalar gradient. The insets represent the local scaling exponent. For the sake of clarity only case F1, F3 and F5 are presented. (a) ξ0 " 0, (b) ξ0 " ξrms. The small scale similarity variables have been used for normalization.

Figure 28 :

 28 Figure 28: Statistical errors for case F5. Error bars for each quantity x'y R,E correspond to 'rms{ ? M . (a) xpδφq 2 y R,E and its local scaling exponent (inset), (b) the different terms of the budget.

  . Hence, the transfer term with respect to spatial position X is zero. Assuming a divergence free flow yields xpδφq 2 pσt∇ x ¨uuqy E,R " 0. Statistical homogeneity leads to 2xσtS d |∇φ|uy E,R " xS d |∇φ|y E,R and xφ ´Sd |∇ x φ| `yE " xφ `Sd |∇ x φ| ´yE . With " x|∇ x φ|y R is the surface density of the iso-scalar surface, i.e. the area of the iso-scalar surface divided by V box . x'y s denotes the area weighted average.

	orientations of the vector r (Thiesset et al. 2021):
		x'y Ω "	1 4π	ij	' sin θdθdϕ,	(2.16)
				Ω
	these simplifications, Eq. (2.12) then writes as:
	B t xpδφq 2 y R loooomoooon	" ´∇r ' xpδuqpδφq 2 y R looooooooooomooooooooooon	loooooooooooooooooomoooooooooooooooooon `2xS d y s Σ ´4xφ `Sd |∇ x φ| ´|y R	. (2.15)
	Unsteady	Transfer´r			S d ´Term
	where Σ At this stage, Eq. (2.15) depends on time t and the separation vector r (a 4D
	space). The problem can further be reduced supposing isotropy, i.e. statistical
	invariance by rotation of the separation vector r, i.e. x'y R prq " x'y R p|r|q. In
	case of anisotropic flows, one can apply an angular average, noted x'y Ω , over all

  " K T `KC measures the relative time increase of surface density. K T denotes the tangential strain rate that may include compressibility effect, and K C " ´2xS d Hy s is the curvature component of the stretch rate.

		al. 2019):
	B t Σ `∇x ' xuy s Σ " pK T `KC qΣ,	(2.23)
	where the stretch rate K

Thiesset et al. (2021) considered the case where S d " 0 (a material surface

Table 1 :

 1 Physical parameters and typical one-point statistics of the DNS database used in the present work

	0.00
	0.00
	1.00

Table 2 :

 2 Expression of the different components of the displacement speed in forced turbulence with a mean scalar gradient G ξ .

	Total S d " S s d	`Sd d
	= Source	S s d "	uy G ξ |∇xξ|
	+ Diffusion	S d d " S c d	`Sn d " ∇x'D∇xξ |∇xξ|
			= Tangential diffusion	S c d " 2DH
			+ Normal diffusion	S n d " n'∇xpn'D∇ξq |∇xξ|

These surface conditioned statistics are noted xS i d |Hy s (i " tc, n, su denotes the contribution to the displacement speed). Such conditional statistics are relevant notably for providing insights into the curvature component of the stretch rate K C , which reads as the mean product of the displacement speed S d by the mean curvature H.

We first analyze the effect of R λ on S d and its components. Surface conditioned

  M. Gauding, F. Thiesset, E. Varea, L. Danaila which may also vary with R λ and Sc.We first characterize the influence of the Péclet number on K T , K d C and K s C .K T is inferred from Eq. (2.24), i.e. from the slope of the transfer term in the limit of small scales. Similarly, Eq. (2.25) reveals that K d C and K s C can be obtained from the slope of the S d d and S s d terms, respectively, when r approaches zero.

	R "	xξ 2 y x ξ y	x y xky	(4.8)

  allows one to easily grasp that for small values of the separation |r|, xpδφq 2 y R is proportional to the area of the ξ 0 iso-surface, while for large scales, xpδφq 2 y R depends on the volume of the excursion set. For intermediate scales, xpδφq 2 y R becomes a morphological descriptor where the scale r plays the role of a morphological parameter. Appendix B. Numerical method for computing xφ ´Sd |∇ x φ| `yR

Acknowledgements

We are thankful to A. Poux, research engineer at CORIA laboratory for his contribution to the pyarcher project. We are also grateful to Marc Massot for

Declaration of Interests

The authors report no conflict of interest. No noticeable differences are observed, and both curves compare well with the asymptotic theoretical limits at either large and small scales. This proves that our method for computing xφ ´Sd |∇ x φ| `yR either from the real or approximated level-set is accurate.

Appendix C. Validation of the partial angular average

We checked that results obtained by operating a partial angular average over only the three directions r x , r y , r z were similar to the ones issued from the full 3D angular average. We considered here the case F0, Sc " 1.0, ξ 0 " 0 (Fig. 19(a)) and the T0 case (Fig. 19(b)).

A careful analysis of Fig. 19 reveals that the budget is accurately closed even if one employs the partial angular average. Some very slight differences are perceptible which are due to a small anisotropy. Note that, in absence of mean scalar gradient, the anisotropy should be interpreted as a statistical effect. Indeed, the number of structures formed by the scalar excursion set ξpxq ą ξ 0 is typically around 50-100 which might not be enough for the two-point statistics of isovolumes to be fully converged. Increasing the number of independent simulations will very likely improve the statistical convergence and the degree of isotropy (a similar situation was encountered in Thiesset et al. (2020)). This effect is however considered to be marginal in the present analysis. The small values of the budget residuals shown in Fig. 19 is a further evidence of the appropriateness of our postprocessing procedure for computing the triple correlation xφ ´Sd |∇ x φ| `yR . The transfer and unsteady term also compare favorably well with their asymptotic theoretical expressions at small and large scales. To conclude Fig. 19 provides the validation of altogether the theory (Eq. (2.15)), the DNS data and the postprocessing procedures.

Appendix D. Large-scale similarity

Here, we present the appropriate normalization of two-point statistics in the large-scale limit.

In Fig. 20, the second-order structure function xpδφq 2 y E,R,Ω is normalized by its asymptotic large-scale value 2xφy E,R p1 ´xφy E,R q while the separation r is normalized using the radius of gyration R g . We observe that the different curves Appendix F. Anisotropy effects due to the mean scalar gradient

The present numerical configuration leads to statistical anisotropy due to the presence of a mean scalar gradient G ξ . Since the latter is active in the y direction, two-point statistics are invariant by rotation around the r y axis. Here, we infer anisotropy from the variations of xpδφq 2 y R,E along the different orientations of the separation vector r.

In Fig. 27, we plot xpδφq 2 y R,E along the direction parallel and perpendicular to the mean scalar gradient. The parallel direction corresponds to r y while the contributions along r x and r z were averaged because of axisymmetry. We observe that the curves corresponding to the parallel and perpendicular directions collapse at small up to intermediate scales, meaning that the second-order statistics of φ are locally isotropic. The previously discussed scaling with respect to H rms (and Σ) holds. The scaling exponent remains the same irrespective of the orientation of the separation vector r. This means that the estimation of the fractal dimension of the present iso-scalar surfaces is quite robust and is not affected by the mean scalar gradient. Differences between parallel and perpendicular directions (viz. anisotropy) are perceptible only at the end of the scaling range and at large scales. It is seen that the scaling range is systematically wider in the direction perpendicular to the mean scalar gradient.

Appendix G. Estimation of statistical errors

We here address the question of statistical errors. We attribute such errors to a lack of statistical convergence. For the forced turbulence cases explored here, any spatially averaged quantity x'y R is also averaged over M independent snapshots which constitute our ensemble E.