Marco Dentz

Peter K Kang

Tanguy Le Borgne

Continuous Time Random Walks for Non-Local Radial Solute Transport

Keywords: Continuous Time Random Walks, Multirate Mass Transfer, Radial Transport, Random Walk Particle Tracking, Stochastic Modeling, Non-Local Transport

This paper derives and analyzes continuous time random walk (CTRW) models in radial flow geometries for the quantification of non-local solute transport induced by heterogeneous flow distributions and by mobile-immobile mass transfer processes. To this end we derive a general CTRW framework in radial coordinates starting from the random walk equations for radial particle positions and times. The particle density, or solute concentration is governed by a non-local radial advection-dispersion equation (ADE). Unlike in CTRWs for uniform flow scenarios, particle transition times here depend on the radial particle position, which renders the CTRW non-stationary. As a consequence, the memory kernel characterizing the non-local ADE, is radially dependent. Based on this general formulation, we derive radial CTRW implementations that (i) emulate non-local radial transport due to heteroge-

Introduction

Solute transport in heterogeneous porous media displays behaviors that cannot be captured by transport models based on an equivalent advection dispersion equation (ADE) parameterized by (constant) effective transport parameters. Such behaviors range from the non-linear evolution of solute dispersion to power-law tails in solute breakthrough curves [START_REF] Bouchaud | Anomalous diffusion in disordered media: Statistical mechanisms, models and physical applications[END_REF][START_REF] Berkowitz | Modeling non-Fickian transport in geological formations as a continuous time random walk[END_REF]. The last three decades have seen intense research to quantify these behaviors in terms of effective transport models that can be obtained by moment equation ap-proaches [START_REF] Neuman | Stochastic-theory of field scale Fickian dispersion in anisotropic porous media[END_REF], and projector formalisms [START_REF] Cushman | Nonequilibrium statistical mechanics of preasymptotic dispersion[END_REF], for example, and include time and space fractional ADEs [START_REF] Benson | Application of a fractional advection-dispersion equation[END_REF][START_REF] Benson | Radial fractional-order dispersion through fractured rock[END_REF], multirate mass transfer (MRMT) models [START_REF] Haggerty | Multiple-rate mass transfer for modeling diffusion and surface reactions in media with pore-scale heterogeneity[END_REF][START_REF] Carrera | On matrix diffusion: formulations, solution methods, and qualitative effects[END_REF], as well as continuous time random walks [START_REF] Berkowitz | Anomalous transport in random fracture networks[END_REF][START_REF] Dentz | Time behavior of solute transport in heterogeneous media: transition from anomalous to normal transport[END_REF], see also the reviews in [START_REF] Berkowitz | Modeling non-Fickian transport in geological formations as a continuous time random walk[END_REF][START_REF] Neuman | Perspective on theories of anomalous transport in heterogeneous media[END_REF][START_REF] Dentz | Mixing, spreading and reaction in heterogeneous media: A brief review[END_REF].

In this paper, we focus on the CTRW approach to modeling non-Fickian solute transport in heterogeneous media. Classical random walks model particle movements by using variable spatial steps which are taken within constant time increments at equidistant times [START_REF] Delay | Simulating solute transport in porous or fractured formations using random walk particle tracking: A review[END_REF][START_REF] Salamon | A review and numerical assessment of the random walk particle tracking methods[END_REF]. A CTRW, in contrast, models particle movements in a heterogeneous medium effectively as a random walk in which both space and time increments are variable. The spatial transitions may reflect the geometry of the underlying medium and flow heterogeneity, while particle transition and waiting times reflect persistent particle velocities over given transition distances, or particle retention due to adsorption or diffusion into immobile zones, for example [START_REF] Berkowitz | Anomalous transport in random fracture networks[END_REF][START_REF] Le Borgne | Lagrangian statistical model for transport in highly heterogeneous velocity fields[END_REF][START_REF] Dentz | Effective transport dynamics in porous media with heterogeneous retardation properties[END_REF][START_REF] Dentz | Distribution-versus correlation-induced anomalous transport in quenched random velocity fields[END_REF][START_REF] Kang | Spatial Markov model of anomalous transport through random lattice networks[END_REF].

The medium heterogeneity is mapped into the probability distribution density (PDF) of characteristic particle transition times. The evolution of the particle density, or, equivalently the solute concentration is governed by a temporally non-local ADE whose memory kernel is given in terms of the PDF of transition times [START_REF] Berkowitz | Physical pictures of transport in heterogeneous media: Advection-dispersion, random-walk, and fractional derivative formulations[END_REF][START_REF] Dentz | Time behavior of solute transport in heterogeneous media: transition from anomalous to normal transport[END_REF].

The MRMT approach is phenomenologically similar to the CTRW modeling framework as it models the impact of medium heterogeneity on large scale transport through a distribution of typical solute retention times in immobile regions. In fact, it can be shown [START_REF] Schumer | Fractal mobile/immobile solute transport[END_REF][START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF] that one model can under certain conditions be mapped onto the other. The latter amounts essentially to identifying the relation between the PDF of particle transition times and particle retention times in immobile regions [START_REF] Margolin | Continuous time random walk and multirate mass transfer modeling of sorption[END_REF][START_REF] Benson | A simple and efficient random walk solution of multi-rate mobile/immobile mass transport equations[END_REF][START_REF] Dentz | Diffusion and trapping in heterogeneous media: An inhomogeneous continuous time random walk approach[END_REF]. As pointed out above, the CTRW model is a random walk approach in that particle movements are governed by random walk equations for the space and time coordinates. Therefore the solution of CTRW and equivalent models is directly accessible to numerical solution through random walk particle tracking simulations [START_REF] Dentz | Time behavior of solute transport in heterogeneous media: transition from anomalous to normal transport[END_REF]. This provides an avenue for the efficient simulation of transport in the presence of mobile-immobile processes [START_REF] Le Borgne | Non-Fickian dispersion in porous media: 2. model validation from measurements at different scales[END_REF][START_REF] Benson | A simple and efficient random walk solution of multi-rate mobile/immobile mass transport equations[END_REF][START_REF] Dentz | Diffusion and trapping in heterogeneous media: An inhomogeneous continuous time random walk approach[END_REF], for example, and for temporally non-local transport in general [START_REF] Zhang | Lagrangian simulation of multidimensional anomalous transport at the made site[END_REF].

Many formulations of the above models are for transport situations under uniform mean flow. Thus, for the interpretation of tracer tests under forced conditions they are only of limited applicability because the nonstationarity of the underlying flow field is not accounted for. Haggerty et al. [START_REF] Haggerty | Tracer tests in a fractured dolomite 2. analysis of mass transfer in single-well injectionwithdrawal tests[END_REF] use a Eulerian radial MRMT implementation to interpret radial single-well injection-withdrawal (SWIW) tracer tests in fractured dolomite.

Le Borgne and Gouze [START_REF] Le Borgne | Non-Fickian dispersion in porous media: 2. model validation from measurements at different scales[END_REF] used a CTRW based random walk implementation of MRMT to model tracer breakthrough data from SWIW tracer tests.

Benson et al. [START_REF] Benson | Radial fractional-order dispersion through fractured rock[END_REF] developed a fractional-order dispersion model in radial coordinates to model tracer tests under forced conditions. A general issue when interpreting field tracer data is to decipher the origin of the observed nonlocal transport behavior, which may range from mobile-immobile diffusive mass transfer processes to highly heterogeneous advective transport [START_REF] Becker | Interpreting tracer breakthrough tailing from different forced-gradient tracer experiment configurations in fractured bedrock[END_REF][START_REF] Kang | Impact of velocity correlation and distribution on transport in fractured media: field evidence and theoretical model[END_REF].

In the latter case, non-Fickian transport may be caused by a broad distribution of flow and transport velocities; the distribution of particle transit times depends on the flow rate and heterogeneity in the flow properties. In the former case, anomalous transport features are due to mass transfer be-tween mobiel and immobile zones; particle transition times may depend on the retention properties and geometries of the immobile regions.

Testing these different hypothesis requires non-local transport models, that integrate both diffusive and advective mass transfer processes in nonuniform flow conditions.

In this paper, we develop a general CTRW approach that allows for the modeling of non-local solute transport under radial conditions. The derivation from the space-time random walk equations gives directly the particle tracking method for its numerical solution. We present three non-local CTRW based radial transport implementations, for the modeling of heterogeneous advection, mobile-immobile mass transfer (MRMT), and the combination of both. To this end we review in Section 2 briefly the random walk formulation of general radial advective-dispersive transport. Section 3 then derives the general radial CTRW framework and defines the specific CTRW models. The model breakthrough curves then are analyzed in Section 4 using numerical random walk simulations and explicit analytical expressions for the asymptotic breakthrough behavior developed in Appendix B. In particular, we discuss the expected differences in non-Fickian transport behaviors induced by purely advective processes, purely diffusive processes, and the combination of these processes.

Radial Random Walks

The classical advection-dispersion equation (ADE) for the solute concentration c(r, t) in radial coordinates can be written as

∂c(r, t) ∂t + 1 r ∂ ∂r v(r)rc(r, t) - 1 r ∂ ∂r rD(r) ∂c ∂r = 0, (1)
where v(r) and D(r) are the radially dependent transport velocity and dispersion coefficient; r denotes the radial distance, t denotes time. We set the constant porosity equal to one, which is equivalent to rescaling time. The equivalent random walk particle tracking formulation is obtained by rewriting (1) in mass conservative form. Therefore, we define the conserved radial concentration as

p(r, t) = 2πrc(r, t). (2)
Notice that p(r, t) denotes the concentration per unit radial distance. Inserting the latter into (1) and rearranging terms we obtain the radial Fokker-Planck equation

∂p(r, t) ∂t + ∂ ∂r v(r) + D(r) r + D (r) p(r, t) - ∂ 2 ∂r 2 D(r)p(r, t) = 0, (3)
where D (r) denotes the derivative of D(r) with respect to r. The equivalent Langevin equation is given by [START_REF] Risken | The Fokker-Planck Equation[END_REF]

dr(t) dt = v[r(t)] + D[r(t)] r(t) + D [r(t)] + 2D[r(t)]ξ r (t), (4)
where ξ r (t) is a Gaussian white noise of zero mean and the correlation function ξ r (t)ξ r (t ) = δ(t -t ). Here and in the following, we employ the Ito interpretation [START_REF] Risken | The Fokker-Planck Equation[END_REF] of the Langevin equation (4). The particle density is given in terms of the radial trajectories as p(r, t) = δ[r-r(t)] , and by virtue of (2), we obtain for the concentration distribution

c(r, t) = δ[r -r(t)] 2πr . (5)
In the following, we will consider the case of [START_REF] Bear | Dynamics of fluids in porous media[END_REF]

v(r) = k v r , D(r) = αk v r , (6)
where α is dispersivity, and k v = Q/(2π) with Q the flow rate. Notice that more general radial dependences of flow velocity and dispersion can be considered within the approaches developed in the following. Here, we focus on the choice [START_REF] Benson | Radial fractional-order dispersion through fractured rock[END_REF]. With these definitions, the Langevin equation ( 4) simplifies to

