J.-N Roux

Pre-peak deformation of model granular materials: a DEM study

Initial states of granular assemblies, which determine the response to quasi-statically applied deviator stresses in mechanical tests, are usually classified according to their density and structural anisotropy. Discrete simulations (DEM) of model granular materials (spherical bead packs) reveal however that dense assemblies with isotropic fabric can be prepared with very different coordination numbers, ranging from 4 to 6 under low confining stress. When subsequently subjected to triaxial compression tests equally dense systems differing by their coordination number have the same peak strength, but strongly differ in their pre-peak stressstrain response. Well-coordinated systems have contact networks that tend to resist throughout large deviator intervals, and exhibit a stiff response, with a global stiffness scaling with the contact stiffness (type I strains). Poorly coordinated ones essentially deform in a series of microscopic instabilities in which contact networks get continuously broken and repaired (type II strains). The distinction between strains of types I and II (not to be confused with elastic and plastic strains) has important consequences in terms of sensitivities of macroscopic response to micromechanical features and propensity to develop localization instabilities. The coordination number, which determines to a large extent the type (I or II) of the response, can be related to initial elastic moduli. These observations call for changes in the practice of DEM approaches and refinements in the classification of internal states and stress-strain behaviour of granular materials. Some experimental confrontations are possible.

INTRODUCTION

Discrete element methods [START_REF] Radjaï | Discrete-element modeling of granular materials[END_REF][START_REF] O'sullivan | Particulate Discrete Element Modeling: a Geomechanics perspective[END_REF] have become widespread numerical simulation tools in the investigation of the microscopic, grain-level origins of the macroscopic mechanical properties of granular materials. The essential features of the drained behaviour of sands, as classically probed in triaxial compression tests, are thus reproduced by simulations of model materials such as spherical bead assemblies [START_REF] Thornton | Numerical simulations of deviatoric shear deformation of granular media[END_REF][START_REF] Suiker | Frictional collapse of granular assemblies[END_REF]. In particular, the deviator stress maximum and the dilatant behaviour of dense samples, the monotonic deviator increase with axial strain and the contractant behaviour of loose samples, as well as the critical state approached for large strain whatever the initial density, are retrieved [START_REF] Thornton | Numerical simulations of deviatoric shear deformation of granular media[END_REF][START_REF] Radjaï | Contact dynamics study of 2d granular media : critical states and relevant internal variables[END_REF]). The influence of intergranular (contact) friction was assessed [START_REF] Lemaître | What do dry granular flows tell us about dense non-brownian suspension rheology?[END_REF], and, beyond the simplest model materials, that of rolling resistance at contacts [START_REF] Estrada | Shear strength and force transmission in granular media with rolling resistance[END_REF], wide polydispersities [START_REF] Voivret | Multiscale force networks in highly polydisperse granular media[END_REF], non-spherical grain shapes, either smooth [START_REF] Antony | Evolution of force distribution in three-dimensional granular media[END_REF] or angular [START_REF] Azéma | Force transmission in a packing of pentagonal particles[END_REF], have also been investigated. Most studies, however, focus on maximum de-viator stress or on critical state properties, and classify initial states by their sole density (with possible, yet limited, effects of inherent anisotropy).

The purpose of the present communication is to show that the pre-peak behaviour of dense, isotropic systems, is not only determined by the initial density, or influenced by inherent anisotropy [START_REF] Lam | Effects of initial anisotropic fabric and σ 2 on strength and deformation characteristics of soils[END_REF], but is sensitive to the coordination number, which, according to the sample assembling procedure, may vary independently of the density (Agnolin & Roux 2007a).

MODEL MATERIAL, INITIAL STATE AND ISOTROPIC COMPRESSION

The numerical model, with the material properties and the simulation procedure are described in detail by Agnolin & Roux (2007a). We consider assemblies of elastic and frictional spherical beads of diameter a, with elastic-frictional contact forces. The friction coefficient is µ = 0.3. Contact elasticity abides by Hertz's law for the normal force [START_REF] Johnson | Contact Mechanics[END_REF], relating normal force F N , to deflection h, for a material of Young modulus E and Poisson ratio ν, as

F N = Ẽ√ a 3 h 3/2 , with Ẽ = E 1 -ν 2
or contact stiffness K N to h, or to F N , as

K N = dF N dh = Ẽ√ a 2 h 1/2 = 1 2 (3a) 1/3 Ẽ2/3 F 1/3 N . (1)
A simplified model is implemented for Mindlin-Deresiewicz laws (Agnolin & Roux 2007a) combining tangential elasticity and friction in an energetically consistent way [Elata and Berryman 1996]. Dense numerical samples are assembled under isotropic pressure P in cuboidal cells periodic in all three directions, with appropriate dynamics such that the cell dimensions L x , L y , L z may adapt to prescribed stress levels σ xx , σ yy , σ zz . Starting from loose configurations with no contact, all three principal stress values are set to P , and the cell is compressed until equilibrium is reached with intergranular forces balancing the prescribed stresses. Particles have the elastic properties of glass (Young modulus E = 70 GPa, Poisson ratio ν = 0.3). Under pressure P , a dimensionless stiffness number is defined as

κ = Ẽ P 2/3 , (2)
so that the typical contact deflection h in contacts, due to the Hertz law (1), satisfies

h a ∝ κ -1 . (3)
The results are based on 5 statistically similar samples of 4000 particles for each preparation procedure. As explained by Agnolin & Roux (2007a), dense packings are obtained on suppressing the effects of friction in the assembling stage. Setting µ = 0 results in the maximum solid fraction Φ = 0.64 ("random close packing") of disordered sphere assemblies, and coordination number z * = 6 (excluding a small proportion x 0 ≃ 2% of "rattlers", grains carrying no force). But if the system is strongly vibrated at high Φ, just below the final equilibrium value, and then finally compressed with µ = 0.3, a much lower coordination number is obtained, z * = 4.5 for nearly equal Φ, while x 0 exceeds 10%. Samples numerically prepared according to these two methods, subsequently referred to as A and C, under low pressure P = 10 kPa or κ = 39000) are then isotropically compressed (with friction) to various confinement levels (up to 1 MPa or κ = 1810), before subjecting them to triaxial compression tests in which "axial strain" ǫ a = ǫ 1 is increased while lateral stresses σ 2 = σ 3 are kept constant. Tab. 1 gives values of solid fraction Φ, coordination number z * (for non-rattler grains) and rattler fraction x 0 in states A and C under isotropic pressure. For reference, we also included similar data on initial state B, prepared by direct compression like A, but with a small friction coefficient µ = 0.02, and on initial state D, assembled with the final friction coefficient value µ = 0.3. It should be remarked that B is less dense than C, but better coordinated. Under growing pressure, density and coordination numbers increase, but moderately so in the pressure range, up to 1MPa, for which simulation results are presented below. From 10 kPa (κ = 39000) to 1 MPa (κ = 1810), z * raises to about 6.2 in A samples, to about 5.1 in C ones (and x 0 decreases to about 6%), while Φ increases by less than 0.02 in all cases (Agnolin & Roux 2007b).

STRESS-STRAIN CURVES AND EFFECT OF CONFINING STRESS

Figure 1: Normalized deviator q/σ 3 and volumetric strain ǫ v vs. axial strain ǫ a for 5 different C-type samples prepared in the same conditions.

We represent deviator stress, normalized by constant lateral stress, q/σ 3 = (σ 1 -σ 3 )/σ 3 , and volumetric strain ǫ v versus axial strain ǫ a . Fig. 1 checks for reproducibility of the results on comparing different samples prepared in the same conditions. (Other checks, not shown here in figures for lack of space, show that the implemented DEM procedure indeed probes a quasi-static, rate-independent reponse, and that creep effects vanish for small enough strain rateds ǫa ). Fig. 2 then compares deviator stress, plotted up to its maximum, and volumetric strain curves obtained from initial states A and C for lateral stress σ 3 = 100 kPa.