dr(t) dt = k v r(t) + 2αk v r(t) ξ r (t). ( 7
)
The temporally discretized version of the radial Langevin equation is given by

r n+1 = r n + k v ∆t r n + 2αk v ∆t r n ξ n , (8)
where r n = r(t n ), t n = n∆t, and ξ n is a Gaussian random variable with zero mean and unit variance.

Radial Continuous Time Random Walks

The radial random walk particle tracking formulations developed in the following are based on the generalization of the radial random walk process (8) in terms of the continuous time random walk

r n+1 = r n + + √ 2α ξ n (9a) t n+1 = t n + τ n (r), (9b)
where is a constant transition length, and τ (r) a radially dependent, independently distributed random transition time with the probability density function (PDF) ψ τ (τ, r). Notice that the classical formulation (8) is obtained by setting τ n (r) = r/k v in [START_REF] Berkowitz | Anomalous transport in random fracture networks[END_REF]. The distribution of the spatial transition lengths ∆r = + √ 2α ξ n is denoted by ψ r (∆r). The mean and mean square displacements are given by ∆r = and ∆r 2 = 2α , where we disregard contributions of order 2 . Notice that the transition length are chosen such that α.

A straightforward application of the general CTRW framework [START_REF] Berkowitz | Modeling non-Fickian transport in geological formations as a continuous time random walk[END_REF] gives for the radial particle density p(r, t) the equation

p(r, t) = t 0 dt R(r, t ) ∞ t-t dτ ψ τ (τ, r), (10)
where R(r, t) denotes the probability per time that the particle has just arrived at the radius r. This equation can be read as follows: The particle density at the position r at a time t is given by the probability that a particle arrives there at an earlier time t times the probability that the next transition takes longer than t-t . The density R(r, t) satisfies the mass balance equation

R(r, t) = δ(r -r 0 )δ(t) + ∞ 0 dr t 0 dt R(r , t )ψ τ (t -t , r )ψ r (r -r ), ( 11
)
where r 0 is the initial radial particle position. Equations ( 10) and ( 11) can be combined into the radial generalized Master equation

dp(r, t) dt = ∞ 0 dr t 0 dt ψ r (r -r )τ k (r ) -1 M (t -t , r )p(r , t ) - t 0 dt τ k (r) -1 M (t -t , r)p(r, t), (12)
where the memory kernel M (t -t , r) is defined in Laplace space by M (λ, r) = λτ k (r) ψτ (λ, r)

1 -ψτ (λ, r) . ( 13
)
The Laplace transform is defined in [START_REF] Abramowitz | Handbook of Mathematical Functions[END_REF]. Laplace transformed quantities in the following are marked by a hat, the Laplace variable is denoted by λ.

We consider here a ψ r (r) that is sharply peaked about its mean value .

In this case, the generalized Master equation ( 12) can be localized in space, and a Taylor expansion of the integrand on the right side of [START_REF] Dentz | Mixing, spreading and reaction in heterogeneous media: A brief review[END_REF] gives the non-local radial Fokker-Planck equation

∂p(r, t) ∂t + t 0 dt ∂ ∂r k v r - ∂ 2 ∂r 2 αk v r M (t -t , r)p(r, t ) = 0. ( 14
)
We disregard contributions of order 2 . Substituting the radial concentration (2) into this equation, we obtain the non-local radial advectiondispersion equation

∂c(r, t) ∂t + t 0 dt k v r ∂ ∂r - αk v r ∂ 2 ∂r 2 M (t -t , r)c(r, t ) = 0. ( 15
)
Radial non-local partial differential equations such as ( 14) and ( 15) can be solved subject to given initial and boundary conditions, either numerically or semi-analytically using common numerical schemes [START_REF] Benson | Radial fractional-order dispersion through fractured rock[END_REF][START_REF] Silva | A general real-time formulation for multi-rate mass transfer problems[END_REF] and analytical methods. In this paper, we use random walk particle tracking to solve for the transport behavior described by these equations.

In the following sections, we present three different models of increasing complexity for the time increment τ (r). The first one emulates non-local transport due to advective heterogeneity, the second model solves radial transport under multirate mass transfer, the third model provides a random walk approach combining heterogeneous advection and multirate mass transfer into immobile zones in radial coordinates.

Heterogeneous Advection

We first consider a radial CTRW that represents non-Fickian transport originating from broad distributions of advective transit times. In radial flow conditions, the mean transit time increases linearly with the radial distance.

To represent this important property, we propose to scale the random transit time τ (r) with the characteristic radial transition time τ k (r) as

τ (r) = τ k (r)η, τ k (r) = r k v . ( 16
)
The time scale τ k (r) denotes the transition time over the space increment r under homogeneous flow conditions. The dimensionless random increments η are independent identically distributed according to ψ η (η). They reflect the non-dimensional fluctuations of the radial transport velocity due to spatial heterogeneity and are related to the inverse radial flow velocity.

Recently, Edery et al. [START_REF] Edery | Origins of anomalous transport in heterogeneous media: Structural and dynamic controls[END_REF] studied the relation of the distribution of hydraulic conductivity and the transition time distribution under uniform flow conditions. The distribution of inverse velocities, and thus transition times, may be obtained from estimates of the hydraulic conductivity distribution. Alternatively it may be modeled by a parametric model [START_REF] Kang | Impact of velocity correlation and distribution on transport in fractured media: field evidence and theoretical model[END_REF], whose parameters are adjusted from the observed breakthrough curves, which in turn may give insight into the flow and medium heterogeneity.

Thus, the transition times τ (r) = τ k (r)η are distributed according to

ψ τ (τ, r) = 1 τ k (r) ψ η [τ /τ k (r)]. ( 17
)
The memory kernel M (t, r) defined in the previous section can now be written as

M (r, λ) = Mη [λτ k (r)] ≡ λτ k (r) ψη [λτ k (r)] 1 -ψη [λτ k (r)] , (18)
where we used that the Laplace transform of ( 17) is given by ψ(λ, r) =

ψ η [λτ k (r)
]. Thus, it follows that the memory kernel in real time has the scal-

ing form M (t, r) = τ k (r) -1 M η [t/τ k (r)].
The generalized advection-dispersion equation ( 15) can then be written as

∂c(r, t) ∂t + t 0 dt k v r ∂ ∂r - αk v r ∂ 2 ∂r 2 M η [(t -t )/τ k (r)] τ k (r) c(r, t ) = 0. ( 19
)
Notice the exponential distribution

ψ η (η) = exp(-η) gives M η (η) = δ(η)
in ( 18) and thus the generalized radial advection-dispersion equation [START_REF] Berkowitz | Physical pictures of transport in heterogeneous media: Advection-dispersion, random-walk, and fractional derivative formulations[END_REF] reduces to the radial advection-dispersion equation ( 1) in a homogeneous medium.

Multirate Mass Transfer

The MRMT model considers solute transport under mass transfer between a single mobile zone and a series of immobile zones [START_REF] Haggerty | Multiple-rate mass transfer for modeling diffusion and surface reactions in media with pore-scale heterogeneity[END_REF][START_REF] Harvey | Temporal moment-generating equations: Modeling transport and mass transfer in heterogeneous aquifers[END_REF][START_REF] Carrera | On matrix diffusion: formulations, solution methods, and qualitative effects[END_REF][START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF][START_REF] Schumer | Fractal mobile/immobile solute transport[END_REF][START_REF] Silva | A general real-time formulation for multi-rate mass transfer problems[END_REF].

Here we derive a radial random walk approach that simulates mass transfer between mobile and immobile regions based on the radial CTRW ( 9) and the CTRW model presented in the previous section. An important difference with the advective CTRW model described in the previous section is that the distribution of trapping times does not depend on the radial position.

Mass transfer between mobile and immobile regions is modeled by a compound Poisson process for the transition time τ (r) following the works of [START_REF] Margolin | Continuous time random walk and multirate mass transfer modeling of sorption[END_REF], [START_REF] Benson | A simple and efficient random walk solution of multi-rate mobile/immobile mass transport equations[END_REF] and [START_REF] Dentz | Diffusion and trapping in heterogeneous media: An inhomogeneous continuous time random walk approach[END_REF]. The particle transition time τ (r) is split into a mobile time τ m (r), which is the time the particle needs to traverse the distance in the mobile portion of the medium, and a series of immobile times τ im,i , which measure the times the particles spent in the immobile portion of the medium.

The number of times n im a particle gets trapped during a mobile transition of duration τ m (r) is given by a Poisson distribution characterized by the mean γτ m (r), where γ is the rate by which particles get trapped in immobile zones,

P im [n, τ m (r)] = [γτ m (r)] n n! exp [-γτ m (r)] . (20)
Thus, the time increment τ (r) in (9b) associated to the spatial transition (9a) is given by

τ (r) = τ m (r) + n im i=1 τ im,i . (21)
The distribution of immobile times is denoted by p im (τ ). For simplicity, we consider here the situation that initially all particles are mobile. The mobile times are given by τ m (r) = τ k (r)η as in [START_REF] Dentz | Effective transport dynamics in porous media with heterogeneous retardation properties[END_REF], with an exponentially distributed η. Notice that, as pointed out at the end of the previous section, an exponential η models transport in a homogeneous medium. This means, τ m (r) is exponentially distributed with mean τ k (r). Thus the transition time PDF can be expressed in terms of its Laplace transform as (see Appendix

A) ψτ (λ, r) = 1 1 + λτ k (r) + γτ k (r)[1 -pim (λ)] . ( 22
)
Inserting this expression into (13) for the Laplace transform M (λ, r) of the memory kernel gives the compact expression

M (λ, r) ≡ M (λ) = 1 1 + γλ -1 [1 -pim (λ)] . (23)
Note that the memory kernel is independent on the radial distance, which reflects the difference between mobile particle transitions, which depends on the local flow velocity, and particle retention in immobile zones, which depends on the distribution of diffusion times, for example. Using this expression in the Laplace transform of [START_REF] Le Borgne | Lagrangian statistical model for transport in highly heterogeneous velocity fields[END_REF] gives for the radial concentration

λĉ(r, λ) + k v r ∂ ∂r - αk v r ∂ 2 ∂r 2
ĉ(r, λ)

1 + γλ -1 [1 -pim (λ)] = c(r, t = 0). ( 24
)
We define now the mobile solute concentration by its Laplace transforms as

ĉm (r, λ) = ĉ(r, λ) 1 + γλ -1 [1 -pim (λ)] . (25)
Thus, we obtain from [START_REF] Dentz | Diffusion and trapping in heterogeneous media: An inhomogeneous continuous time random walk approach[END_REF] for ĉm (r, λ)

λĉ m (r, λ) + λ γλ -1 [1 -pim (λ)]ĉ m (r, λ) + k v r ∂ ∂r - αk v r ∂ 2 ∂r 2 ĉ(r, λ) = c m (r, t = 0). ( 26
)
Note that we assume for simplicity that initially all particles are mobile. An initial presence of particles in the immobile zones would give rise to a source term in (26) [START_REF] Tecklenburg | A non-local two-phase flow model for immiscible displacement in highly heterogeneous porous media and its parametrization[END_REF]. We furthermore define the density of immobile particles by the expression in curly brackets on the left side of [START_REF] Zhang | Lagrangian simulation of multidimensional anomalous transport at the made site[END_REF]. It reads in time space as

c im (r, t) = γ t 0 dt c m (r, t ) ∞ t-t dτ p im (τ ). ( 27
)
The right hand side expresses the density of immobile particles by the probability per time that mobile particles get trapped at a given time t , γc m (r, t ), times the probability that the residence time in the immobile region is larger than t-t . Thus, at asymptotically long times, the ratio of the time averaged immobile and mobile concentrations is given by

lim t→∞ c im (r, t) c m (r, t) = lim λ→0 γλ -1 [1 -p * im (λ)] = γ τ im , (28)
where τ im is the mean immobile time and the overline denotes the time average c(r, t) = t -1 t 0 dt c(r, t ). We now define the memory function ϕ(t) as

ϕ(t) = 1 τ im ∞ t dτ p im (τ ). ( 29
)
Thus, the governing equation ( 26) of the mobile solute concentration can be written in time space as

∂c m (r, t) ∂t + β t 0 dt ϕ(t -t )c m (r, t ) + k v r ∂ ∂r - αk v r ∂ 2 ∂r 2 c(r, t) = 0. ( 30
)
where we defined β = γ τ im . Equation [START_REF] Risken | The Fokker-Planck Equation[END_REF] describes transport under multirate mass transfer in radial flow [START_REF] Haggerty | Tracer tests in a fractured dolomite 2. analysis of mass transfer in single-well injectionwithdrawal tests[END_REF][START_REF] Le Borgne | Non-Fickian dispersion in porous media: 2. model validation from measurements at different scales[END_REF]. The memory function ϕ(t) encodes the mass transfer mechanism [START_REF] Haggerty | Multiple-rate mass transfer for modeling diffusion and surface reactions in media with pore-scale heterogeneity[END_REF][START_REF] Harvey | Temporal moment-generating equations: Modeling transport and mass transfer in heterogeneous aquifers[END_REF][START_REF] Carrera | On matrix diffusion: formulations, solution methods, and qualitative effects[END_REF][START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF][START_REF] Dentz | Effective non-local reaction kinetics for transport in physically and chemically heterogeneous media[END_REF] between the mobile and immobile regions. For linear first-order mass exchange it reflects the distribution of transfer rates between mobile and immobile regions [START_REF] Haggerty | Multiple-rate mass transfer for modeling diffusion and surface reactions in media with pore-scale heterogeneity[END_REF]. For diffusive mass transfer, it encodes the geometries and the characteristic diffusion scales of the immobile regions [START_REF] Noetinger | Up-scaling of double porosity fractured media using continuous-time random walks methods[END_REF][START_REF] Gouze | Non-fickian dispersion in porous media explained by heterogeneous microscale matrix diffusion[END_REF]. For transport through highly heterogeneous porous and fractured media, the memory function may be related semi-analytically or empirically to the medium heterogeneity [START_REF] Willmann | Transport upscaling in heterogeneous aquifers: What physical parameters control memory functions?[END_REF][START_REF] Zhang | Linking aquifer spatial properties and non-Fickian transport in mobile-immobile like alluvial settings[END_REF].

The memory function is here defined in terms of the distribution of residence times in the immobile regions given by [START_REF] Kang | Impact of velocity correlation and distribution on transport in fractured media: field evidence and theoretical model[END_REF]. Reversely, the distribution of residence times p im (τ ) can be obtained for a given memory function ϕ(t)

according to p im (τ ) = -τ im dϕ(τ ) dτ . (31)
Expressions ( 29) and ( 31) establish the relation between the distribution of residence times and the memory function.

Combination of Heterogeneous Advection and Multirate Mass Transfer

We consider now the random walk implementation of the combination of CTRW, which accounts for heterogeneous transport in the mobile zone, and multirate mass transfer between the mobile and immobile regions. Using the approach presented in the previous section, the mobile transition time τ m (r) in ( 21) is now given by the general relationship τ m (r) = τ k (r)η. The random variable η is distributed according to a general ψ η (η) as outlined in Section 3.1. Accordingly, the PDF for the mobile transitions is given by [START_REF] Dentz | Distribution-versus correlation-induced anomalous transport in quenched random velocity fields[END_REF] and reads as

ψ m (τ, r) = 1 τ k (r) ψ η [τ /τ k (r)]. (32)
Thus, the transition time PDF for the general compound [START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF] process reads in terms of its Laplace transform as (see Appendix A)

ψτ (λ, r) = ψm (λ + γ[1 -pim (λ)], r) . (33)
Specifically, by using the Laplace transform of (32), ψm (λ, t) = ψη [λτ k (t)],

the Laplace transform of ψ τ (τ, r) can be written as

ψτ (λ, r) = ψη (λτ k (r) + γτ k (r)[1 -pim (λ)], r) . (34)
Note that unlike for the mobile immobile model discussed in the previous section, the transition time PDF renders the memory kernel (13) dependent on the radial distance. The total trapping time during a mobile transition is related to the number of trapping events, whose mean is given by the mobile transition time times the trapping rate. While the individual trapping times are independent on advective heterogeneity, their collective is related to advective heterogeneity through the number of trapping event. This interrelation is reflected in the radial dependence of [START_REF] Edery | Origins of anomalous transport in heterogeneous media: Structural and dynamic controls[END_REF].

For times t γ -1 , or equivalently λ γ, the transition time PDF [START_REF] Edery | Origins of anomalous transport in heterogeneous media: Structural and dynamic controls[END_REF] can be approximated by the one for purely advective heterogeneity. Thus for t γ -1 , transport is dominated by advective heterogeneity, which is evident because at t γ -1 the number of trapping events is very small. If the PDF of dimensionless mobile transition times ψ η (η) has the finite mean transition time η < ∞, the transition time PDF (34) can be approximated for small arguments as

ψτ (λ, r) = 1 -τ k (r) η λ 1 + λ -1 γ[1 -pim (λ)] . (35)
Notice that this is valid as long as λτ k (r) + γτ k (r)[1 -pim (λ)] 1. Inserting [START_REF] Harvey | Temporal moment-generating equations: Modeling transport and mass transfer in heterogeneous aquifers[END_REF] into [START_REF] Delay | Simulating solute transport in porous or fractured formations using random walk particle tracking: A review[END_REF] gives in leading order for the memory kernel

M (λ, r) ≡ M (λ) = η -1 1 1 + γλ -1 [1 -pim (λ)] . (36)
Thus, in this limit, the model behaves as the MRMT model introduced in the previous section. These asymptotic behaviors of the different models are discussed in more detail in the next section, which studies solute breakthrough in the presented CTRW models.

Breakthrough Curves

We discuss here the behavior of solute breakthrough curves in the radial non-local models developed in the previous section for an instantaneous pulse tracer injection at radius r = r 0 . Notice that the response to a pulse has a fundamental character because responses to other injection conditions can be obtained by superposition of pulse responses. Furthermore, in field tracer experiments the duration of the tracer pulse is typically much shorter than the duration of the experiment and may be approximated as instantaneous [START_REF] Becker | Interpreting tracer breakthrough tailing from different forced-gradient tracer experiment configurations in fractured bedrock[END_REF][START_REF] Benson | Radial fractional-order dispersion through fractured rock[END_REF][START_REF] Kang | Impact of velocity correlation and distribution on transport in fractured media: field evidence and theoretical model[END_REF].

We use numerical random walk particle tracking simulations to solve for the breakthrough behavior in each of the models and derive explicit analytical expressions for their asymptotic behaviors. The numerical simulations are based on the recursion relations (9) for the particle positions and times.

Specifically, the initial position at n = 0 random walk steps is set to r 0 = and the initial time is set to t 0 = 0.

The solute breakthrough curve is identical to the distribution f (t, r) of first passage times of solute particles at the radius r. The first passage time τ a (r) is defined by

τ a (r) = t nr (37)
where n r = min(n|r n ≥ r) denotes the minimum number of steps needed to pass the radius r. The first passage time PDF then is given by

f (t, r) = ∞ n=0 f 0 (n, r)p(t, n), (38)
where f 0 (n, r) denotes the distribution of the minimum number of steps needed to exceed r, and p(t, n) the distribution of times after n random walk steps in (9b).

Heterogeneous Advection

The CTRW approach representing heterogeneous advection models the transition time by the time increment [START_REF] Dentz | Effective transport dynamics in porous media with heterogeneous retardation properties[END_REF], which is determined by the dimensionless increment η. In the following, we consider for ψ η (η) the Pareto distribution

ψ η (η) = βη -1-β , η ≥ 1, (39)
which models a broad distribution of transport velocities. This type of pure power-law behavior may be observed in an intermediate time regime. Asymptotically one would expect that the transition time PDF is truncated on a scale corresponding to a largest heterogeneity scale, for example, see also the discussion in [START_REF] Dentz | Time behavior of solute transport in heterogeneous media: transition from anomalous to normal transport[END_REF].

For a power-law distribution ψ η (η) of dimensionless transition times, we derive in Appendix B.1 the following scaling forms of the breakthrough curves for times larger than τ k (r)

f (t, r) = 1 θ(r) f 01 t θ(r) , 0 < β < 1 (40) f (t, r) = 1 θ(r) f 02 t -τ a (r) θ(r) , 1 < β < 2, ( 41
)
where we defined

θ(r) =   r/ i=1 τ k (r i ) β   1/β , τ a (r) = 1 β -1 r/ i=1 τ k (r i ). ( 42
)
We furthermore obtain for the asymptotic behavior of the scaling functions f 01 (x) and f 02 (x) for x 1 the power-law decay ∼ x -1-β such that we obtain for t θ(r) the following behavior for the breakthrough curves

f (t, r) ∼ 1 θ(r) t θ(r) -1-β . ( 43
)
For 0 < β < 1, mean and variance of the arrival time do not exist, for 1 < β < 2, the mean exists and is given by τ a (r) , while the variance diverges. The mean arrival time depends on distribution of transition times through β -1 in the denominator. The second term in the expression for τ a (r) is the mean arrival time for a homogeneous model characterized by an exponential PDF of dimensionless transition times η. For 0 < β < 1, the characteristic time θ(r) scales peak width and position, as can be deduced from the scaling form [START_REF] Willmann | Transport upscaling in heterogeneous aquifers: What physical parameters control memory functions?[END_REF]. For 1 < β < 2, it is a measure for the peak width, as indicated by [START_REF] Zhang | Linking aquifer spatial properties and non-Fickian transport in mobile-immobile like alluvial settings[END_REF]. Using τ k (r) = r/k v , we obtain for θ(r) and τ a (r)