As expected, given their nearly identical initial density, both systems share the same maximum deviator value, about q/σ 3 = 1.4, within statistical uncertainties (see Fig. 1), corresponding to internal friction angle ϕ ≃ 24 • . Volumetric curves also exhibit the same slope (equal dilatancy), with a more pronounced initial contractant regime in system C. But the axial strain corresponding to a given level of internal friction mobilisation is strikingly different, with, e.g., deviator q = 1.1σ 3 being reached for strains about one order of magnitude larger in the C case (∼ 10 -2 ) than in systems with initial state A (∼ 10 -3 ). This behaviour is better understood on comparing triaxial tests carried out on one single sample prepared under varying confining pressure. Fig. 3 shows that stress-strain curves on one sample of type A, for σ 3 ranging from 10 kPa to 1 MPa, strongly depend on σ 3 . The larger the confining stress, the softer the macroscopic response expressed in terms of ratio q/σ 3 . To explain such observations, one may note that material deformations, at the microscopic level, i.e., contact deflections, relative to grain diameter a, scale as κ -1 a. If macroscopic strains are averages of microscopic ones in the contact regions, then they should scale as κ -1 a too. This is checked in Fig. 4, which shows that, at least within quite a large initial deviator range, up to about q = 0.9σ 3 , stress-strain curves pertaining to different confinement pressures coincide if plotted with stress ratio q/σ 3 and rescaled strains κǫ. Fig. 5 and Fig. 6 respectively display stress ratio q/σ 3 and volumetric strain ǫ v versus axial strain in one C-type sample under the same confining stress values as in Figs. 34(initial pressure P = σ 3 ).

In sharp contrast with the results obtained with Atype systems,stress-strain curves are now reproduced for the different confinement levels without rescaling strains: macroscopic strains do not scale with contact deflections in that case.

Behaviours of A and C-type samples in the prepeak deviator range are characteristic of two different physical origins of macroscopic strains. For a material specimen prepared in state A, strains result from local deformations in the intergranular contact regions -we refer to such strains as type I strains. With C samples, macroscopic strains are of geometric origin, Figure 6: Same as Fig. 3 for volumetric strain in one C sample, the inset magnifying the very small strain region.

stemming from contact network failure and rearrangement (type II strains). Such strains are also present in A samples in a certain deviator range closer to the maximum (say, above q = σ 3 ). One can also observe type I strain response over a significant deviator interval in both systems A and C upon unloading [START_REF] Roux | Quasistatic rheology and the origins of strain[END_REF], in which both kinds of strains were characterized in some detail with various methods, on two-dimensional numerical simulations of disk assemblies).

An interval of loading parameters in which type I response is observed corresponds to the stability range of a contact network, whence the much larger loading range corresponding to such behaviour in well-coordinated systems such as A. (Actually, a closer scrutiny of results on system C would reveal a small range of type I response, up to about q/σ 3 = 0.1). On unloading, the evolution of tangential contact forces tends to be reversed, thus suppressing full friction mobilisation in many contacts and stabilising the network.

Within a type I strain interval, the system is microscopically stable, each one of its successive configurations can be made arbitrarily close to an equilibrium configuration upon reducing the strain rate. As a consequence, the macroscopic behaviour should be devoid of instability. A strain response of type II, on the other hand, involves microscopic instabilities, as broken contact networks are not repaired until new contacts form, which implies some finite strain interval during which a momentary surge of kinetic energy is observed. This may result in a macroscopic instability which is not, of course, bound to happen. While microscopic stability entails macroscopic stability, microscopic instabilities, which occur with a frequency that increases with the system size [START_REF] Roux | Quasistatic rheology and the origins of strain[END_REF] are compatible with well-behaved constitutive laws expressing smooth stress-strain dependencies on the macroscopic scale. The absence of microscopic instability characterising type I strains, or, equivalently, the enduring capability of a given contact network to support varying deviator levels, has been related by [START_REF] Mcnamara | Quasirigidity: some uniqueness issues[END_REF] to a positive second-order work criterion, expressed in terms of the elastoplastic stifffness matrix of the contact network.

INITIAL ELASTIC PROPERTIES

In experimental practice, the coordination number is usually not accessible to direct measurements. Yet we now show that elastic moduli, which can be measured [START_REF] Shibuya | Elastic deformation properties of geomaterials[END_REF][START_REF] Di Benedetto | Sand behaviour in very small to medium strain domains[END_REF][START_REF] Kuwano | On the applicability of crossanisotropic elasticity to granular materials at very small strains[END_REF], provide indirect information on contact nework coordination. In Figs. 5 and 6 the insets (subplots), with blown-up strain scales, focus on the initial, very small strain response. The straight lines then correspond to the elastic behaviour of the initial contact network, with slope E * , the effective Young modulus of the granular system, in Fig. 5, and-(1 -2ν * ), with ν * its effective Poisson ratio, in Fig. 6. These moduli are obtained on computing the response of the initial contact network, assuming each contact responds elastically (to see how these moduli are calculated, on building appropriate stiffness matrices, see Agnolin & Roux (2007c)). We thus observe in simulations that the moduli give the tangent slopes to the stress-strain curves. In a very small range, for strains typically of order 10 -5 (just like in experiments [START_REF] Di Benedetto | Sand behaviour in very small to medium strain domains[END_REF])), the incremental behaviour about the initial equilibrium state can be regarded as elastic in good approximation (Agnolin & Roux 2007c). The elastic reange should not be confused with the type I strain regime, as the lat- ter may obviously extend to much larger stress and strain intervals (Fig. 4). Figs. 7 and7 respectively plot bulk (B) and shear (G) moduli of all samples of Tab.1, once isotropically compressed, versus pressure P . Those graphs make it very clear that moduli are primarily sensitive to coordination number, as C and D results are close to each other, despite different densities. Results for state B, not shown in Figs. 7 and 8 are actually close to A ones. Note, also, that those moduli (especially shear modulus G) of well coordinated and poorly coordinated packings exhibit different pressure dependences -different apparent pressure exponents, see the slopes in the logarithmic plots of Figs. 7 and8.

More detailed analyses of elastic properties of isotropic spherical particle assemblies are given by Agnolin & Roux (2007c) (such as a discussion of "anomalous" shear moduli in poorly coordinated isotropic assemblies) along with some experimental comparisons for glass bead systems.

CONCLUSIONS

Numerical simulations clearly reveal that the initial coordination number in a dense isotropically equilibrated granular material is not determined by the density, and strongly influences the stress-strain response in the pre-peak range. This is related to the ability of well-coordinated contact networks to support large deviator stress increases without failure, whence macroscopic strains of type I. While systems prepared with large coordination numbers exhibit type I strains, i.e., deform because the contacts are deformable, poorly coordinated ones macroscopically deform because contacts networks get continuously broken and repaired. Initial coordination numbers appear to be directly related to elastic moduli.

These observations have important practical and fundamental implications.

On the practical side, these results call for studies of initial state coordination numbers, depending on assembling procedures, which should imply more systematic comparisons between experiments and simulations, preferentially for numerically tractable model materials. It seems that one should eventually refine the classification of initial states in terms of density and inherent fabric anisotropy only, as the coordination number is another source of variability. Results of triaxial tests on sands under varying confining stresses often resemble Figs. 5 and 6 more than Fig. 3, which certainly raises doubts on the relevance of the numerical procedure of compression without friction, most often resorted to in practice because of its simplicity, to produce dense samples.

Figure 3 :

 3 Figure 3: Stress-strain curves for initial pressure P = σ 3 , in one A sample, taking 5 different values growing from 10 kPa to 1 MPa by factors of √ 10.

Figure 4 :

 4 Figure 4: Same as Fig. 3, with rescaled strain coordinates (P 0 = 100 kPa).

Figure 5 :

 5 Figure 5: Same as Fig. 3 for normalized deviator stress in one C sample. Inset shows detail of curve near origin (note the blow-up scale).

Figure 7 :

 7 Figure 7: Bulk modulus B vs. pressure P , on logarithmic scales, in isotropic compression of samples A, C, D of Tab. 1.