θ(r) ≈ 2 k v 1 1 + β 1/β r 1+1/β , τ a (r) ≈ 1 β -1 r 2 2k v . ( 44
)
Note that the late time scaling (43) is the same as for a uniform flow scenario that is characterized by a power-law distribution of transition times [START_REF] Dentz | Time behavior of solute transport in heterogeneous media: transition from anomalous to normal transport[END_REF].

The radial geometry is reflected in the dependence of the characteristic time scales (42) on r. increasing β, the peak arrival time and the width of the breakthrough peak both decrease as indicated by (44). For times t τ k (r) we clearly observe the power-law behavior (43). As outlined above, the late time power-law behavior is the same as the one observed for uniform form. The scaling of peak arrival and width, however, are impacted on by both the radial geometry and the transition time PDF, as given in (42).

Multirate Mass Transfer

The CTRW model for multirate mobile-immobile mass exchange is characterized by the transition times [START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF]. The mobile times are given by τ m (r) = τ k (r)η, where η here is distributed according to the exponential PDF

ψ η (η) = exp(-η). ( 45
)
The distribution of immobile times p im (τ ) is given by the Pareto distribution

p im (τ ) = 1 τ c τ τ c -1-δ , τ ≥ τ c . ( 46
)
Notice that the memory function ϕ(t) for diffusive mass transfer between a mobile and homogeneous immobile zones behaves as ϕ(t) ∼ t -1/2 [START_REF] Haggerty | Multiple-rate mass transfer for modeling diffusion and surface reactions in media with pore-scale heterogeneity[END_REF][START_REF] Carrera | On matrix diffusion: formulations, solution methods, and qualitative effects[END_REF].

Thus, we obtain from ( 31) that the corresponding PDF of immobile times behaves as p im (τ ) ∼ t -3/2 , which corresponds to δ = 1/2 in (46). For heterogeneous immobile zones, one obtains in general different behaviors for the memory function ϕ(t) [START_REF] Gouze | Non-fickian dispersion in porous media explained by heterogeneous microscale matrix diffusion[END_REF] and therefore for the PDF of immobile times. The asymptotic behavior of the breakthrough curves for times t τ c is derived in Appendix B.2. It is given by

f (t, r) ∼ 1 θ c (r) t θ c (r) -1-δ , θ c (r) = τ c   γ r/ i=1 τ k (i )   1/δ . ( 47
)
Recall that γ is the rate for trapping in the immobile zones. The characteristic time θ c (r) is a measure for peak width and the time for the onset of the power-law tail behavior f (t, r) ∼ t -1-δ . Using the explicit τ k (r) = r/k v gives for θ c (r)

θ c (r) ≈ τ c γr 3 2k v 1/δ (48)
Figure 2 shows the behavior of the breakthrough curves for δ = 1/2 and varying dispersivity α and trapping rate γ. As for the case of heterogeneous advection, increase in dispersivity α leads to an increase of the width of the breakthrough peak, as illustrated in Figure 2a. The peak arrival time remains unchanged. Figure 2b shows solute breakthrough curves for varying trapping rates γ. As γ increases, the weight in the power-law tail increases as well as the peak width and therefore also the time of the onset of the power-law behavior. Note that increasing the trapping rate increases the proportion of trapped particles at any time, which explains the increased weight in the tail of the breakthrough curve. At the same time, increasing γ decreases the time that particles spend mobile, and thus particularly the average time until which particles get trapped for the first time. Thus, as particles notice the presence of immobile zones earlier, the peak width increases, and the breakthrough curve breaks off earlier from the behavior without immobile zones.

Heterogeneous Advection with Multirate Mass Transfer

The CTRW model combining heterogeneous advection in the mobile zone with mass exchange between mobile and immobile regions is based on the transition times [START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF]. The mobile time increment is again given by τ m (r) = δ-1 for τ k (r) τim . We consider the case 0 < δ < 1 such that this condition can only be achieved for γτ c 1, which implies a low trapping rate. Thus, in this time regime, most particles have not yet encountered a trapping event and thus, the breakthrough curve behaves as in the case of heterogeneous advection given by (43), f (t, r) ∼ t -1-β .

In the long time regime t τim , we need to distinguish between the cases 0 < β < 1 and 1 < β < 2. For 0 < β < 1, the long time breakthrough behavior is given by

f (t, r) ∼ t -1-βδ . (49)
Note that 0 < β < 1 indicates quite strong tailing in the mobile zone.

Thus both the exponents characteristic for the mobile and immobile zones determine the long time breakthrough behavior. For 1 < β < 2 the long time breakthrough behavior is dominated by particle retention due to trapping in the immobile regions. The breakthrough curve behaves asymptotically as given by (47).

These behaviors are illustrated in Figure 3. Figure 3a shows the behaviors in the intermediate and asymptotic long-time regimes for β = 3/2 and δ = 1/2. We see a relatively short intermediate regime (which increases as the trapping rate γ decreases) dominated by the advective heterogeneity and a the long time regime dominated by particle retention in the immobile regions.

A very similar behavior has been observed in the field from push pull tracer tests [START_REF] Le Borgne | Non-Fickian dispersion in porous media: 2. model validation from measurements at different scales[END_REF] pointing to a combined effect of advective and diffusive processes on non-Fickian transport. Figure 3b illustrates the breakthrough curve for β = 0.9 and δ = 0.75. The intermediate time regime again is dominated by the advective heterogeneity while the asymptotic tailing behavior is determined by both the advective and trapping exponents.

Summary and Conclusions

The CTRW approach provides a versatile framework for the modeling of non-Fickian solute transport in heterogeneous media. In this paper we develop a general CTRW approach for transport under radial flow conditions starting from the random walk equations for the radial and temporal particle coordinates. In contrast to CTRW formulations under uniform flow conditions here the random transition times form a non-stationary stochastic process, which reflects the radial dependence, i.e., non-stationarity, of the flow conditions. Thus, the evolution of the solute concentration is governed by a non-local ADE that is characterized by a radially dependent memory kernel. Within this general framework, we present CTRW models that may account for the impact of advective heterogeneity on large scale transport, that implement multirate mass transfer between a mobile immobile regions, and that combine non-local mobile transport due to heterogeneous advection and mass transfer into immobile zones. The three models are characterized by the specific forms of the stochastic process of particle times.

For heterogeneous advection, the transition time reflects the inverse (heterogeneous) flow velocity, which in average depends on the total flow rate and on the radial distance. This is accounted for by a transition time that is modeled as the product of a non-dimensional time increment and an advective transition time that is proportional to the radial distance and the inverse flow rate. Thus, transition times decrease with increasing flow rate as expected for purely advective heterogeneity. Flow heterogeneity is accounted for by the distribution of the dimensionless time increment [START_REF] Kang | Impact of velocity correlation and distribution on transport in fractured media: field evidence and theoretical model[END_REF]. The memory kernel of the non-local ADE for solute concentration depends here on the radial distance.

Radial transport under multirate mass transfer between mobile and immobile regions is implemented by separating the transition time into the time the particle is mobile and the sum of immobile times, which is represented by a compound Poisson process [START_REF] Margolin | Continuous time random walk and multirate mass transfer modeling of sorption[END_REF][START_REF] Benson | A simple and efficient random walk solution of multi-rate mobile/immobile mass transport equations[END_REF] whose mean is given by the trapping rate times the local mobile trapping time. The mobile time is modeled by an exponentially distributed dimensionless time increment, which models Fickian mobile transport. The mean number of trapping events during a mobile transition is given by the trapping rate and the mean mobile time. The memory kernel here is independent of the radial distance, which reflects the fact that the individual trapping events do not depend on mobile advection, which is consistent with radial MRMT formulation that assume homogeneous advection in the mobile zone [START_REF] Haggerty | Tracer tests in a fractured dolomite 2. analysis of mass transfer in single-well injectionwithdrawal tests[END_REF]. We note that ψ(λ, r i ) = 1 + ∆ ψ(λ, r i ), where ∆ ψ(λ, r i ) decreases with decreasing λ, see below. Thus, for small λ, we can approximate

p(λ, n) = exp n-1 i=0 ∆ψ(λ, r i ) , (B.8)
The asymptotic behavior will be determined starting from expression (B.6)

for the density of a single radial transition time.

Appendix B.1. Heterogeneous Advection Model

For the CTRW model defined through [START_REF] Dentz | Effective transport dynamics in porous media with heterogeneous retardation properties[END_REF], the single transition time τ i (r i ) is defined by

τ i (r i ) = τ k (r i )η i . (B.9)
The random variable η i is distributed according to a power-law such that In the regime (α im γτ c ) In the long time regime t τ c (α im γτ c ) 1 δ-1 , we need to distinguish between the cases 0 < β < 1 and 1 < β < 2. In the former case, the leading order of ψτ (λ, r i ) for λ (α im γτ c ) In the case 1 < β < 2, the leading order of ψ(λ, r i ) for λ (α im γτ c )

Figure 1 Figure 1 :

 11 Figure 1 shows solute breakthrough curves for different values of dispersivity α and exponents β of the Pareto distribution (39). The width of the peak of the breakthrough curves increases with increasing α, as shown in Figure 1a. The peak arrival time remains unaffected by α.Figure 1 shows the behavior of the breakthrough curves for different exponents β. With

Figure 2 : 1 and γ = 10 - 2 , 5 × 10 - 2 , 2 . 5 × 10 - 1 .

 21102510225101 Figure 2: Breakthrough curves for the MRMT model (21) with (a) α = 10 -2 , 4 × 10 -2 , 10 -1 , δ = 1/2, τ c = 10 -1 and γ = 10 -1 , and (b) α = 2 × 10 -2 , δ = 1/2, τ c = 10 -1and γ = 10 -2 , 5 × 10 -2 , 2.5 × 10 -1 . The green lines show the power-law ∼ t -1-β behavior for β = 1/2.

τ

  k (r)η, where η now is modeled by the Pareto distribution (39) characterized by the exponent β. The PDF of immobile times is again given by the Pareto distribution (46) characterized by the exponent δ. We focus here on situations, for which β > δ. Notice that the exponents β and δ encode the width of the distributions of characteristic advection scales and trapping time scales, respectively. The width of characteristic retention time scales is typically smaller than the one of characteristic advection time scales in the mobile zone, and thus we set β > δ. Appendix B.3 develops the asymptotic breakthrough behavior for this model. One can distinguish two time regimes with distinct temporal behavior. An intermediate time regime is set by the advection scale τ k (r) and the immobile time scale τim = τ c (α im γτ c ) 1

Figure 3 :

 3 Figure 3: Breakthrough curves for the combined model (21) with (a) α = 2×10 -2 , β = 3/2, δ = 1/2, τ c = 1 and (from top to bottom) γ = 10 -2 , 10 -3 , 10 -4 , (b) α = 2 × 10 -2 , β = 0.9, δ = 3/4, τ c = 1 and γ = 5 × 10 -2 . The green lines indicates the power-law ∼ t -1-β for (a) β = 3/2 and (b) β = 0.9. The blue line in (a) indicates the power-law ∼ t -1-δ for δ = 1/2. The blue line in (b) indicates the power-law ∼ t -1-βδ for β = 0.9 and δ = 3/4.