Figure 8 :

 8 Figure 8: Shear modulus G vs. pressure P in isotropic compression of samples A, C, D of Tab. 1.

Table 1 :

 1 Isotropic states, with µ = 0.3 (κ ≃ 39000 for A and C, κ ≃ 181000 for B and D), obtained with different assembling procedures.

	Procedure	Φ	z *	x 0 (%)
	A (µ 0 = 0.) B (µ 0 = 0.02) C (vibrated) D (µ 0 = µ = 0.3) 0.5923 ± 0.0006 4.546 ± 0.009 0.6370 ± 0.0002 6.074 ± 0.0015 0.6271 ± 0.0002 5.80 ± 0.007 0.635 ± 0.002 4.56 ± 0.03	1.3 ± 0.2 1.95 ± 0.02 13.3 ± 0.5 11.1 ± 0.4

Figure 2: Deviator stress and volumetric strain curves averaged over available A-type and C-type samples.


Auguste Gires
email: augustegires@enpc.fr

Jean-Baptiste Abbes

Igor Da

Silva Rocha Paz

Ioulia Tchiguirinskaia

Daniel Schertzer

Multifractal characterisation of a simulated surface flow: a case study with Multi-Hydro in Jouy-en-Josas, France

In this paper we suggest to innovatively use scaling laws and more specifically Universal Multifractals (UM) to analyse simulated surface runoff and compare the retrieved scaling features with the rainfall ones. The methodology is tested on a 3 km 2 semi-urbanised with a steep slope study area located in the Paris area along the Bièvre River. First Multi-Hydro, a fully distributed model is validated on this catchment for four rainfall events measured with the help of a C-band radar. The uncertainty associated with small scale unmeasured rainfall, i.e. occurring below the 1km x 1km x 5min observation scale, is quantified with the help of stochastic downscaled rainfall fields. It is rather significant for simulated flow and more limited on overland water depth for these rainfall events. Overland depth is found to exhibit a scaling behaviour over small scales (10 m -80 m) which can be related to fractal features of the sewer network. No direct and obvious dependency between the overland depth multifractal features (quality of the scaling and UM parameters) and the rainfall ones was found.

1) Introduction

The combined effects of a growing urbanisation -approximately 80% of Europe's population will live in cities by 2020 (EEA, 2014) -and potential increase of extreme events as a consequence of climate change (IPCC, 2013) expose more and more people to surface pluvial flooding. [START_REF] Pitt | The Pitt Review : Learning lessons from the 2007 floods[END_REF] carried out a review on flood events in the United Kingdom and showed that two thirds of the flood related damages were caused by surface water flooding.

Urban flooding has become a growing concern in Europe, hence a significant number of European research projects address this issue, along with national counterparts. The purpose of these projects is to increase the resilience of urban areas through improvement of both real time management of extreme events and long term planning. We can cite FP7 SMARTesT (http://floodresilience.eu/), CORFU (http://www.corfu-fp7.eu/), Climate KIC Blue Green Dream (www.bgd.org.uk) or the INTERREG IV RainGain project (http://www.raingain.eu) among others.

There is a need to improve the understanding of urban surface flow. Indeed, there is a growing interest for 2D models in urban environment for both operational and research applications [START_REF] Bolle | Hydraulic modelling of the two-directional interaction between sewer and river systems, Urban Drainage Modelling and Water Sensitive Urban Design[END_REF][START_REF] Carr | Linking of 2D and Pipe hydraulic models at fine spatial scales, Urban Drainage Modelling and Water Sensitive Urban Design[END_REF][START_REF] Chen | The urban inundation model with bidirectional flow interaction between 2D overland surface and 1D sewer networks[END_REF]Deltares, 2013;DHI, 2011;[START_REF] Giangola-Murzyn | Multi-Hydro : an open access coupling model for urban hydrology[END_REF][START_REF] Innovyze | InfoWorks CS v13.0.6. Innovyze[END_REF]Innovyze, , 2103;;[START_REF] Phillips | 1D and 2D Modelling of Urban Drainage Systems using XP-SWMM and TUFLOW[END_REF]XP Solutions, 2012). Such models aim at actually modelling processes in a physically based manner, while the most commonly used semi-distributed models take them into account through tailored lumped models. In case of overflow they simply consider a volume output from the sewer system and deduce a local water depth, but the dynamical behaviour of the water added on the ground is not addressed. Basically, urban surface flow is not commonly perceived as a geophysical process and is therefore not addressed with geophysical tools capable of grasping its intrinsic complexity visible across all scales. Indeed, it results from the non-linear interactions between the highly spatially and temporally variable rainfall field, the topography and the strongly inhomogeneous land use cover.

In this paper we suggest to use multifractal tools, which are commonly used in geophysics to characterise and simulate fields extremely variable over a wide range of scales; such as wind turbulence, rainfall, river flow or topography (see [START_REF] Schertzer | Multifractals, generalized scale invariance and complexity in geophysics[END_REF] for review). Such tools have seldom been used in an urban context. [START_REF] Gires | Multifractal analysis of a semi-distributed urban hydrological model[END_REF]Gires et al. ( , 2014b) ) used them to downscale rainfall to quantify the uncertainty associated with small scale rainfall variability, or to characterise the variability across scales of simulated flow in conduits (sewer). To the knowledge of the authors it has never been used to study either surface runoff flow (urban drainage) or surface flow in general including stream rivers. Investigating the potential multifractal features of surface flow and notably whether it inherits rainfall features is the main purpose of this paper and constitute its main novelty. In addition, this case study will also be used to quantify the uncertainty associated with small scale rainfall variability, not only on the simulated flow which has already been done on other catchments, but also on the surface flow.

Given the lack of measurements of distributed data of surface runoff, outputs of a numerical model are analysed. The model used is Multi-Hydro (El Tabach et al., 2009 for an initial version and Giangola-Murzyn, 2014 for a recent one) developed at the Ecole des Ponts ParisTech. It is implemented on a 3.017 km² peri-urban catchment in Jouy-en-Josas (South-East of Paris), which exhibits steep slopes and both forest and urbanised areas. Achieving such an analysis is relevant only if a distributed rainfall field is used as model input. Météo-France radar mosaics with a resolution of 1 km in space and 5 min time (Tabary, 2007;Tabary et al., 2007) for four events that occurred between 2009 and 2011 are used. When needed, the rainfall field is downscaled both in space and time from the raw radar data, in order to simulate the improvement that could be made with higher radar resolution.

The model and the study area data for its implementation are presented in details in section 2. The multifractal framework and analysis methods are presented in section 3.

Results are discussed in section 4 and 5. More precisely, the validation of the model and quantification of the uncertainty associated with small scale unmeasured rainfall variability on both simulated sewer flow and maximum water depth is carried out in section 4. Multifractal characterization of overland water depth is addressed in section 5. Main conclusions are highlighted in section 6.

2) Model and catchment

2.1) The Multi-Hydro model

Multi-Hydro is a multi-module model whose goal is to model and predict the impacts of rainfall events in urban and peri-urban areas. In this paper, there is an emphasis on heavy rainfall events. Following the approach of various recent developments of hydrological models [START_REF] Djordjevic | An approach to simulation of dual drainage[END_REF][START_REF] Fletcher | Understanding, management and modelling of urban hydrology and its consequences for receiving waters: A state of the art[END_REF][START_REF] Hsu | Inundation simulation for urban drainage basin with storm sewer system[END_REF][START_REF] Jankowfsky | Understanding and modelling of hydrological processes in small periurban catchments using an object-oriented and modular distributed approach. Application to the Chaudanne and Mercier sub-catchments (Yzeron catchment[END_REF][START_REF] Rodriguez | A distributed hydrological model for urbanized areas -model development and application to case studies[END_REF]; it makes different modules interact, each of them echoing a portion of the water cycle in urban areas (surface runoff, infiltration, ground water flow, sewer flow).

Each of the modules integrated in Multi-Hydro relies on open-source software packages that have already been widely used and validated by the scientific community. The surface module is based on TREX (Two dimensional Runoff, Erosion and eXport model, [START_REF] Velleux | TREX Watershed Modelling Framework User's Manual: Model Theory and Description[END_REF] which solves fluid mechanics equations for surface flow (diffusive wave approximation of 2D Saint-Venant, see p. 6-7 of the TREX user manual) and infiltration (simplification of Green and Ampt equation). The sewer or drainage module, which is based on SWMM developed by the US Environmental Agency (Storm Water Management Model, [START_REF] Rossman | Storm Water Management Model, User's Manual[END_REF], is a 1D-model dealing with sewer flows through numerical solutions of Saint-Venant 1D equations in pipes. The interactions between the surface and sewer flow is handled through the gully pixels. These interactions (input or output of water) between the surface and sewer flow are carried out every 3 min. When there is no overflow, gully pixels are considered to have an infinite infiltration rate, and the water passing through them is directly inputted into the corresponding node of the sewer model. This way of modelling implies that a large transport capacity is assumed for gully, especially with 10 m pixel size as in this paper (see below). Future developments of Multi-Hydro will enable to improve the model with regards to this coarse assumption. They could notably rely on the experimental and computational studies of gully inflow capacity, including 3D CFD studies, which analyse phases in the flow, inlet capacity, reverse flow when the piezometric level in the sewer is beyond the ground level [START_REF] Despotovic | Inefficiency of storm water inlets as a source of urban floods[END_REF]Djorjevic et al., 2005). In case of sewer overflow through a node, the corresponding gully pixel is converted into a road pixel and the water exiting the node is inputted on this pixel (considered as a source in TREX). There is also a module handling ground water flow which was not included in this study to limit computation time.

In order to run Multi-Hydro, data needs to be shaped in a standard format. Commonly available Geographical Information System (GIS) data, such as land use and topography provided in France by IGN (the French agency producing geographical information) are inputted to MH-AssimTool [START_REF] Richard | Gis data Assimilation interface for distributed hydrological models[END_REF]. This software formats the inputs with the desired resolution and makes Multi-Hydro a transportable model, rather easy to implement on a new catchment. Once a resolution is chosen, one has to affect an elevation and a land use class to each pixel. The elevation is obtained by an interpolation of the raw available data.

With regards to the land use, a priority order has been determined to assign a unique land use class for a given pixel according to the hydrological importance of the given class instead of the surface represented by this class: if a gully is located on a pixel, the entire pixel will be considered as a gully. This process is repeated in the following order for this case study: roads, houses, forest, grass, and water surface. See [START_REF] Ichiba | Scale effect challenges in urban hydrology highlighted with a distributed hydrological model[END_REF] for a comparison with other possible strategies.

In this paper, the model was implemented with pixels of size 10 m x 10 m. Given the obtained results discussed below it was not found necessary to run it at higher resolution which makes computation time too long. For an in-depth analysis of the relation between the selected pixel size and simulated flow, which is not the purpose of this paper, refer to Ichiba (2016). Multi-hydro was not calibrated, i.e. standard values for the parameters describing a land use class are used (hydraulic conductivity, capillary suction, moisture deficit, Manning's coefficient, depth of interception). Raw or downscaled radar data are used as input of the model.

2.2) Presentation of the study area

The catchment studied in this article is located in Jouy-en-Josas (Yvelines County, South-west of Paris). It occupies a 3.017 km² area, mainly on the left bank of the Bièvre River. A small portion of the right bank near the river bed is also included. The remaining portion of the right bank is drained to a small river that flows into the Bièvre River downstream the outlet of the studied catchment. The Bièvre River is a tributary of the Seine River which it meets in Paris. It flows through increasingly urbanised areas along its 33 km path. This has led to strongly modify its natural bed, both in underground pipes which are integrated in the storm water sewer system, or in a highly artificial open air bed. An effort is currently undertaken to restore its "natural" aspect.

A striking feature of this catchment is that, unlike the previous ones studied with Multi-Hydro (Giangola-Murzyn et al., 2014;Gires et al., 2014a), it exhibits steep slopes.

There is a difference of approximately 100 m between the plateau in the north of the catchment, and the outlet of the catchment (Fig. 1). The downhill portion strengthens overland runoff, and the combination of pluvial and fluvial processes on the river bank has led to severe flooding in 1973 and 1982. Some details are available on the SIAVB (Syndicat Intercommunal d'Assainissement de la Vallée de la Bièvre, the local authority in charge of urban drainage of the area) website http://www.siavb.fr/gestion_des_crues.aspx. Urbanisation and imperviousness are concentrated along the river bank, and on a housing estate along one major North-South road. The remaining of this semi-urban catchment is mainly made of forests. The sewer system is a separate one, and the storm water is routed into the Bièvre River.

Following the severe flooding, the SIAVB has created 15 storage basins (integrated in the landscape) along the Bièvre River to mitigate flooding risks. One, the Bassin des Bas Près, is located just upstream the Jouy-en-Josas catchment. The outlet of this basin is equipped with flow and height gauges operated in real time. There is a second measuring point of water depth, few meters upstream the outlet of the catchment, at the "Pont de Pierre" (Fig. 1). This gauge has been installed to monitor the river level and to protect a music school by triggering a warning system in case of elevated height. Given the position of the two measuring points, Multi-hydro will only be validated on the area drained by the sewer network represented in green in Fig. 1. The forest corresponds approximately to 60% of the catchment (~ 2 km 2 ).

Although it is only possible to validate the implementation of the model on a portion of the catchment, the whole area is modelled to ensure the accuracy of flow over the areas actually used for validation. The river is part of the storm water sewer system in Jouy-en-Josas and is modelled as a pipe in Multi-Hydro drainage module. Indeed, through the city, the river bed is highly artificial or even underground. The long and West -East oriented pipe located in the South of the Basin (Fig. 1, left) is actually the Bièvre River.

2.3) Fractal dimensions of the impervious surfaces and of the sewer system

The studied catchment is located in a semi-urbanised area. The impervious surfaces are highly relevant for hydrology since they basically correspond to areas where runoff is quickly active during a storm event. Thanks to the determination of land use per pixel in MH-AssimTool, the evaluation of the impervious areas can be done in an apparent simple way by calculating the number of pixels of roads, buildings and gullies (since the water falling on gully pixels is immediately routed to the sewer network, they are considered as impervious).

This impervious surface depends on the resolution at which it is computed. Indeed, an imperviousness of 55%, 50%, 42%, 32% and 25% is obtained with pixels of size 20, 15, 10, 5, 2 m respectively. This is due to the priority order set in the data assimilation tool that affects a land use for each pixel. This order prioritizes impervious areas (Fig. 2). Obviously these values strongly depend on the approach implemented to affect a land use class to a pixel. As previously mentioned, comparison with other approaches can be found in Ichiba et al. (2017). Investigations on the possibility of having different pixel size according to the land use should also be envisaged in the future, in order to for instance refine the pixels for roads and gullies and coarser them for forests. Coming back to the imperviousness percentages found in this paper, it is possible to use the notion of fractal dimension, which is scale invariant, to explain these figures. The fractal dimension D F of a geometrical set (here the impervious pixels) is obtained with the help of the following equation:

F D N λ λ ≈ Eq. 1
where N λ is the number of impervious pixels, and λ is the resolution defined as the ratio between the outer scale L of the phenomenon and the observation scale l (

l L = λ
). It characterizes the space occupied by a geometrical set in a scale invariant way. The symbol ≈ denotes an asymptotic convergence and absorbs slowly varying prefactors.

For the studied catchment, it appears that the impervious areas exhibit a fractal dimension. Indeed Eq. 1 is plotted in log-log for the geometrical set consisting of the impervious pixels at the 2-m resolution (imperviousness of 25 %), and a straight line is retrieved on the whole range of scales, i.e. 2m-2048m (Fig. is simply a consequence of the priority order set for affecting a land use class to a pixel (impervious classes are prioritised over pervious ones). This confirms the fact that even though the represented imperviousness varies with scale, a feature (the fractal dimension) is conserved and provides a quantification of the level of urbanisation. We find D F equal to 1.73 for this catchment. In a previous study Gires et al. (2014a), found that for a highly urbanised area in Seine-Saint-Denis (North-East of Paris), the fractal dimension was of 1.85 from on scales ranging from 1 m to 1024 m. Given that this catchment is less urbanised, it was expected to obtain a smaller fractal dimension.

The same study was performed on the sewer system (Fig. 3.b). In this case, the geometrical set studied is the "rasterised" sewer system. If a pixel is crossed by a conduit belonging to the storm water sewer network, then it is considered as part of the sewer system..

Two scaling regimes can be identified: from 10 m to 80 m the fractal dimension is 1.03 and from 80 m to 1280 m it is 1.76. For small scales, the dimension is close to 1, and it simply reflects the 1D intrinsic nature of the sewer system. For large scales, the structure of the network becomes apparent, and exhibits a scaling behaviour. For large scales the value is slightly smaller than the 1.85 found on the Seine-Saint-Denis catchment in Gires et al.

(2014a) which is consistent with the fact that this one is less urbanised. The similarity between both fractal dimensions (imperviousness and large scale sewers) indicates that it is a relevant way of quantifying a level of urbanisation for the area. See [START_REF] Gires | Fractal analysis of urban catchments and their representation in semi-distributed models: imperviousness and sewer system[END_REF] for an extension of this approach to 10 areas in 5 European countries.

2.4) Rainfall data

Four rainfall events, which occurred between 2009 and 2011, are studied in this paper.

Simulations are performed using Météo-France radar mosaic, which provides a spatially distributed data with a resolution of 1 km x 1 km x 5 min (the closest radar is the C-band one of Trappes located 15 km West). For three events the data recorded with the help of a rain gauge operated by the SIAVB located a few hundred meters south of the catchment at the "Bassin des Bas Près" is also available. Because of (i) the standard 0.2 mm discretization issue of the tipping bucket rain gauge (data is number of tips equal to 0.2 mm) which prevents it from providing reliable intensity, (ii) the gap between the observation scales of the two measuring devices (see Gires et al., 2014b, for an in-depth analysis of this issue) and (iii) the fact that the rain gauge is furthermore outside of the catchment; it is not possible to use the rain gauge data for other purpose than a rough check of the accuracy of radar data. It is done by comparing the cumulative volumes of rainfall for each studied event which are displayed in Table 1 along with their main features. Gires et al. (2014b) used data from dense network of point measurement devices (rain gauges or disdrometers) distributed over 1 km 2 and showed that the cumulative depth differences between devices could reach more than 40 % for individual rainfall events (of the same order of magnitude as the one discussed here). They showed with the help of numerical simulations that similar values were found simply taking into account small scale rainfall variability. Here the maximum observed differences are 34%, which suggests that the agreement between the two devices is acceptable, i.e. smaller than expected uncertainty simply due to the scale gap between the two measuring devices. Authors did not have access to longer time series of both radar and rain gauge to perform a more indepth evaluation of the radar versus rain gauge measurement for this specific point, which would be the topic of another study. The temporal evolutions of the radar rain rate averaged over the catchment are displayed in Fig. 4. These events were selected because they are heavy ones. However they are not extreme ones, indeed over durations of 1 h and 4 h, only the 14

July 2010 event has a return period greater than 1 year (data from a rain gauge located in the Paris area that was available to the authors was used to obtain these estimates). For the July event, the return period is of about 1 year for a duration of 1 h and of about 2 years for a duration of 4 h.

3) Methods