Finally, we propose

  a CTRW model that combines heterogeneous advection in the mobile continuum with multirate mass transfer into immobile continua. The mobile transition time now is modeled by a non-exponential dimensionless time increment. The immobile time increment is again a compound Poisson process. While the individual trapping events are independent on the heterogeneous advection, the collective of trapping times depends on the radial distance through the mean number of trapping events. Unlike for the case of a homogeneous mobile advection, heterogeneity leads to an interrelation of the PDFs of mobile transition times and retention times in immobile regions, and thus to a radial dependence of the memory kernel.The derivation of the general framework from the random walk equations for the particle coordinates in space and time provides directly the random walk particle tracking methods for the numerical solution of radial MRMT models, and temporally non-local radial transport formulations in general.Notice that the developed non-local transport models assume that dispersion is proportional to the flow velocity, as well as a constant injection or withdrawal rate of the flow at the well. It is straightforward to relax the first assumption and introduce a more general velocity dependence of dispersion in the developed models. The assumption of a constant flow rate may be relaxed to account for scenarios encountered in single well injection withdrawal tests, i.e., piecewise constant flow rates.The solute breakthrough behaviors for the different non-local radial transport models are studied in detail by using random walk particle tracking simulations and explicit analytical expressions for the asymptotic time behaviors of the solute breakthrough curves. A broad distribution of advective transition times is modeled by a Pareto distributions of dimensionless times.The distribution of immobile times is also given by a Pareto distribution, which mimics a broad distribution of diffusion or in general retention times in the immobile zones. We find distinct power-law tail behavior in all cases, which are similar to the ones that can be observed for uniform flow. However, the non-stationarity of the random time increment leads to a non-trivial radial dependence of the breakthrough curves. The model combining heterogeneous advection and MRMT can display an intermediate time regime in which heterogeneous advection dominates, before the onset of trapping, and a late time regime that is governed by the retention properties in the im-mobile regions. The different behaviors in the two time-regimes allow for the discrimination and identification of heterogeneous advection and mobileimmobile mass transfer as drivers of anomalous transport. The CTRW model accounting for advective heterogeneity is controlled by the flow rate, while the CTRW model for MRMT is controlled by the trapping rate and the properties of the immobile zone. In conclusion, this paper provides a radial CTRW framework that allows for the systematic interpretation of tracer data from radial forced gradient test, through the implementation of different mechanisms that lead to nonlocal radial transport. where ψ(λ, r i ) is given by ψ(λ, r i ) = exp [-λτ i (r i )] . (B.6) Equation (B.5) can also be written as p(λ, n) = exp n-1 i=0 ln[ ψ(λ, r i )] , (B.7)