3.1) Multifractal framework

The Multifractal framework is used for several purposes throughout this paper to characterize the variability across scales of fields, and is therefore presented here in a generic way. Only basic properties are discussed here, and interested readers are referred to the recent review by [START_REF] Schertzer | Multifractals, generalized scale invariance and complexity in geophysics[END_REF] for more details. The general assumption of multifractal fields is that they are generated by an underlying scale invariant multiplicative cascade process. In such process, a structure at a given scale is divided into smaller structures at smaller scale and the value of a child structure is equal to the value of the parent structure multiplied of a random increment. The process is scale invariant in the sense that the way structures are divided into sub-structures and the probability distribution of the random multiplicative increments are the same at all scales. A consequence is that statistical properties of such fields are conserved across scales. More precisely let us denote ε λ a field at resolution λ (=L/l, where l is the observation scale and L the outer scale of the phenomenon as for the fractal dimension definition). The probability of exceeding a given threshold (λ γ ), defined with the help of the scale invariant notion of singularity γ (the thresholds depend on the observation scale, but not the singularity),

) ( ) Pr( γ γ λ λ λ ε c - ≈ ≥ Eq. 2,
and the moment of order q, K(q) q λ λ ε ≈ Eq. 3, exhibit a power law relation with regards to the resolution at which they are computed. As for Eq. 1, the symbol ≈ denotes an asymptotic convergence and absorbs slowly varying prefactors. Equations 2 and 3 define respectively the codimension function c(γ) and the moment scaling function K(q), which both fully characterize the variability across scales of the field. c(γ) and K(q) contain the same information and are related by a Legendre transform [START_REF] Parisi | A multifractal model of intermittency[END_REF]. Eq. 2 can be understood from the simpler notion of fractal dimension (Eq. 1). Indeed, an intuitive interpretation of a multifractal field is that the geometrical sets made of each portion of the field greater than given thresholds are fractal and characterized by fractal dimensions. To be mathematically more rigorous the notion of threshold is replaced by the scale invariant one of singularity.

By generalizing the central limit theorem [START_REF] Schertzer | Physical modelling and analysis of rain and clouds by anisotropic scaling and multiplicative processes[END_REF] showed that any conservative scale-invariant multiplicative processes converge toward Universal Multifractals (in a similar way as re-normalized sum of identical and independent random variables converge toward normal distribution as long as their variance is defined). For Universal Multifractals (UM), i.e. this limit behaviour, K(q) and c(γ) functions are defined with the help of only two relevant parameters with a physical interpretation. They are known as UM parameters C 1 and α:

-C 1 is the mean intermittency which measures the average sparseness of the field. C 1 =0 for a homogeneous field.

α is the multifractality index ( 2 0 ≤ ≤ α

) and measures how fast the intermittency evolves when considering level of activity slightly different from the average one.

Great values of α and C 1 corresponds to strong extreme. A common tool to assess the extremes of a field is the scale invariant notion of maximum probable singularity γ s observable [START_REF] Hubert | Multifractals and extreme rainfall events[END_REF][START_REF] Douglas | Probable maximum precipitation estimation using multifractals: Application in the eastern United States[END_REF][START_REF] Royer | Multifractal analysis of the evolution of simulated precipitation over France in a climate scenario[END_REF]Gires et al., 2014a). It is defined for a unique sample by

d c s = ) (γ Eq. 4
Where d is the dimension of the embedding space, i.e. 1 = d for time series and 2 = d for maps.

The power spectrum (Fourrier transform of the auto-correlation function) of such multifractal field exhibits a scaling relation with wave number k:

β - ≈ k k E ) ( Eq. 5
where β is the spectral slope.

3.2) Uncertainty associated with small scale rainfall

The purpose of this section is to explain the approach implemented to quantify the uncertainty associated with small scales rainfall variability, i.e. which is occurring below the 1 km x 5 min scale currently provided by the C-band radar operating in this area. The same methodology as in [START_REF] Gires | Multifractal analysis of a semi-distributed urban hydrological model[END_REF]Gires et al. ( , 2014a) ) where

radar PF t ,
is the time of the peak flow simulated with the raw radar data (PF radar )..

3.3) Multifractal analysis of overland water depth maps

There is no distributed data available for overland water depth over large areas, but it is possible to study the fields obtained with the help of numerical simulations with spatially distributed rainfall as input. Maps of water depth during runoff at the end of each 3 min Multi-Hydro loop are studied.

Technically, in this paper an area of 128 x 128 pixels of size 10 m x 10 m is extracted from the map of the catchment to carry out the analysis. Both ensemble analysis (i.e. considering all successive maps as independent realisations of the same process and upscaling them individually before taking the mean in Eq. 2 and 3) and individual time step analysis (i.e. to obtain temporal evolutions of the various parameters) are performed. Finally, analyses are done in 2D on the maps but also in 1D on the columns or the lines of pixels over the catchment, in a North-South direction and in an East-West direction respectively (Fig. 5). The purpose of this is to monitor a possible influence of the slope over the generated runoff scaling properties. 4) Implementation of the Multi-Hydro model on the Jouy-en-Josas catchment

4.1) Validation with raw radar data

The validation of the model is achieved by comparing the water height measured at the Pont-de-Pierre gauge with the simulated one. Before going on authors would like to highlight that a proper validation on this case study is not possible given the available data, and will therefore limit this section to checking that the model approximately behaves well. The main reasons for this problem are:

-Only one measuring point is available for the whole catchment taking into account approximately an area of 2 km 2 .

-The uncertainty associated with this water level gauge is high. Indeed, it is not operated for accurate hydraulic measurement but to trigger an alarm to evacuate a music school located nearby. The main issue is that the shape of the river bed cross section at this point is not available. The width was estimated at around 1.80 m, using aerial photography from IGN and an approximate measure from few meters away. In order to correctly model the pipe, we used Multi-Hydro and tested various types of conduits. Finally, we chose to model the Bièvre as a circular pipe, with free surface of 2 m diameter, which is close to the approximate measurement. This choice is only an approximation which does not take into account the variations in time of this shape due the fact that the bottom of the river bed is not flat and contains moving rocks and changing vegetation.

-There is a lack of available data on initial soil saturation which is one of main sources of uncertainty and can biased runoff (see [START_REF] Shah | Modelling the effects of spatial variability in rainfall on catchment response. 1. Formulation and calibration of a stochastic rainfall field model[END_REF]Zehe et al., 2005) especially at the beginning of the event. In this paper, dry conditions were considered at the beginning of each event. A sensitivity test was conducted by considering a saturated soil at the beginning. A slight increase (few percent) of simulated flow was noted only during approximately the first hour (not shown here). Having longer rainfall time series would enable to simulate the catchment's behaviour some time before the event and limit the uncertainties associated with this issue.

-The uncertainties on the water input in the Bièvre River at the outlet the Bas-Près storage basin upstream the catchment are not quantified.

-Obviously there are some uncertainties on the radar rainfall measurement itself.

The simulation and measurement at the "Pont de Pierre" point for the selected rainfall events are displayed in Fig. 6. For the 09-02-2009 event we observe a clear overestimation at the beginning of the event. For the 14-07-2010 event Multi-Hydro with the radar rainfall data reproduces well the two main peaks, but overestimates the first local maximum of rainfall intensity and misses the second one. The 15-08-2010 event shows a greater variability in the first half of the simulation (variations are more pronounced on the model than on the measurements) but reproduces well the last peak. Finally, for the 15-12-2011 event, the Multi-Hydro model reproduces well the first peak, but the flow decreases more rapidly than the observations. Given the available data on a limited number of events it is difficult to attribute the observed discrepancies to one or several of the previously mentioned sources of uncertainty.

Proper validation would indeed require the analysis of much longer time series and more accurate measurements with better position of sensors. Nevertheless, the obtained results do not highlight strikingly wrong behaviour of simulated water heights in conduit, and enable to partially reproduce observations. Finally, it seems that for some events the simulated flows might be too noisy compared with observed water levels. This should not affect the UM analysis that follows because the analyses carried out in this paper are spatial ones, i.e. maps are studied and not time series so the potential effect should be limited. Keeping in mind the previously mentioned limitations, results suggest that it remains relevant to use this implementation of Multi-Hydro with a rather coarse 10 m resolution for testing its sensitivity to small scale rainfall variability and analysing surface runoff with the help of multifractals.

The authors acknowledge that further investigations on other catchments with more accurately validated models would be needed to fully confirm the findings discussed after.

4.2) Uncertainty associated with small scales rainfall variability

The envelop curves Q 0.05 , Q 0.25 , Q 0.75 , and Q 0.95 are displayed in Fig. 7 for the 09-02-2009 event for 5 conduits selected from upstream to downstream, which enables to analyse the effect of the position of the conduit within the network. Link #4 corresponds to the Pontde-Pierre measurement, and #5 to the outlet of the catchment. As it can be seen in Fig. 7, link #4 and #5 are located along the Bièvre River, and they take into account the significant base flow in the river coming from upstream the Jouy-en-Josas catchment. It means that they are obviously less sensitive to local rainfall variability. Similar curves were also generated for water height (not shown) at the Pont-de-Pierre. The computed uncertainty is small and certainly does not explain the discrepancies between simulations and measurements noticed in Fig. 6, which are hence not simply due to effects of small scale rainfall variability.

CV' 95 and CV' 75 values computed for the selected conduits (Fig. 7) and the four events are displayed in Table 2. As expected they decrease while considering more and more downstream conduits. There is a sharp decrease in CV' when the Bièvre River is reached because the base flow of the river dampens the effect of local small scale rainfall variability occurring over the 3 km 2 catchment, but the uncertainty only associated with this effect remains of roughly 10 % at the outlet whatever the event. The values for up-stream and midstream pipes are great for all events, even for CV' 75 which highlights a significant impact of small scale rainfall variability on the simulated flow. The variability observed in the simulated flow is basically due to the disparities in the simulated downscaled rainfall fields which are transferred through the hydrological model (see [START_REF] Gires | Influence of the zerorainfall on the assessment of the multifractal parameters[END_REF] for more detailed analysis of this issue). Small scale rainfall data is needed to understand better, and plan better, some local flooding due to sewer overflows which have been reported in some areas, notably the street parallel to the Bièvre River bed in the city (just North of it), There does not seem to have a straightforward relation between the computed uncertainty and the strength of the event (in terms of maximum rainfall peak intensity over 5 min). Indeed, the tendency that could be observed on the 09-02-2010, 15-08-2010 and 15-12-2011 (not a linear one as for example the peak rainfalls are equal to approximately 7 and 24 mm.h -1 for respectively the 15-08-2010 and 15-12-2011 event while the computed uncertainties are close) is not confirmed by the results for the 14-07-2010 event (see Tab. 2). Finally, these values are comparable to the ones that were obtained on a 1.5 km 2 highly urbanised catchment located 40 Km North-East on the other side of the Paris area in Gires et al. (2014a). For this catchment, CV' 95 values were ranging from 21 to 56%, 26 to 94% and 22 to 50% from downstream to upstream for the same 09-02-2009, 15-08-2010 and 15-12-2011 events respectively (at a different location). The values are slightly smaller for this catchment and this is likely to be due to lower level of imperviousness resulting in a smaller portion of rainfall becoming immediately active.