ψ(α im γτ c ) 1 δ- 1 . (B. 35 ) 1 δ- 1 1 δ- 1 .

 11351111 η (η) ∼ η -1-β , (B.[START_REF] Dentz | Time behavior of solute transport in heterogeneous media: transition from anomalous to normal transport[END_REF] for large η and 0 < β < 2. Thus, (B.6) can be written asψ(λ, r i ) = exp [-λτ k (r i )η i ] = ψη [λτ k (r i )]. (B.11)Notice that the Laplace transform of the power-law density (B.10) can be expanded for small λ as [e.g., 10]ψη (λ) = 1 -α 11 λ β , 0 < β < 1 (B.12) ψη (λ) = 1 -α 12 λ + α 22 λ β , 1 < β < 2, (B.13)where the parameters α 11 , α 12 and α 22 depend on the details of the underlying distribution ψ η (η). For the Pareto distribution[START_REF] Gouze | Non-fickian dispersion in porous media explained by heterogeneous microscale matrix diffusion[END_REF], they are given byα 11 = Γ(1 -β), (B.14) (λ, r i ) can be approximated for λτ k (r i ) 1 as ψ(λ, r i ) = 1 -α 11 [λτ k (r i )] β , 0 < β < 1 (B.16a) ψ(λ, r i ) = 1 -α 12 λτ k (r i ) + α 22 [λτ k (r i )] β , 1 < β < 2. (B.16b)Inserting the latter into (B.8) givesp(λ, n) = exp -α 11 (λθ n ) β , 0 < β < 1 (B.17) p(λ, n) = exp -λ τ a,n + α 22 (λθ n ) β ,the mean arrival time, which is only defined for 1 < β < 2. It is impacted on by the distribution of dimensionless transition times through α 11 . For α 11 = 1, it is equal to the mean arrival time of the homogeneous Using now expansion (B.16) for ψτ (λ, r) gives ψτ (λ, ri ) = 1 -α 1 λ + γα im (λτ c ) δ τ k (r i ) β , 0 < β < 1 (B.33a) ψτ (λ, r i ) = 1 -α 1 λ + γα im (λτ c ) δ τ k (r i ) + α 2 λ + γα im (λτ c ) δ τ k (r i ) β , 1 < β < 2. (B.33b) From these expressions, we can identify two characteristic times scales that separate time regimes with different behaviors. The first one is set by the condition λ (α im γτ c ) The second time scale is set by the condition that λτ k (r i ) 1, which implies t τ k (r i ). (B.36) These scales set the time regimes τ k (r i ) t τ c (α im γτ c ) , and t τ c (α im γτ c )

1 1-δ τ - 1 c 1 δ- 1

 1111 λ τ k (r) -1 , the transition time density ψ(λ, r i ) is given in leading order by ψτ (λ, r i ) = 1 -α 1 [λτ k (r i )] β , 0 < β < 1 (B.37a) ψτ (λ, r i ) = 1 -α 1 λτ k (r i ) + α 2 [λτ k (r i )] β , 1 < β < 2, (B.37b)which is identical to (B.16). Thus, the breakthrough curves for τ k (r) t τ c (α im γτ c ) behave as in (B.26).

1 1-δ τ - 1 c

 11 is given by ψτ (λ,r i ) = 1 -α 1 α β im [γτ k (r i )] β (λτ c ) βδ . (B.38)Inserting this expression into (B.5) and expanding up to leading order gives for p(λ, n)p(λ, n) = 1 -α 1 α β im (λτ c ) βδ n i=1 [γτ k (r i )] β . (B.39)Thus, p(t, n) behaves asymptotically as p(t, n) ∼ t -1-βδ . (B.40)

1 1-δ τ - 1 c

 11 is given by ψ(λ, r i ) = 1 -α 1 α im (λτ c ) δ γτ k (r i ).(B.41)Inserting this expression into (B.5) and expanding up to leading order gives for p(λ, n) similar as in the previous sectionp(λ, n) = 1 -α 1 α im (λτ c ) δ γ n i=1 τ k (r i ). (B.42)Thus, the asymptotic behavior of p(t, n) is given by (B.31).

Acknowledgements: MD acknowledges the support of the European Research Council (ERC) through the project MHetScale (contract number 617511).

We thank Dr. Antoine Aubenau and two anonymous reviewers for their comments.

Appendix A. Derivation of the Transition Time Distribution for the Compound Time Process

We derive here the probability density function (PDF) ψ(t, r) for the general compound process [START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF] of the time increment τ (r). The density of the mobile time τ m (r) = τ k (r)η is denoted by ψ m (τ, r), and can be expressed in terms of the PDF ψ η (η) as

The density of the compound process τ (r) can be written in general as

where the angular brackets denote the ensemble average over all particles.

Inserting [START_REF] Dentz | Transport behavior of a passive solute in continuous time random walks and multirate mass transfer[END_REF] and performing the average over the τ m (r) gives

Performing the Laplace transform in τ and performing the average of the

Performing the average of the Poisson variable n im results in

The infinite series can be summed up to an exponential such that ψτ (λ, r)

The integration over τ m can be performed explicitly by noting that the integral is equal to the Laplace transform of ψ m (τ m , r). Thus we obtain directly [START_REF] Silva | A general real-time formulation for multi-rate mass transfer problems[END_REF]. Using the Laplace transform of the exponential distribution

gives directly [START_REF] Margolin | Continuous time random walk and multirate mass transfer modeling of sorption[END_REF].

Appendix B. Asymptotic Behavior

We consider now the asymptotic long time behavior of the breakthrough curves for the CTRW models presented in Sections 3.1, 3.2 and 3.3. The distribution f (t, r) of arrival times at a radius r is given by

where f 0 (n, r) is the distribution of the number of steps needed to arrive at the radius r in the random walk (9a), and p(t, n) the distribution of times needed to make n steps in the random walk (9b). For large n it can be approximated by an inverse Gaussian distribution that is sharply peaked about n = r/ . Thus, we approximate the breakthrough curve by

In the following, we derive the asymptotic behavior of p(t, n) for large times t. Notice that p(t, n) can be defined by

where t n is given by (9b) and can be written as

Thus, the Laplace transform of p(t, n) can be written as

model. From (B.17) and (B.18), we can deduce the scaling forms

where the Laplace transforms of f 01 (t) and f 02 (t) are defined by

We obtain the long-time behavior of the breakthrough curves by expanding the latter expressions up to leading order in λ 1, which gives

Thus, f 01 (t) and f 02 (t) behave at long times as ∼ t -1-β and thus p(t, n) behaves for t θ n as

Expression ( 43) is obtained by setting n = r/ , r i = i , and θ(r) = θ r/ .

Appendix B.2. MRMT Model

We employ for the trapping time distribution p im (τ ) the power-law

for t τ c , τ c a characteristic immobile time and 0 < δ < 1. Thus, its Laplace transform for λτ c 1 can be written as

Inserting the latter into [START_REF] Margolin | Continuous time random walk and multirate mass transfer modeling of sorption[END_REF] for the transition time distribution and expanding the resulting expression up to the leading order, we obtain ψτ (λ,

Using this expression in (B.5) and again expanding to leading order yields for p(λ, n)

where α im is given by the details of p im (τ ). Thus, the asymptotic behavior of p(t, n) is given by

Expression (47) is obtained by setting r i = i , n = r and θ c,r = θ c (r).

Appendix B.3. Heterogeneous Advection with MRMT

We consider now the case that the mobile transition time is given by (B.9)


K A Ïs

Mourad Choulli
email: mourad.choulli@univ-lorraine.fr

LOGARITHMIC STABILITY IN DETERMINING A BOUNDARY COEFFICIENT IN AN IBVP FOR THE WAVE EQUATION

Keywords: inverse problem, wave equation, boundary damping coefficient, logarithmic stability, boundary measurements. MSC: 35R30

we introduced a method combining together an observability inequality and a spectral decomposition to get a logarithmic stability estimate for the inverse problem of determining both the potential and the damping coefficient in a dissipative wave equation from boundary measurements. The present work deals with an adaptation of that method to obtain a logarithmic stability estimate for the inverse problem of determining a boundary damping coefficient from boundary measurements. As in our preceding work, the different boundary measurements are generated by varying one of the initial conditions.

Introduction

We are concerned with an inverse problem for the wave equation when the spatial domain is the square Ω = (0, 1) × (0, 1). To this end we consider the following initial-boundary value problem (abbreviated to IBVP in the sequel) :

(1.1)

       ∂ 2 t u -∆u = 0 in Q = Ω × (0, τ ), u = 0 on Σ 0 = Γ 0 × (0, τ ), ∂ ν u + a∂ t u = 0 on Σ 1 = Γ 1 × (0, τ ), u(•, 0) = u 0 , ∂ t u(•, 0) = u 1 .
Here Γ 0 = ((0, 1) × {1}) ∪ ({1} × (0, 1)), Γ 1 = ((0, 1) × {0}) ∪ ({0} × (0, 1))

1
and ∂ ν = ν • ∇ is the derivative along ν, the unit normal vector pointing outward of Ω. We note that ν is everywhere defined except at the vertices of Ω and we denote by Γ = Γ 0 ∪ Γ 1 . The boundary coefficient a is usually called the boundary damping coefficient.

In the rest of this text we identify a |(0,1)×{0} by a 1 = a 1 (x), x ∈ (0, 1) and a |{0}×(0,1) by a 2 = a 2 (y), y ∈ (0, 1). In that case it is natural to identify a, defined on Γ 1 , by the pair (a 1 , a 2 ).

1.1. The IBVP. We fix 1/2 < α ≤ 1 and we assume that a ∈ A , where

A = {b = (b 1 , b 2 ) ∈ C α ([0, 1]) 2 , b 1 (0) = b 2 (0), b j ≥ 0}.
This assumption guarantees that the multiplication operator by a j , j = 1, 2, defines a bounded operator on H 1/2 ((0, 1)). The proof of this fact will be proved in Appendix A.

Let V = {u ∈ H 1 (Ω); u = 0 on Γ 0 } and we consider on V × L 2 (Ω) the linear unbounded operator A given by A a = (w, ∆v), D(A a ) = {(v, w) ∈ V × V ; ∆v ∈ L 2 (Ω) and ∂ ν v = -aw on Γ 1 }. One can prove that A a is a m-dissipative operator on the Hilbert space V × L 2 (Ω) (for the reader's convenience we detail the proof in Appendix B). Therefore, A a is the generator of a strongly continuous semigroup of contractions e tAa . Hence, for each (u 0 , u 1 ), the IBVP (1.1) possesses a unique solution denoted by u a = u a (u 0 , u 1 ) so that

(u a , ∂ t u a ) ∈ C([0, ∞); D(A a )) ∩ C 1 ([0, ∞), V × L 2 (Ω)). 1.2. Main result. For 0 < m ≤ M , we set A m,M = {b = (b 1 , b 2 ) ∈ A ∩ H 1 (0, 1) 2 ; m ≤ b j , b j 2 H 1 (0,1) ≤ M }. Let U 0 given by U 0 = {v ∈ V ; ∆v ∈ L 2 (Ω) and ∂ ν v = 0 on Γ 1 }. We observe that U 0 × {0} ⊂ D(A a ), for any a ∈ A . Let C a ∈ B(D(A a ); L 2 (Σ 1 )) defined by C a (u 0 , u 1 ) = ∂ ν u a (u 0 , u 1 ) |Γ1 .
We define the initial to boundary operator

Λ a : u 0 ∈ U 0 -→ C a (u 0 , 0) ∈ L 2 (Σ 1 ).
Clearly C a ∈ B(D(A a ); L 2 (Σ 1 )) implies that Λ a ∈ B(U 0 ; L 2 (Σ 1 )), when U 0 is identified to a subspace of D(A a ) endowed with the graph norm of A a . Precisely the norm in U 0 is the following one

u 0 U0 = u 0 2 V + ∆u 0 2 L 2 (Ω) 1/2 .
Henceforth, for simplicity sake, the norm of Λ a -Λ 0 in B(U 0 ; L 2 (Σ 1 )) will denoted by Λ a -Λ 0 .

Theorem 1.1. There exists τ 0 > 0 so that for any τ > τ 0 , we find a constant c > 0 depending only on τ such that

(1.2) a -0 L 2 ((0,1)) 2 ≤ cM ln m -1 Λ a -Λ 0 -1/2 + m -1 Λ a -Λ 0 L 2 (Σ1) ,
for each a ∈ A m,M .

We point out that our choice of the domain Ω is motivated by the fact the spectral analysis of the laplacian under mixed boundary condition is very simple in that case. However this choice has the inconvenient that the square domain Ω is no longer smooth. So we need to prove an observability inequality associated to this non smooth domain. This is done by adapting the existing results. We note that the key point in establishing this observability inequality relies on a Rellich type identity for the domain Ω.

The inverse problem we discuss in the present paper remains largely open for an arbitrary (smooth) domain as well as for the stability around a non zero damping coefficient. Uniqueness and directional Lipschitz stability, around the origin, was established by the authors in [START_REF] Ammari | Determining a boundary coefficient in a dissipative wave equation: uniqueness and directional Lipschitz stability[END_REF].

The determination of a potential and/or the sound speed coefficient in a wave equation from the so-called Dirichlet-to-Neumann map was extensively studied these last decades. We refer to the comments in [START_REF] Ammari | Logarithmic stability in determining two coefficients in a dissipative wave equation[END_REF] for more details.

Preliminaries

2.1. Extension lemma. We decompose Γ 1 as follows Γ 1 = Γ 1,1 ∪ Γ 1,2 , where Γ 1,1 = (0, 1) × {0} and Γ 1,2 = {0} × (0, 1). Similarly, we write Γ 0 = Γ 0,1 ∪ Γ 0,2 , with Γ 0,1 = {1} × (0, 1) and Γ 0,2 = (0, 1) × {1}.

Let (g 1 , g 2 ) ∈ L 2 ((0, 1)) 2 . We say that the pair (g 1 , g 2 ) satisfies the compatibility condition of the first order at the vertex (0, 0) if (2.1)

1 0 |g 1 (t) -g 2 (t)| 2 dt t < ∞.
Similarly, we can define the compatibility condition of the first order at the other vertices of Ω.

We need also to introduce compatibility conditions of the second order. Let (f j , g j ) ∈ H 1 ((0, 1)) × L 2 ((0, 1)), j = 1, 2. We say that the pair [(f 1 , g 1 ), (f 2 , g 2 )] satisfies the compatibility conditions of second order at the vertex (0, 0) when

(2.2) f 1 (0) = f 2 (0), 1 0 |f ′ 1 (t) -g 2 (t)| 2 dt t < ∞ and 1 0 |g 1 (t) -f ′ 2 (t)| 2 dt t < ∞.
The compatibility conditions of the second order at the other vertices of Ω are defined in the same manner.

The following theorem is a special case of [4, Theorem 1.5.2.8, page 50].

Theorem 2.1.

(1) The mapping

w -→ (w |Γ0,1 , w |Γ0,2 , w |Γ1,1 , w |Γ1,2 ) = (g 1 , . . . , g 4 ),
defined on D(Ω) is extended from H 1 (Ω) onto the subspace of H 1/2 ((0, 1)) 4 consisting in functions (g 1 , . . . , g 4 ) so that the compatibility condition of the first order is satisfied at each vertex of Ω in a natural way with the pairs (g j , g k ).

(2) The mapping

w → (w |Γ0,1 , ∂ x w |Γ0,1 , w |Γ0,2 , ∂ y w |Γ0,2 w |Γ1,1 , -∂ y w |Γ1,1 , w |Γ1,2 , -∂ x w |Γ1,2 ) = ((f 1 , g 1 ), . . . (f 4 , g 4 ))
defined on D(Ω) is extended from H 2 (Ω) onto the subspace of [H 3/2 ((0, 1)) × H 1/2 ((0, 1))] 4 of functions ((f 1 , g 1 ), . . . (f 4 , g 4 )) so that the compatibility conditions of the second order are satisfied at each vertex of Ω in a natural way with the pairs [(f j , g j ), (f k , g k )].

Lemma 2.1. (Extension lemma) Let g j ∈ H 1/2 ((0, 1)), j = 1, 2, so that (g 1 , g 2 ), (g 1 , 0) and (g 2 , 0) satisfy the first order compatibility condition respectively at the vertices (0, 0), (1, 0) and (0, 1). Then there exists u ∈ H 2 (Ω) so that u = 0 on Γ 0 and ∂ ν u = g j on Γ 1,j , j = 1, 2.

Proof. (i) We define f 1 (t) = t 0 g 2 (s)ds and f 2 (t) = t 0 g 1 (s)ds. Then (f 1 , g 1 ) and (f 2 , g 2 ) satisfy the compatibility conditions of the second order at the vertex (0, 0).

(ii) Let g 1 ∈ H 1/2 ((0, 1)) be such that 1 0 | g1(t)| 2 t dt < ∞. Let f 1 (t) = t 0 g 2 (s)ds.
Hence, it is straightforward to check that ( f 1 , g 1 ) and (0, g 2 ) satisfy the compatibility conditions of the second order at (0, 0).

(iii) From steps (i) and (ii) we derive that the pairs [(f 1 , g 1 ), (f 2 , g 2 )], [(f 1 , g 1 ), (0, g 2 )] and [(0, g 1 ), (f 2 , g 2 )] satisfy the second order compatibility conditions respectively at the vertices (0, 0), (1, 0) and (0, 1). We see that unfortunately the pair [(0, g 1 ), (0, g 2 )] doesn't satisfy necessarily the compatibility conditions of the second order at the vertex (1, 1). We pick χ ∈ C ∞ (R) so that χ = 1 in a neighborhood of 0 and χ = 0 in a neighborhood of 1. Then [(0, χg 1 ), (0, χg 2 )] satisfies the compatibility condition of the second order at the vertex [START_REF] Alves | Solving inverse source problems using observability. Applications to the Euler-Bernoulli plate equation[END_REF][START_REF] Alves | Solving inverse source problems using observability. Applications to the Euler-Bernoulli plate equation[END_REF]. Since this construction is of local character at each vertex, the cutoff function at the vertex (1, 1) doesn't modify the construction at the other vertices. In other words, the compatibility conditions of the second order are preserved at the other vertices. We complete the proof by applying Theorem 2.1.

Corollary 2.1. Let a = (a 1 , a 2 ) ∈ A and g j ∈ H 1/2 ((0, 1)), j = 1, 2, so that (g 1 , g 2 ), (g 1 , 0) and (g 2 , 0) satisfy the first order compatibility condition respectively at the vertices (0, 0), (1, 0) and (0, 1). Then there exists u ∈ H 2 (Ω) so that u = 0 on Γ 0 and ∂ ν u = a j g j on Γ 1,j , j = 1, 2.

Proof. It is sufficient to prove that (a 1 g 1 , a 2 , g 2 ) and (a j g j , 0), j = 1, 2, satisfy the first order compatibility condition at (0, 0) with a 1 (0) = a 2 (0) for the first pair and without any condition on a j for the second pair.

Using a 1 (0) = a 2 (0), we get

t -1 |a 1 (t) -a 2 (t)| 2 ≤ 2t -1 |a 1 (t) -a 1 (0)| 2 + 2t -1 |a 2 (t) -a 2 (0)| 2 ≤ 2t -1+2α ([a 1 ] 2 α + [a 2 ] 2 α ) ≤ 2([a 1 ] 2 α + [a 2 ] 2 α
). This estimate together with the following one

|a 1 (t)g 1 (t) -a 2 (t)g 2 (t)| 2 ≤ 2|a 1 (t) -a 2 (t)| 2 |g 1 (t)| 2 + 2|a 2 (t)| 2 |g 1 (t) -g 2 (t)| 2 yield 1 0 |a 1 (t)g 1 (t) -a 2 (t)g 2 (t)| 2 dt t ≤ 4([a 1 ] 2 α + [a 2 ] 2 α ) f L 2 ((0,1)) + 2 a 2 L ∞ ((0,1)) 1 0 |g 1 (t) -g 2 (t)| 2 dt t .
Hence

1 0 |g 1 (t) -g 2 (t)| 2 dt t < ∞ =⇒ 1 0 |a 1 (t)g 1 (t) -a 2 (t)g 2 (t)| 2 dt t < ∞.
If (g j , 0) satisfies the first compatibility at the vertex (0, 0). Then

1 0 |g j (t)| 2 dt t < ∞. Therefore 1 0 |a j g j (t)| 2 dt t ≤ a j 2 L ∞ ((0,1)) 1 0 |g j (t)| 2 dt t < ∞.
Thus (a j g j , 0) satisfies also the first compatibility at the vertex (0, 0).

Observability inequality.

We discuss briefly how we can adapt the existing results to get an observability inequality corresponding to our IBVP. We first note that

Γ 0 ⊂ {x ∈ Γ; m(x) • ν(x) < 0}, Γ 1 ⊂ {x ∈ Γ; m(x) • ν(x) > 0},
where m(x) = x -x 0 , x ∈ R 2 , and x 0 = (α, α) with α > 1.

The following Rellich identity is a particular case of identity [5, (3.5), page 227]: for each 3/2 < s < 2 and ϕ ∈ H s (Ω) satisfying ∆ϕ ∈ L 2 (Ω),

(2.3) 2 Ω ∆ϕ(m • ∇ϕ)dx = 2 Γ ∂ ν ϕ(m • ∇ϕ)dσ - Γ (m • ν)|∇ϕ| 2 dσ. Lemma 2.2. Let (v, w) ∈ D(A a ). Then 2 Ω ∆v(m • ∇v)dx = 2 Γ ∂ ν v(m • ∇v)dσ - Γ (m • ν)|∇v| 2 dσ.
Proof. Let (v, w) ∈ D(A a ). By Corollary 2.1, there exists v ∈ H 2 (Ω) so that v = 0 on Γ 0 and ∂ ν v = -aw on Γ 1 . In light of the fact that z = v -v is such that ∆z ∈ L 2 (Ω), z = 0 on Γ 0 and ∂ ν z = 0 on Γ 1 , we get z ∈ H s (Ω) for some 3/2 < s < 2 by [5, Theorem 5.2, page 237]. Therefore v ∈ H s (Ω). We complete the proof by applying Rellich identity (2.3).

Lemma 2.2 at hand, we can mimic the proof of [7, Theorem 7.6.1, page 252] in order to obtain the following theorem:

Theorem 2.2. We assume that a ≥ δ on Γ 1 , for some δ > 0. There exist M ≥ 1 and ω > 0, depending only on δ, so that

e tAa (v, w) V ×L 2 (Ω) ≤ M e -ωt (v, w) V ×L 2 (Ω) , (v, w) ∈ D(A a ), t ≥ 0.
An immediate consequence of Theorem 2.2 is the following observability inequality.

Corollary 2.2. We fix 0 < δ 0 < δ 1 . Then there exist τ 0 > 0 and κ, depending only on δ 0 and δ 1 so that for any τ ≥ τ 0 and a ∈ A satisfying δ 0 ≤ a ≤ δ 1 on Γ 1 ,

(u 0 , u 1 ) V ×L 2 (Ω) ≤ κ C a (u 0 , u 1 ) L 2 (Σ1) .
Moreover, C a is admissible for e tAa and (C a , A a ) is exactly observable.

We omit the proof of this corollary. It is quite similar to that of [7, Corollary 7.6.5, page 256].

The inverse problem

3.1. An abstract framework for the inverse source problem. In the present subsection we consider an inverse source problem for an abstract evolution equation. The result of this subsection is the main ingredient in the proof of Theorem 1.1.

Let H be a Hilbert space and A : D(A) ⊂ H → H be the generator of continuous semigroup (T (t)). An operator C ∈ B(D(A), Y ), Y is a Hilbert space which is identified with its dual space, is called an admissible observation for (T (t)) if for some (and hence for all) τ > 0, the operator Ψ ∈ B(D(A), L 2 ((0, τ ), Y )) given by (Ψx

)(t) = CT (t)x, t ∈ [0, τ ], x ∈ D(A), has a bounded extension to H.
We introduce the notion of exact observability for the system

z ′ (t) = Az(t), z(0) = x, (3.1)
y(t) = Cz(t), (3.2)
where C is an admissible observation for T (t). Following the usual definition, the pair (A, C) is said exactly observable at time τ > 0 if there is a constant κ such that the solution (z, y) of (3.1) and (3.2) satisfies

τ 0 y(t) 2 Y dt ≥ κ 2 x 2 H , x ∈ D(A).
Or equivalently

(3.3) τ 0 (Ψx)(t) 2 Y dt ≥ κ 2 x 2 H , x ∈ D(A).
Let λ ∈ H 1 ((0, τ )) such that λ(0) = 0. We consider the Cauchy problem (3.4) z ′ (t) = Az(t) + λ(t)x, z(0) = 0 and we set

(3.5) y(t) = Cz(t), t ∈ [0, τ ].
We fix β in the resolvent set of A. Let H 1 be the space D(A) equipped with the norm x 1 = (β -A)x and denote by H -1 the completion of H with respect to the norm x -1 = (β -A) -1 x . As it is observed in [1, Proposition 4.2, page 1644] and its proof, when x ∈ H -1 (which is the dual space of H 1 with respect to the pivot space H) and λ ∈ H 1 ((0, T )), then, according to the classical extrapolation theory of semigroups, the Cauchy problem (3.4) has a unique solution z ∈ C([0, τ ]; H). Additionally y given in (3.5) belongs to L 2 ((0, τ ), Y ).

When x ∈ H, we have by Duhamel's formula

(3.6) y(t) = t 0 λ(t -s)CT (s)xds = t 0 λ(t -s)(Ψx)(s)ds. Let H 1 ℓ ((0, τ ), Y ) = u ∈ H 1 ((0, τ ), Y ); u(0) = 0 . We define the operator S : L 2 ((0, τ ), Y ) -→ H 1 ℓ ((0, τ ), Y ) by (3.7) (Sh)(t) = t 0 λ(t -s)h(s)ds.
If E = SΨ, then (3.6) takes the form y(t) = (Ex)(t).

Let Z = (β -A * ) -1 (X + C * Y ).
Theorem 3.1. We assume that (A, C) is exactly observable at time τ . Then (i) E is one-to-one from H onto H 1 ℓ ((0, τ ), Y ). (ii) E can be extended to an isomorphism, denoted by E, from Z ′ onto L 2 ((0, τ ); Y ). (iii) There exists a constant κ, independent on λ, so that

(3.8) x Z ′ ≤ κ|λ(0)|e λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 τ Ex L 2 ((0,τ ),Y ) .
Proof. (i) and (ii) are contained in [1, Theorem 4.3, page 1645]. We need only to prove (iii). To do this, we start by observing that

S * : L 2 ((0, τ ), Y ) → H 1 r ((0, τ ); Y ) = u ∈ H 1 ((0, τ ), Y ); u(τ ) = 0 , the adjoint of S, is given by S * h(t) = τ t λ(s -t)h(s)ds, h ∈ H 1 r ((0, τ ); Y ).
We fix h ∈ H 1 r ((0, τ ); Y ) and we set k = S * h. Then

k ′ (t) = λ(0)h(t) - τ t λ ′ (s -t)h(s)ds. Hence [|λ(0)| h(t) ] 2 ≤ τ t |λ ′ (s -t)| |λ(0)| [|λ(0)| h(s) ]ds + k ′ (t) 2 ≤ 2 τ t |λ ′ (s -t)| |λ(0)| [|λ(0)| h(s) ]ds 2 + 2 k ′ (t) 2 ≤ 2 λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 t 0 [|λ(0)| h(s) ] 2 ds + 2 k ′ (t) 2 .
The last estimate is obtained by applying Cauchy-Schwarz's inequality.

A simple application of Gronwall's lemma entails

[|λ(0)| h(t) ] 2 ≤ 2e 2 λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 τ k ′ (t) 2 .
Therefore,

h L 2 ((0,τ );Y ) ≤ √ 2 |λ(0)| e λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 τ k ′ L 2 ((0,τ );Y ) .
This inequality yields

(3.9) h L 2 ((0,τ );Y ) ≤ √ 2 |λ(0)| e λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 τ S * h H 1 r ((0,τ );Y ) .
The adjoint of S * , acting as a bounded operator from [H r ((0, 1); Y )] ′ into L 2 ((0, τ ); Y ), gives an extension of S. We denote by S this operator. By [1, Proposition 4.1, page 1644] S * defines an isomorphism from [H r ((0, 1); Y )] ′ onto L 2 ((0, τ ); Y ). In light of the fact that S B([Hr((0,1);Y )] ′ ;L 2 ((0,τ );Y )) = S * B(L 2 ((0,τ );Y );Hr((0,1);Y )) , (3.9) implies

(3.10) |λ(0)| √ 2 e - λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 τ ≤ S B([Hr((0,1);Y )] ′ ;L 2 ((0,τ );Y )) .
On the other hand, according to [1, Proposition 2.13, page 1641], Ψ possesses a unique bounded extension, denoted by Ψ from Z ′ into [H r ((0, 1); Y )] ′ and there exists a constant c > 0 so that

(3.11) Ψ B(Z ′ ;[Hr((0,1);Y )] ′ ) ≥ c.
Consequently, E = S Ψ gives a unique extension of E to an isomorphism from Z ′ onto L 2 ((0, τ ); Y ). We end up the proof by noting that (3.8) is a consequence of (3.9) and (3.11).

3.2.

An inverse source problem for an IBVP for the wave equation. In the present subsection we are going to apply the result of the preceding subsection to H = V × L 2 (Ω), H 1 = D(A a ) equipped with its graph norm and Y = L 2 (Γ 1 ). We consider the the IBVP (3.12)

       ∂ 2 t u -∆u = λ(t)w in Q, u = 0 on Σ 0 , ∂ ν u + a∂ t u = 0 on Σ 1 , u(•, 0) = 0, ∂ t u(•, 0) = 0.
Let (0, w) ∈ H -1 and λ ∈ H 1 ((0, τ )). From the comments in the preceding subsection, (3.12) has a unique solution u w so that (u w , ∂

t u w ) ∈ C([0, τ ]; V × L 2 (Ω)) and ∂ ν u w|Γ 1 ∈ L 2 (Σ 1 ).
We consider the inverse problem consisting in the determination of w, so that (0, w) ∈ H -1 , appearing in the IBVP (3.12) from the boundary measurement ∂ ν u w |Σ1 . Here the function λ is assumed to be known.

Taking into account that {0} × V ′ ⊂ H -1 , where V ′ is the dual space of V , we obtain as a consequence of Corollary 2.1: Proposition 3.1. There exists a constant C > 0 so that for any λ ∈ H 1 ((0, τ )) and w ∈ V ′ ,

(3.13) w V ′ ≤ C|λ(0)|e λ ′ 2 L 2 ((0,τ )) |λ(0)| 2 τ ∂ ν u L 2 (Σ1) .
3.3. Proof of Theorem 1.1. We start by observing that u a is also the unique solution of

Ω u ′′ (t)vdx = Ω ∇u(t) • ∇vdx -Γ1 au ′ (t)v, for all v ∈ V. u(0) = u 0 , u ′ (0) = u 1 .
Let u = u a -u 0 . Then u is the solution of the following problem

(3.14) Ω u ′′ (t)vdx = Ω ∇u(t) • ∇vdx -Γ1 au ′ (t)v -Γ1 au ′ 0 (t)v, for all v ∈ V. u(0) = 0, u ′ (0) = 0.
For k, ℓ ∈ Z, we set

λ kℓ = [(k + 1/2) 2 + (ℓ + 1/2) 2 ]π 2 φ kℓ (x, y) = 2 cos((k + 1/2)πx) cos((ℓ + 1/2)πy).
We check in a straightforward manner that u 0 = cos( √ λ kℓ t)φ kℓ when (u 0 , u 1 ) = (φ kℓ , 0). In the sequel k, ℓ are arbitrarily fixed. We set λ(t) = cos( √ λ kℓ t) and we define w a ∈ V ′ by

w a (v) = -λ kℓ Γ1 aφ kℓ v.
In that case (3.14) becomes

Ω u ′′ (t)vdx = Ω ∇u(t) • ∇vdx -Γ1 au ′ (t)v + λ(t)w a (v), for all v ∈ V. u(0) = 0, u ′ (0) = 0.
Consequently, u is the solution of (3.12) with w = w a . Applying Proposition 3.1, we find

(3.15) w a V ′ ≤ Ce λ kℓ τ 2 ∂ ν u L 2 (Σ1) . But (3.16) a 1 (0) Γ1 (aφ kℓ ) 2 dσ = 1 √ λ kℓ |w a ((a 1 ⊗ a 2 )φ kℓ )| ≤ 1 √ λ kℓ w a V ′ (a 1 ⊗ a 2 )φ kℓ V ,
where we used a 1 (0) = a 2 (0), and

(3.17) (a 1 ⊗ a 2 )φ kℓ V ≤ C 0 λ kl a 1 ⊗ a 2 H 1 (Ω) .
Here C 0 is a constant independent on a and φ kℓ .

We note (a

1 ⊗ a 2 )φ kℓ ∈ V even if a 1 ⊗ a 2 ∈ V .
Now a combination of (3.15), (3.16) and (3.17) yields

a 1 (0) a 1 φ k 2 L 2 ((0,1)) + a 2 φ ℓ 2 L 2 ((0,1)) ≤ C a 1 H 1 (0,1) a 2 H 1 (0,1) e λ kℓ τ 2 /2 ∂ ν u L 2 (Σ1) ,
where φ k (s) = √ 2 cos((k + 1/2)πs). This and the fact that m ≤ a j (0) and a j H 1 ((0,1)) ≤ M imply

a 1 φ k 2 L 2 ((0,1)) + a 2 φ ℓ 2 L 2 ((0,1)) ≤ C M 2 m e λ kℓ τ 2 /2 ∂ ν u L 2 (Σ1) .
Hence, where j = 1 or 2,

a j φ k 2 L 2 ((0,1)) ≤ C M 2 m e k 2 τ 2 π 2 ∂ ν u L 2 (Σ1) .
Let

a k j = 1 0 a j (x)φ k (x)dx, j = 1, 2.
Since

|a k j | = 1 0 a j (x)φ k (x)dx ≤ a j φ k L 1 ((0,1)) ≤ a j φ k L 2 ((0,1)) ,
we get

(a k j ) 2 ≤ C M 2 m e k 2 τ 2 π 2 ∂ ν u L 2 (Σ1) .
On the other hand

∂ ν u L 2 (Σ1) = Λ a (φ kl ) -Λ 0 (φ kl ) L 2 (Σ) ≤ Ck 2 Λ a -Λ 0 . Hence (3.18) (a k j ) 2 ≤ C M 2 m e k 2 (τ 2 π 2 +1) Λ a -Λ 0 .
Let q = M 2 m and α = τ 2 π 2 + 2. We obtain in a straightforward manner from (3.18)

|k|≤N (a k j ) 2 ≤ Cqe αN 2 Λ a -Λ 0 .
Consequently,

a j 2 L 2 ((0,1)) ≤ |k|≤N (a k j ) 2 + 1 N 2 |k|>N k 2 (a k j ) 2 ≤ C qe αN 2 Λ a -Λ 0 + a j 2 H 1 ((0,1)) N 2 ≤ C qe αN 2 Λ a -Λ 0 + M 2 N 2 ≤ CM 2 1 m e αN 2 Λ a -Λ 0 + 1 N 2 .
That is

(3.19) a j 2 L 2 ((0,1)) ≤ CM 2 1 m e αN 2 Λ a -Λ 0 + 1 N 2 .
Assume that Λ a -Λ 0 ≤ δ = me -α . Let then N 0 ≥ 1 be the greatest integer so that

C m e αN 2 0 Λ a -Λ 0 ≤ 1 N 2 0 . Using 1 m e α(N0+1) 2 Λ a -Λ 0 ≤ 1 (N 0 + 1) 2 , we find (2N 0 ) 2 ≥ (N 0 + 1) 2 ≥ 1 α + 1 ln m Λ a -Λ 0 . This estimate in (3.19) with N = N 0 gives (3.20) a j L 2 ((0,1)) ≤ 2C √ α + 1M ln m -1 Λ a -Λ 0 -1/2 .
When Λ a -Λ 0 ≥ δ, we have

(3.21) a j L 2 ((0,1)) ≤ M δ Λ a -Λ 0 .
In light of (3.20) and (3.21), we find a constants c > 0, that can depend only on τ , so that

a j L 2 ((0,1)) ≤ cM ln m -1 Λ a -Λ 0 -1/2 + m -1 Λ a -Λ 0 .
Appendix A

We prove the following lemma Lemma A.1. Let 1/2 < α ≤ 1 and a ∈ C α ([0, 1]). Then the mapping f → af defines a bounded operator on H 1/2 ((0, 1)).

Proof. We recall that H 1/2 ((0, 1)) consists in functions f ∈ L 2 ((0, 1)) with finite norm

f H 1/2 ((0,1)) = f 2 L 2 ((0,1)) + 1 0 1 0 |f (x) -f (y)| 2 |x -y| 2 dxdy 1/2 . Let a ∈ C α ([0, 1]). We have |a(x)f (x) -a(y)f (y)| 2 |x -y| 2 ≤ a 2 L ∞ (0,1) |f (x) -f (y)| 2 |x -y| 2 + |f (y)| 2 [a] 2 α |x -y| 2(1-α) , where [a] α = sup{|a(x) -a(y)||x -y| -α ; x, y ∈ [0, 1], x = y}.
Using that 1/2 < α ≤ 1, we find that x → |x -y| -2(1-α) ∈ L 1 ((0, 1)), y ∈ [0, 1], and

1 0 dx |x -y| 2(1-α) ≤ 1 2α -1 , y ∈ [0, 1].
Hence af ∈ H 1/2 ((0, 1)) with af H 1/2 ((0,1)) ≤ 1 2α -1 a C α ([0,1]) f H 1/2 ((0,1)) . In other words, A a is dissipative. We complete the proof by showing that A a is onto implying that A a is m-dissipative. To this end we are going to show that for each (f, g) ∈ V × L 2 (Ω), the problem w = f, -∆v = g. In light of the fact ψ → Ω |∇ψ| 2 dx 1/2 defines an equivalent norm on V , we can apply Lax-milgram's lemma. We get that there exists a unique v ∈ V satisfying

Ω ∇v • ∇ψdx = Ω gψdx - Γ1 awψdσ, ψ ∈ V.
From this identity, we deduce in a standard way that -∆v = g and ∂ ν v = -aw on Γ 1 . The proof is then complete

2 Γ1a 2 Γ1a |w 1 | 2

 2212 Herea C α ([0,1]) = a L ∞ ((0,1)) + [a] α .Appendix BWe give the proof of the following lemma Lemma B.1. Let a ∈ A and A a be the unbounded operator defined on V × L 2 (Ω) byA a = (w, ∆v), D(A a ) = {(v, w) ∈ V × V ; ∆v ∈ L 2 (Ω) and ∂ ν v = -aw on Γ 1 }. Then A a is m-dissipative. Proof. Let •, • be scalar product in V × L 2 (Ω). That is (v 1 , w 1 ), (v 2 , w 2 ) = Ω ∇v 1 • ∇v 2 dx + Ω w 1 w 2 dx, (v j , w j ) ∈ V × L 2 (Ω), j = 1, 2.For (v 1 , w 1 ) ∈ D(A a ), we haveA a (v 1 , w 1 ), (v 1 , w 1 ) = (w 1 , ∆v 1 ), (v 1 , w 1 ) (B.1) = Ω ∇w 1 • ∇v 1 dx + Ω ∆v 1 w 1 dxApplying twice Green's formula, we getΩ ∇w 1 • ∇v 1 dx = -Ω w 1 ∆v 1 dx + Γ1 w 1 ∂ ν v 1 dσ, (B.2) Ω ∆v 1 w 1 dx = -Ω ∇v 1 • ∇w 1 dx -Γ1 aw 1 w 1 dσ. (B.3)We take the sum side by side of identities (B.2) and (B.3). Using that ∂ ν v 1 = -aw 1 on Γ 1 we obtainΩ ∇w 1 • ∇v 1 dx + Ω ∆v 1 w 1 dx = -Ω w 1 ∆v 1 dx -Ω ∇v 1 • ∇w 1 dx -|w 1 | 2 dσ = -(v 1 , w 1 ), A a (v 1 , w 1 ) -dσ.This and (B.1) yieldℜ A a (v 1 , w 1 ), (v 2 , w 2 ) = -Γ1 a |w 1 | 2 dσ ≤ 0.

  has a unique solution (v, w) ∈ D(A a ).