In this paper, the uncertainty is computed not only on the simulated flow, but also on the water depth in streets. As for the flow, for each realisation of downscaled rainfall field, the maximum water depth over the whole simulation is retrieved for each pixel. A sample is shown in Fig. 8.a for the 15-12-2011 event. The known hot spots are visible, although with too high values. For example, the modelled maximum water depth reaches more than 15 cm in the street along the Bièvre River bank in the city and the parallel street just north of it (already mentioned in the previous paragraph). Although some flooding is regularly reported by citizens to the SIAVB for these streets, such height was not reported for this event. In the urbanized portion of the catchment the street network is visible on the maximum water depth map, meaning the maximum values of water depth maps are reached on the corresponding pixels. Lower values are found on the on roads/streets located on the steep portion of the catchment because water moves faster in these areas. Same patterns and numerical values are obtained for other realisations of the same event. Similar plots are obtained for the other events with lower depths for the 09-02-2009 and 15-08-2010 (for which a lower cumulative rainfall depth was recorded) and greater depths for the 14-07-2010 event. Then, as for the flow analysis previously carried out, the uncertainty on this maximum water depth is computed with the help of the 5 and 95% quantiles for each pixel and a pseudo-coefficient of variation. Illustrations of the quantiles maps are shown in Fig. 8.b and 8.c for the 15-12-2011 event. Similar patterns are observed on the two maps, notably for the hotspots previously mentioned which are visible on both maps. Maps of CV' 95 for maximum depth are displayed in Fig. 9 for the four rainfall events. It appears that the uncertainty is lower for the areas where the greatest maximum depths are found (i.e. on roads) and is also lower for the heaviest rainfall events. It reaches only few percents on the hottest points. The values (Fig. 9) are anyway much smaller than those found for sewer flow (Table 2 and Fig. 7). This apparent contradiction is likely to be due to the fact that most of the rain water is properly handled by sewers and overflows are limited for these events. It means that for these events disparities in local amounts will not be visible on ground levels, whereas they are indeed in sewer flows and water depths. Further investigations with heavier rainfall events should be carried out to confirm or not this interpretation. The areas with the greatest uncertainty are found in gardens for the weakest event (09-02-2009), and correspond to places with a very small maximum depth (smaller that 1mm), meaning that the hydrological relevance is not very high.

5) Multifractal characterization of overland water depth

Multifractal analyses of overland water depth during rainfall event are presented in this paper for the 14-07-2010 and 15-12-2011 events which are the two heaviest ones in terms of maximum rainfall intensity over 5 min (see Table 1).

Figure 10.a displays the spectral analysis of the water depth for the 14-07-2011 event.

Maps of water depth for each time steps during the event are used to carry out 2D ensemble analyses. The quality of the scaling is low, with a coefficient of determination for the linear regression equal to 0.42. The fact that the spectral slope is close to zero (β is found roughly equal to 0.2) indicates that the field is conservative, i.e. its mean is conserved across scales. It is therefore possible to implement directly on the field a Trace Moment (TM) analysis, which consists in assessing the validity of Eq. 4 by plotting it in log-log. Perfect UM fields would lead to straight lines. Figure 10.b shows the TM ensemble analysis performed over all the time steps of the same 14-07-2011 event. Two scaling regimes can be identified: a small scales regime from 10 m to 80 m (right part of Fig. 10.b) and a large scales regime from 80 m to 1280 m (left part of Fig. 10.b). The coefficient of determination r 2 of the linear regression for q=1.5 in Fig. 10.b is taken as an indication of the quality of the scaling. The scaling from small scales (10 m -80 m) is much more robust than for large scales (80 m to 1280 m), as illustrated by the r 2 equal to respectively 0.99 and 0.91. Given the low quality of the scaling for large scales, UM parameter estimates will not be reported and discussed for this regime because they are not reliable. Furthermore, small scales are crucial for surface runoff because it is at these scales that it is generated into the drainage system. The location of this break at approximately 80 m indicates a possible physical interpretation. Indeed, it is the same location as the break in the fractal analysis of the sewer network and corresponds roughly to the interdistance between roads. This would mean that this break is driven by the influence of the collection of water by sewer network. The more robust scaling behaviour for surface flow is found for the scales for which the sewer network does not behave yet as network but as isolated linear pipes. Before going on, it should be mentioned that numerous pixels have very small depth (see Fig. 8 for an illustration), for which the model uncertainties might be great.

These zeros values or spurious ones close to zero will affect the scaling analysis for small moments (typically q < 0.5) through a multifractal phase transition (see [START_REF] Gires | Influence of the zerorainfall on the assessment of the multifractal parameters[END_REF] a detailed analysis of this issue). Here the influence of this bias does not extend to moments close to 1 around which the estimates of UM parameters are carried out, meaning that they are not affected by this issue.

Although intrinsically less robust since scaling properties are statistical ones requiring numerous data to be properly observed, TM analyses were also carried out independently on each sampling time step of Multi-Hydro (3 min in this paper). The purpose is to see whether there is an impact of the current rainfall rate on it. Figure 11 displays for the 14-07-2010 event the temporal evolution of both the rainfall rate and the r 2 for q=1.5 in the TM analysis for the two regimes identified in the ensemble analysis, i.e. small (10 m -80 m) and large scales (80 m -1280 m). For this event, two rainfall peaks are observed, and they both result in a sudden loss of the scaling quality, more pronounced for large scales than small ones. For the first peak (yellow bars on Fig. 11) the decrease of r 2 lasts approximately 20 min, while it lasts only few minutes for the second peak (red bars on Fig 11). In both cases the quality of the scaling behaviour improves again over few tens of minutes. The physical meaning of such loss is not clear, but could be due to a bad representation of the surface flow process during intense rainfall (it might take some time to retrieve a realistic surface flow simulation following a sudden change in rainfall input), a bias in the geometrical repartition, or an intrinsic feature of the process. For the latter, a possibility is that during intense rainfall period, the surface flow exhibits more directly the rainfall features than its intrinsic ones which are retrieved once the flow process has "adapted" to the new conditions. This would explain both the loss of scaling quality and why scaling properties closer to rainfall ones are observed during these short periods. Analysis with a higher resolution model would be needed to further investigate this issue, which would also enable to have access to a wider range of small scales.

Similar features are retrieved for the other studied event (15-12-2011). Finally, it should also be mentioned that similar results are also found when performing the analysis on the North-South or West-East 1D-samples, which means that the preferential slope of the catchment (North-South) does not seem to have an influence on the scaling features of the simulated water depth. In terms of scaling quality, very similar results are also found with raw radar data, or downscaled rainfall fields suggesting a limited impact of small scale rainfall variability on these features. The same downscaling process as in section 2 is used.

UM parameters retrieved on the maximum water depth were computed for small scales, and are displayed in Table 3 for the two events (14-07-2011 and 15-12-2012) and for simulations with raw radar data and also a realisation of downscaled rainfall field with α=1.8 and C 1 =0.1 (other realisations yield very similar results). The temporal evolutions of α and C 1 for the 14-07-2011 event are shown in Figures 12.

It appears that the UM parameters are also affected by the "jumps" that were noticed on r 2 in Fig. 11. Indeed after an intense period, sharp increase of α and decrease of C 1 are noticed. These pronounced variations mean that the values obtained with ensemble analyses should not be over-interpreted. Nevertheless few comments can be made. First the values of C 1 are much greater than the ones reported for rainfall (typically 0.1-0.3 at small scale) meaning that significant levels of water depth are much more concentrated than the rainfall field, which reflects the influence of the physical processes associated with surface flow on the transferred field, notably the flow concentration. The most relevant one is the topography that routes water through specific paths and tends to concentrate it. Second the values of UM parameters are quite different between the two events. These differences are much greater than the ones observed on the rainfall fields (see [START_REF] Ichiba | X-band radar data and predictive management in urban hydrology[END_REF], for a detailed analysis of these storms) at small scales. This suggests that the large scales rainfall pattern has a strong influence on the retrieved parameters. Indeed, the topography and small scale rainfall features are the same between the two simulations; the only difference is the large scale rainfall features. Thirdly the values of γ s are rather similar for both events (the differences between α and C 1 tend to compensate themselves).

The temporal evolutions of the UM parameters obtained by inputting raw and downscaled rainfall data are very similar. The differences are slightly more pronounced on the values computed on ensemble analysis but as previously said this should not be overinterpreted given the strong variations visible in the temporal analysis. This similarity highlights the low influence of small scale rainfall variability on the retrieved parameters which seems to be more dependent on features associated with surface flow process itself or large scale rainfall.

In order to test the sensitivity of the results to small scales rainfall features, synthetic rainfall fields with various sets of known parameters are used as input to Multi-Hydro simulations. More precisely the pseudo-events tested last 30 min with an average intensity of 10 mm/h. Three pairs (α; C 1 ) of parameters are tested: (1.8; 0.1), (1.8; 0.05), (1.4; 0.1). Figure 13 displays the temporal evolutions of the rain rates, r 2 , α and C 1 for water depth for the three synthetic rainfall events.

The temporal evolution shows the same general tendency as the one observed with the real events. A loss of scaling quality is observed during the event itself, and it improves afterwards. α and C 1 have a constant behaviour during the rainfall, while they decrease and increase respectively after the rainfall has stopped. The comparison of the UM parameters for the overland maximum water depth shows that they do not seem to depend on the small scale rainfall variability in this case. α is constant around 1.4 while C 1 is constant around 0.6 during the rainfall. γ s is again constant around 1.7 on average. The rainfall UM parameters do not seem to modify the structure of the overland flow, and its geometrical distribution. Successive simulations with the same parameters for synthetic rainfall yielded same results. A physical explanation of the C 1 parameter could be that during the rainfall, the surface flow is more homogenous due to a ubiquitous input of water. UM parameters on water depth are thus closer to the rainfall ones (small C 1 ). However after the rain has stopped, the disparities of simulated water depth are increased due to predominant pathways (roads) or topographic depressions where the water can accumulate. The greater C 1 after the event could reflect this fact. The smaller values of α mean that the disparities among the areas where water remains tend to decrease after the rainfall event.

This study seems to highlight the fact that UM parameters α and C 1 for water depth are rather relying on the large scale structure of the rainfall and on the catchment features, while the maximum observable singularity γ s is conserved for all events. Further studies could infirm or confirm the fact that γ s depends on the studied catchment. The temporal evolutions of the UM parameters also deeply rely on the rainfall rate. Synthetic events with block structures enabled to stand out rather simple general tendencies. They become more complex with real rainfall, when the intensity has a higher temporal variability.

The temporal evolutions of the UM parameters also enable to quantify a catchment response time. Due to the sampling time step of the simulations, the uncertainty associated with the response is of 3 min. Still, it can be noted that in urban catchments (or semi-urban here), the response time of water depth UM parameters to the beginning of a rainfall or to an important peak of intensity is almost non-existent. This is due to the presence of impervious area over which rainfall directly transfer into surface runoff.

6) Conclusions

The Multi-Hydro model was implemented on the Jouy-en-Josas catchment in the Paris area. This 3 km 2 semi-urbanised catchment exhibits sharp slopes, and a dense area along the Bièvre River bed. It has often been damaged by major pluvial and fluvial flooding, before the construction of storage basins along the river path. The model was validated on this new catchment on four rainfall events with the help of the data from a height gauge near the outlet.

Rainfall radar data with a resolution of 1km x 1km x 5min was used.

Then ensembles of downscaled rainfall fields were used to quantify the sensitivity of the model outputs to small scales unmeasured rainfall variability, i.e. occurring below the resolution of the available raw radar data. It appears that it is rather significant on flow simulated in conduits with pseudo coefficients of variations ranging from 90 % upstream to 10% downstream. This confirms previous results obtained on a 1.5 km 2 flat highly urbanised catchment also in the Paris area. The methodology was extended here to simulated water depth, and it was found that the sensitivity was much lower than for conduits' flow. This is likely to be due to the fact that the sewer system is mainly able to cope with the storm water for these events limiting the amount of surface runoff.

After using them to downscale the radar data, Universal Multifractals are used in an innovative way to characterize the surface flow process -through simulated water depth for each 10 m x 10 m pixel over 3 min time steps-during rainfall events. UM parameters α and C 1 , and the composite parameter γ s are evaluated on the outputs of Multi-Hydro. Two scaling regimes are identified for this field and estimates are only reliable for small scales, i.e. 10m -80m, and related to the fractal feature of the sewer system which exhibits a scale break at the same scale. There is a loss of the quality of the scaling during intense rainfall periods and UM parameters get closer to rainfall ones. A possible interpretation is that during this short period, a mixture of the scaling behaviour of both surface flow and rainfall is observed. After the event scaling is improved and features more specific to surface flow processes are retrieved with a field strongly concentrated and variability among the wet areas dampened (C 1 greater than 1 and α smaller than 1). Small scale rainfall features do not seem to strongly influence the results which depend more on large scales rainfall spatio-temporal patterns for these events which do not trigger much sewer overflow.

The conclusions found with the help of this innovative methodology are not as straightforward as the authors would have hoped. Further investigations with other rainfall events, other catchments, notably with denser monitoring network including in-sewer measurements, should be carried out to strengthen the results. Higher resolution models should also be tested to extend the range of available scales for the small scales regime to obtain more reliable estimates of scaling features. Such new analysis would enable to generalize the behaviour of the scaling and of the UM parameters which describes the surface flows, and eventually to link them to other geometrical features of the catchment, such as the fractal dimension of its impervious surface, of the roads (which are the preferential path for surface flows) or of the sewer system. This paper should be seen as a promising first step that hints at innovative techniques relying on scale invariance properties to analyse how the rainfall extremes are either dampened or enhanced by hydrological models and also to quantify the extremes at very high spatial resolution (typically 1 m) without having to run the model at these resolutions which would require too much time especially for real time applications.

  3.a). This a basic feature of the catchment. The fact that the points corresponding to the catchment representation at 20, 15, 10, 5, 2 m obtained with MH-AssimTool are along this straight line (circled cross on Fig. 3.a)

  is implemented, and only basic ideas are explained here. Firstly, an ensemble of downscaled rainfall fields is generated, then each realisation is inputted into the numerical model and finally the disparities within the ensemble of outputs, which reflect the studied uncertainty, are analysed and quantified. 100 sample ensembles are used. The downscaling technique relies on the Universal Multifractal framework. It basically consists in stochastically continuing a space-time cascade process that has been validated on the available range of scales. The resolution of the downscaled rainfall field is 12 m in space and 20 s in time starting from the original 1 km and 5 min of the available radar data. The process has been validated down to such small scales(Gires et al., 2014b).The disparities among the simulated ensembles are quantified with the help of quantile analysis. Let us first illustrate this with the flow output, but the same is done for maximum water depth at each pixel. For each time step the 5, 25, 75 and 95 % quantiles are computed, and give the envelop curves Q 0.05 , Q 0.25 , Q 0.75 , and Q 0.95 , respectively. The width between these curves characterizes the uncertainty interval on simulated flow. It is quantified with the help of two pseudo-coefficients of variation computed as:

Figure 1 :

 1 Figure captions:

Figure 2 :

 2 Figure 2: Map of the land use obtained with the help of MH-AssimTool over the Jouy-en-

Figure 3 :Figure 4 :

 34 Figure 3: (a) Evaluation of the fractal dimension of the impervious area for the studied

Figure 6 :Figure 7 :

 67 Figure 6: Water height simulated with the help of Multi-Hydro using raw radar data as rainfall

Figure 8 :

 8 Figure 8: For the 15-12-2011 event. (a) Map of the computed maximum water depth for a

Figure 9 :

 9 Figure 9: Map of CV' 95 (in %) for the maximum water depth for the 09-02-2009 (a), 15-08-

Figure 11 :

 11 Figure 11: For the 14-07-2010 event: (a) Temporal evolution of the rain rate; (b) Temporal evolution of the r 2 for q=1.5 in the TM analysis for the two regimes identified in Fig. 10.

Figure 12 :

 12 Figure 12: Temporal evolution of the UM parameters α and C 1 of the maximum water depth

Figure 13 :

 13 Figure 13: For three synthetic rainfall events with different sets of UM parameters; temporal

Acknowledgments

Authors acknowledge partial financial support from the European Interreg IV RainGain project (www.raingain.eu), the Chair "Hydrology for resilient cities" endowed by Veolia, and the Ile-de-France regional R2DS RadX@IdF project. Authors would also like to thank Hervé Cardinal from the SIAVB for providing data and fruitful discussions.

6) References

Tables:

Radar rain depth analysis with either raw radar data or a realisation of downscaled rainfall field (with α=1.8

and C 1 =0.1) as rainfall input for the 14-07-2010 and 15-12-2011 events.