Éric Bonnetier
email: eric.bonnetier@univ-grenoble-alpes.fr

Pierre Etoré
email: pierre.etore@univ-grenoble-alpes.fr

Miguel Martinez
email: miguel.martinez@univ-eiffel.fr

A probabilistic representation of the solution to a 1D evolution equation in a medium with negative index

Keywords: Negative Index Materials, Evolution equations, Spectral representation of semigroups, Skew Brownian motion, Pseudo processes

In this work we investigate a 1D evolution equation involving a divergence form operator where the diffusion coefficient inside the divergence is changing sign, as in models for metamaterials. We focus on the construction of a fundamental solution for the evolution equation, which does not proceed as in the case of standard parabolic PDE's, since the associated second order operator is not elliptic. We show that a spectral representation of the semigroup associated to the equation can be derived, which leads to a first expression of the fundamental solution. We also derive a probabilistic representation in terms of a pseudo Skew Brownian Motion (SBM). This construction generalizes that derived from the killed SBM when the diffusion coefficient is piecewise constant but remains positive. We show that the pseudo SBM can be approached by a rescaled pseudo asymmetric random walk, which allows us to derive several numerical schemes for the resolution of the PDE and we report the associated numerical test results.

Introduction

Over the last two decades, negative index materials (NIM) have drawn considerable attention, due to the spectacular way in which electromagnetic, acoustic or elastic waves may propagate in such media. Composite materials built as mixtures of NIM's and classical dielectric material may indeed show resonant effects when excited at certain frequencies, in which electromagnetic fields may concentrate near the interfaces, providing a way to channel the fields. The amplitude of the fields may be enhanced by several orders of magnitude in the neighborhood of particles with negative permittivity or permeability. The possibility of localizing and concentrating waves has made NIM's a subject of great interest for many applications, for instance in communication and medical imaging.

From the point of view of mathematical modeling, a large part of the work on NIM's has focused on the so-called electrostatic approximation, for a composite medium made of inclusions (or phases) of NIM's embedded into a homogeneous dielectric matrix phase. This means that the time-harmonic Maxwell system is reduced to a diffusion for one of the components of either the electric or magnetic fields, provided the geometry of the device is assumed to have a direction of invariance, and provided the typical dimensions of the inhomogeneities are small with respect to the incident wavelength. In the simplest case, one remains with a transmission equation of the form div(A(x)∇u(x)) = f + boundary conditions [START_REF] Chesnel | T-coercivity for scalar interface problems between dielectrics and metamaterials[END_REF] in a bounded domain D (one could also consider this equation in the whole space), and where the conductivity A(x) takes negative values in the inclusions of NIM's and a positive value in the dielectric medium in which they are embedded.

Because the conductivity changes sign in the domain, the bilinear form associated to the above PDE a(u, v) = D A(x)∇u(x)∇v(x), fails to be coercive in the natural Sobolev space H 1 (D) or in the appropriate subspace H that accounts for the imposed boundary conditions. Thus, one cannot invoke the Lax-Milgram Lemma to show existence of solutions to [START_REF] Chesnel | T-coercivity for scalar interface problems between dielectrics and metamaterials[END_REF]. In many cases though, depending on the geometry of the NIM inhomogeneities, one can show that the form a is T-coercive : there exists an invertible operator T : H -→ H such that (u, v) ∈ H × H -→ a(u, T v) is coercive. In this case, applying the Lax-Milgram Lemma to the latter bilinear form yields well-posedness of the PDE [START_REF] Chesnel | T-coercivity for scalar interface problems between dielectrics and metamaterials[END_REF][START_REF]T-coercivity for the Maxwell problem with sign-changing coefficients[END_REF][START_REF] Nguyen | Limiting absorption principle and well-posedness for the Helmholtz equation with sign changing coefficients[END_REF], except possibly for some values of the negative conductivities, for which the associated operator is not invertible and may even loose its Fredholm character.

T-coercivity is thus a simple way of paliating the lack of coercivity. One of our objectives is to investigate whether T-coercivity also grants that the operator in [START_REF] Chesnel | T-coercivity for scalar interface problems between dielectrics and metamaterials[END_REF] shares other characteristics of elliptic operators. In particular, given the relation of the latter to stochastic processes, we would like to answer the following questions : Does there exists a probabilistic representation of the solutions to (1) akin a Feynman-Kac formula ? What is the nature of the underlying stochastic process ? Can one design Monte Carlo type numerical schemes to approximate the solutions to (1) ?

In this work, we investigate these questions in a simple one-dimensional situation. Since the probabilistic interpretation of an elliptic PDE is strongly related to the parabolic equation whose infinitesimal generator is the associated elliptic operator, we consider a parabolic version of (1) of the form

         A(x)∂tu(t, x) = 1 2 ∂x A(x)∂xu(t, x) , x ∈ I, t > 0, u(0, x) = u0(x) x ∈ I, u(t, ±1) = 0, t > 0, (2)
where I is the interval (-a, a), and where the conductivity A is defined by A(x) = k, x ∈ I -:= (-a, 0) 1, x ∈ I + := (0, a).

In other words, we assume that a dielectric with (positive) conductivity 1 fills in the right part of the interval I, whereas when k < 0, the left part is filled with a negative index material. Note that the time derivative of u is multiplied by the conductivity A(x), and that weak solutions to (2) satisfy    ∂tu(t, x) = 1 2 u ′′ (t, x), x ∈ (-a, 0) ∪ (0, a), t > 0 u(0, x) = u0(x), x ∈ (-a, a) u(t, a)

= u(-t, a) = 0 [START_REF] Allaire | Analyse numérique et optimisation: une introduction à la modélisation mathématique et à la simulation numérique[END_REF] and the transmission conditions at the interface

u(t, 0 -) = u(t, 0 + ), k u ′ (t, 0 -) = u ′ (t, 0 + ). ( 4
)
When k > 0, the coefficient in front of the time derivative of u does not bring significant changes to the usual parabolic setting, and existence and uniqueness of solutions are guaranteed. When k < 0 however, a 'parabolic' version of (1) of the form B(x)∂tu(t, x) = 1 2 ∂x A(x)∂xu(t, x) , with initial datum u0 ∈ L 2 (I), may only have solutions with finite energy (i.e. in the space L 2 (0, ∞, H 1 0 (I))) if for any x ∈ I, B(x) and A(x) have the same sign [START_REF] Bonnetier | Spectral bases for graphs mixing dielectric and metamaterials[END_REF], so here we simply choose B = A. Our objective is to construct a fundamental solution for [START_REF]T-coercivity for the Maxwell problem with sign-changing coefficients[END_REF] which can be interpreted as the transition function of a stochastic process, or rather as we explain below, of a pseudo-process.

We explicitely characterize the associated measure, and inquire how it can be decomposed so as to build probabilistic-like numerical schemes to compute the solutions to (4).

When k > 0, it is well known that the fundamental solution of (2) can be interpreted as the transition function of a stochastic process : this can be shown for example by using the stochastic calculus for Dirichlet forms (see for e.g. [START_REF] Fukushima | Dirichlet forms and symmetric Markov processes[END_REF]). The solution u(t, x) of ( 3)-( 4) can be represented via a Feynman-Kac formula as the expectation E x (u0( Xt)) where X stands for a Skew Brownian motion (SBM) of parameter β = 1-k 1+k , killed when reaching -a or a (we have denoted E x (•) = E(•| X0 = x); see e.g. [START_REF] Karatzas | Brownian motion and stochastic calculus[END_REF], and [START_REF] Harrison | On skew Brownian motion[END_REF] for an account on the Skew Brownian motion). Thus the transition probability density function (transition function, in short) of X provides a fundamental solution of (2).

When k < 0, the skew Brownian motion may not be well defined. However, the above construction formally provides a solution to (3)-( 4) that generalizes the probabilistic representation. With this perspective, we consider a pseudo expectation involving a pseudo Skew Brownian Motion, which we properly define as a pseudo process, i.e., as a family of measurable functions Y = (Yt) defined on a space (Ω, F) endowed with a signed measure P with P(Ω) = 1. The functions Y = (Yt) are called pseudo-random variables Y = (Yt), and (Ω, F, P) is referred to as a pseudo-probability space. See [START_REF] Lachal | From pseudorandom walk to pseudo-brownian motion: First exit time from a one-sided or a two-sided interval[END_REF] and the references therein for an account on pseudo random variables. A pseudo process is defined mainly by its transition function, which allows to compute pseudo transition probabilities.

When I = R (a = +∞), the construction of a pseudo SBM can be successfully carried out, taking advantage of the fact that in the case k > 0, the transition function of the SBM is known ( [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF]). It follows that one can indeed construct a fundamental solution of (2) and define a pseudo SBM in the general case k ∈ R * \ {-1}. Moreover, the (pseudo) probabilistic representation obtained in this way naturally lends itself to numerical approximation, via the convergence of pseudo skew random walks to pseudo skew Brownian motions.

The bounded case, where I = (-a, a), a < +∞, is more involved, since it is not clear how one could define a 'killed' pseudo skew Brownian motion as a pseudo-probabilistic process. This is mainly because 'killing' is a trajectorial operation, and because the trajectories of a pseudo processes do not have a clear meaning. We are however able to extend the representation of the solution u(t, x) of ( 3)-( 4) in terms of the fundamental solution, when k ∈ R * \ {-1}. One preliminary step in this direction, consists in computing the transition function of the killed skew Browian X, by probabilistic arguments in the case k > 0. This work is organized as follows. In Section 2, we construct solutions to (2) in L 2 (0, ∞, H 1 0 (I)) for any k ̸ = -1, with the help of the eigenfunctions of the bilinear form a(u, v) = I A(x)∂xu∂xv which are shown to form a basis of the space H 1 0 (I). Section 3 focuses on the unbounded case when I = R (or a = ∞). We show that the solutions to (2) can be obtained as convolutions of the initial datum with a kernel p(t, x, y) which has a pseudo probabilistic interpretation, and to which we associate a pseudo skew-Brownian motion. We further show that such pseudo skew-Brownian motion can be approximated by pseudo-skew random walks. In Section 4, we consider the bounded case (a < ∞). Finally, in Section 5, we construct several numerical schemes for solving (2) based on these developments and report numerical results obtained with these schemes.

Spectral representation of the semigroup

In this section, we consider the eigenvalue problem : find λ ∈ C and u ∈ H 1 0 (I) such that

-A(x)u ′ (x) = λ 2 A(x)u(x), in I, (5)
Seeking a solution to [START_REF] Étoré | On random walk simulation of one-dimensional diffusion processes with discontinuous coefficients[END_REF] in the form

u(x) = a1 cos(λx) + b1 sin(λx), x ∈ I -, a2 cos(λx) + b2 sin(λx), x ∈ I + , (6)
and expressing the transmission and boundary conditions (4), one is led to solving the linear system

    1 0 -1 0 0 k 0 -1 cos(λa) -sin(λa) 0 0 0 0 cos(λa) sin(λa)         a1 b1 a2 b2     = 0,
from which one obtains the dispersion relation sin(2λa)(k + 1) = 0.

Assuming k ̸ = 1, the associated eigen-elements can be grouped in two families : a set of even functions

f k,n (x) = cos (2n -1)π 2a x , x ∈ I, associated to λn = (2n -1)π 2a
, n ≥ 1, and the set of functions

g k,n (x) =    sin( nπ a x) if x ∈ I -, k sin( nπ a x) if x ∈ I + , associated to µn = nπ a , n ≥ 1.
When k > 0, the bilinear form u, v

→ I A(x)u(x)v(x) dx is a scalar product in L 2 (I), while the form u, v ∈ H 1 0 (I) → a(u, v) = I A(x)u ′ (x)v ′ (x) dx,
is coercive and symmetric in H 1 0 (I), and it is well-known that the solutions to (5) form a basis of L 2 (I) and of H 1 0 (I). The next Proposition shows that this is still the case when k < 0. Proposition 2.1. When k < 0, the functions (f k,n , g k,n ) n≥1 form a Hilbert basis of L 2 (I) and of H 1 0 (I). Proof. Let k < 0. In view of the above remark, the functions (f -k,n , g -k,n ) n≥1 are a basis of L 2 (I), and one has the orthogonal decomposition (with respect to the scalar product associated to -k)

L 2 (I) = H f ⊕ Hg,
where H f (resp. Hg) is the vector space generated by the f -k,n 's (resp. by the g -k,n 's). Consider the mapping T : L 2 (I) -→ L 2 (I) defined by

T u(x) =    u(x) if u ∈ H f , u(x) if u ∈ Hg and x ∈ I -, -u(x)
if u ∈ Hg and x ∈ I + .

As T • T = I, T is an isomorphism, thus the basis (f -k,n , g -k,n ) is transformed into a basis of L 2 (I). It is easy to check that (T f -k,n , T g -k,n ) = (f k,n , g k,n ). The same arguments show that the (f k,n , g k,n )'s form a basis of H 1 0 (I). From now on, we assume that k ∈ R * \ {-1}, and we drop the index k in the notation of the basis functions f k,n and g k,n . Note that the functions fn, gn, n ≥ 1 satisfy the following relations:

λ 2 p I A(x)fp(x)φ(x) = I A(x)f ′ p (x)φ ′ (x) = 0,
for φ = fq, q ̸ = p or φ = gq, q ≥ 1, and similarly

µ 2 p I A(x)gp(x)φ(x) = I A(x)g ′ p (x)φ ′ (x) = 0,
for φ = gq, q ̸ = p or φ = fq, q ≥ 1. This corresponds to orthogonality or pseudo-orthogonality properties, depending on the sign of k.

It follows that any function u0 ∈ L 2 (I) can be written in the form

u0(x) = n≥1 anfn(x) + bngn(x), (8)
where the coefficients an and bn are given by an = I A(x)u0(x)fn(x) dx

I A(x)|fn(x)| 2 dx , bn = I A(x)u0(x)gn(x) dx I A(x)|gn(x)| 2 dx . (9)
In particular, one can check that

I A(x)|fn(x)| 2 dx = a(k + 1) 2 , I A(x)|gn(x)| 2 dx = ak(k + 1) 2 . ( 10
)
Remark 2.2. Note that if k = -1 it is impossible to provide the forthcoming representations [START_REF] Karoui | Les outils stochastiques des marchés financiers : une visite guidée de einstein à black-scholes[END_REF] and [START_REF] Lachal | From pseudorandom walk to pseudo-brownian motion: First exit time from a one-sided or a two-sided interval[END_REF], with the an's and bn's computed by ( 9) and [START_REF] Karatzas | Brownian motion and stochastic calculus[END_REF]. Therefore our assumption k ∈ R * \ {-1}.

We now focus on the evolution problem (2) assuming that u0 ∈ L 2 (I). Decomposing u0 on the basis of eigenfunctions as in [START_REF] Harrison | On skew Brownian motion[END_REF], it is easy to check that

u(t, x) = Ptu0(x) = n≥1 ane -λ 2 n t/2 fn(x) + bne -µ 2 n t/2 gn(x), (11)
is a weak solution to (2) in the sense that u ∈ L 2 (0, ∞, H 1 0 (I)), ∂tu ∈ L 2 (0, ∞, H -1 (I)) and

∀ v ∈ H 1 0 (I), I 2A(x)∂tu(t, x)v(x) + I A(x)∂xu(t, x)∂xv(x) = 0, (12)
and u(0, x) = u0(x), a.e. x ∈ I. With this definition, using the fn, gn's as test functions, the weak solution to (2) is easily seen to be unique. Note that the family (Pt) forms a semigroup on L 2 (I). Finally, the expression (11) can be rewritten as

u(t, x) = I u0(y)p(t, x, y) dy,
where the kernel has the form

p(t, x, y) = n≥1 A(y) 2 a(k + 1) fn(y)fn(x) e -λ 2 n t/2 + 2 ak(k + 1) gn(y)gn(x) e -µ 2 n t/2 . ( 13
)
We give a probabilistic derivation of the fundamental solution for (2) in the next sections.

3 A probabilistic construction of the fundamental solution for the evolution equation on R

In this section, we consider the case when I = R, and obtain the expression of the fundamental solution to (2) following a probabilistic construction. More precisely, we consider the Cauchy problem

   2A(x)∂tu(t, x) = (A(x)u ′ (t, x)) ′ , t > 0, x ∈ R u(0, x) = u0(x), u(t, •) ∈ H 1 (R) and Au ′ (t, •) ∈ H 1 (R) t > 0, (14)
or equivalently the PDE

∂tu(t, x) = 1 2 u ′′ (t, x), x ∈ (-∞, 0) ∪ (0, ∞), t > 0,
with the initial condition u(0, x) = u0(x), x ∈ R and the radiation and transmission conditions

lim |x|→∞ u(t, x) = 0, u(t, 0 -) = u(t, 0 + ) ku ′ (t, 0 -) = u ′ (t, 0 + ). ( 15
)
We first introduce a few notations.

Notations. We denote by g(t, x, y)

= 1 √ 2πt exp -|y-x| 2 2t
the density of a N (x, t). Let b < c ∈ R. The transition density of the Brownian motion killed at the points b or c is given by p (b,c)

W (t, x, y) := ∞ n=-∞ g(t, x, y -2n(c -b)) -g(t, x, 2b -y -2n(c -b)) (16)
(see [START_REF] Nikolaevitch | Handbook of brownian motion : facts and formulae, Probability and its applications[END_REF], Appendix I, N� 6). Note that one has for any y ∈ (b, c),

lim x-→ > b p (b,c) W (t, x, y) = 0 and lim x-→ < c p (b,c) W (t, x, y) = 0 (17)
3.1 The fundamental solution when k > 0

We first assume that k > 0, so that equation ( 2) is parabolic. It is known that in this case, under mild assumptions on the initial condition u0 (for example u0 is continuous and bounded, see [START_REF] Lejay | On the constructions of the skew Brownian motion[END_REF]) the solution to the above Cauchy problem is given by

u(t, x) = E x [u0(Xt)], ( 18
)
where X is the SBM with parameter β = 1-k 1+k ∈ (-1, 1). We refer to the survey [START_REF] Lejay | On the constructions of the skew Brownian motion[END_REF] for the definition and main properties of the SBM. In particular X solves the Stochastic Differential Equation (SDE) with local time dXt = dWt + βdL 0 t (X), [START_REF] Revuz | Continuous martingales and Brownian motion[END_REF] where W denotes a standard Brownian motion driving the SDE, and where L 0 t (X) is the symmetric local time at the point zero and at time t of X. Note that |β| < 1 ensures the existence of X, see e.g. [START_REF] Lejay | On the constructions of the skew Brownian motion[END_REF]. The SBM behaves like a Brownian motion, except at the times when it touches zero, at which its dynamics are biased by the term of local time in [START_REF] Revuz | Continuous martingales and Brownian motion[END_REF]. In particular we have: [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF]] Let β ∈ (-1, 1) and let X be the solution to [START_REF] Revuz | Continuous martingales and Brownian motion[END_REF]. Under P 0 one has:

Lemma 3.1. [Walsh,
i) The process |X| is distributed as a reflecting Brownian motion |W | (starting from zero).

ii) The processes (sign(Xt)) and (|Xt|) are independent.

In addition, for any t > 0 one has P 0 (Xt > 0) = 1+β 2 . From this Lemma and the reflection principle for the Brownian motion, Walsh was able to give an explicit expression of the transition probability density p(t, x, y) of the SBM, in the form [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF]

p(t, x, y) =                    (1 -β)g(t, x, y) if x ≥ 0, y < 0 g(t, x, y) + βg(t, x, -y) if x > 0, y > 0 (1 + β)g(t, x, y) if x ≤ 0, y > 0 g(t, x, y) -βg(t, x, -y) if x < 0, y < 0. (20)
Note that (20) can also be written as

p(t, x, y) =                    (1 -β)g(t, x, y) if x ≥ 0, y < 0 -β[g(t, x, y) -g(t, x, -y)] + (1 + β)g(t, x, y) if x > 0, y > 0 (1 + β)g(t, x, y) if x ≤ 0, y > 0 β[g(t, x, y) -g(t, x, -y)] + (1 -β)g(t, x, y) if x < 0, y < 0. ( 21
)
where g(t, x, y) -g(t, x, -y) =: ǧ(t, x, y), x, y > 0 (resp. x, y < 0) can be interpreted as the transition function of a Brownian motion on (0, ∞) (resp. (-∞, 0)) killed at zero (see [START_REF] Nikolaevitch | Handbook of brownian motion : facts and formulae, Probability and its applications[END_REF] and Remark 4.2). Yet another way to rewrite [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF] is to consider ĝ(t, x, y) := g(t, x, y) + g(t, x, -y), x, y > 0, the transition function of a reflected Brownian motion on [0, ∞) (see [START_REF] Nikolaevitch | Handbook of brownian motion : facts and formulae, Probability and its applications[END_REF]), which yields

p(t, x, y) =                    (1 -α)[ĝ(t, x, -y) -ǧ(t, x, -y)] if x ≥ 0, y < 0 αĝ(t, x, y) + (1 -α)ǧ(t, x, y) if x > 0, y > 0 α[ĝ(t, -x, y) -ǧ(t, -x, y)] if x ≤ 0, y > 0 (1 -α)ĝ(t, -x, -y) + αǧ(t, -x, -y) if x < 0, y < 0, (22)
where we have set

α = (1 + β)/2 = 1/(1 + k). (23)
It follows from [START_REF] Raviart | Introduction à l'analyse numérique des équations aux dérivées partielles[END_REF] that

u(t, x) = E x [u0(Xt)] = R u0(y)p(t, x, y)dy, (24)
so that p(t, x, y) identifies with the fundamental solution of ( 14).

This fact also follows by mere computation: indeed one can easily check that the function p, defined by [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF], satisfies for any y ∈ R

∂tp(t, x, y) = 1 2 ∂ 2 xx p(t, x, y), ∀x ∈ (-∞, 0) ∪ (0, ∞). (25)
Using ∂xg(t, x, y) = y-x t g(t, x, y) one can check the transmission conditions

p(t, 0-, y) = p(t, 0+, y), k∂xp(t, 0-, y) = ∂xp(t, 0+, y), (26)
and the radiation condition is also easy to check:

lim |x|→∞ p(t, x, y) = 0. ( 27
)
Taking derivatives of the integral on the right-hand side of (24) and using ( 25)-( 27) shows that p is indeed the fundamental solution of ( 14).

3.2

The fundamental solution in the case k < 0 and the pseudo SBM

When k ∈ R * -\ {-1}, setting β = 1-k 1+k
, we may again define a function p as in [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF]. It is easy to check that this function solves (25)-( 27), so that

u(t, x) = R u0(y)p(t, x, y) dy,
is a solution to the Cauchy problem [START_REF] Gall | One-dimensional stochastic differential equations involving the local times of the unknown process, Stochastic analysis and applications[END_REF]. Consequently, p(t, x, y) is a fundamental solution to ( 14) also in the case k ∈ R * -\ {-1}. Note that β is not in (-1, 1), so that p(t, x, y) dy is only a signed measure, and ( 19) does not define a SBM (see e.g. [START_REF] Gall | One-dimensional stochastic differential equations involving the local times of the unknown process, Stochastic analysis and applications[END_REF]). Stochastic Differential equations of the type

dXt = σ(Xt)dWt + βdL 0 t (X) (28)
with β / ∈ (-1, 1) have been addressed for example in [START_REF] Kopytko | An example of a stochastic differential equation with the property of weak non-uniqueness of a solution[END_REF]. But this latter work does not allow to take σ ≡ 1, β / ∈ (-1, 1) and to get a weak solution to (28) (indeed to touch β / ∈ (-1, 1) the coefficient σ has to be different from 1). However, we can still define a function p(t, x, y) by ( 20) or ( 21), to which we can associate a pseudo SBM, as we describe below.

By pseudo-random variable, we mean a measurable function defined on a space (so-called pseudoprobability space) endowed with a signed measure with a total mass equal to 1. We observe that in the case k ∈ R * -\ {-1}, the signed asymmetric heat-type kernel p(t, x, y) defined in ( 21) is no longer positive everywhere, however it integrates to 1 w.r.t. dy. In view of (21), this is clearly the case when x = 0. When x > 0, we have

+∞ -∞ p(t, x, y)dy = -β +∞ 0 [g(t, x, y) -g(t, x, -y)]dy + (1 + β) +∞ 0 g(t, x, y)dy + (1 -β) 0 -∞ g(t, x, y)dy = β +∞ 0 g(t, x, -y)dy + +∞ 0 g(t, x, y)dy + (1 -β) 0 -∞ g(t, x, y)dy = +∞ -∞ g(t, x, y)dy + β +∞ 0 g(t, x, -y)dy -β 0 -∞ g(t, x, y)dy = 1 + β - -∞ 0 g(t, x, z)dz - 0 -∞ g(t, x, y)dy = 1.
A similar computation applies to the case x < 0. One may also check that the Chapman-Kolmogorov idendity for the family of transition kernels (p(t, x, y)dy

) x∈R is also preserved in this case k < 0, k ̸ = -1.
Hence, in accordance to the usual Markov rules, we may define the pseudo-skew Brownian motion as the pseudo-Markov process '(Xt) t≥0 ' associated with the signed asymmetric heat-type kernel p(t, x, y)dy defined in (21) -which is the fundamental solution to (14) also in the case k ∈ R * -\ {-1} -by : for t > 0 and 0 = t0 < t1 < . . . tm and x = x0, x1, . . . , xm, y ∈ R, P x (Xt ∈ dy) = p(t, x, y)dy and

P x (Xt 1 ∈ dx1, . . . , Xt m ∈ dxm) = m i=1 p(ti -ti-1, xi-1, xi)dxi.
Note that since pseudo-Markov processes are defined in terms of a signed measure, it is not clear how one could generalize the definition of the skew Brownian motion over all t ≥ 0 in this context. In particular, from a strict probabilistic point of view, the notion of trajectory for pseudo Markov processes indexed by continuous time does not have a clear meaning.

For more results on pseudo-Markov processes, we refer to e.g. [START_REF] Lachal | From pseudorandom walk to pseudo-brownian motion: First exit time from a one-sided or a two-sided interval[END_REF].

We now give a definition, inspired from [START_REF] Lachal | From pseudorandom walk to pseudo-brownian motion: First exit time from a one-sided or a two-sided interval[END_REF], for the convergence of a family of pseudo-processes ((Y ε t ) t≥0 )ε>0 towards a pseudo-process (Yt) t≥0 . Definition 3.2. Let ((Y ε t ) t≥0 )ε>0 denote a family of pseudo-processes and (Yt) t≥0 a pseudo-process. We write

(Y ε t ) t≥0 pseudo-w ------→ ε↘0 (Yt) t≥0 if ∀ℓ ∈ N * , ∀t1, . . . , t ℓ ≥ 0, ∀u1, . . . , u ℓ ∈ R, E exp i ℓ j=1 ujY ε t j ---→ ε↘0 E exp i ℓ j=1 ujYt j . ( 29
)
Remark 3.3. Note that in (29), the left hand side expectation symbol might depend on ε. But in order to avoid cumbersome notations, and as it will cause no ambiguity in the proofs, we simply denote it by E.

Convergence of the scaled pseudo asymmetric random walk to the pseudo SBM

In this section, we assume that k ∈ R * \ {-1}, and define α as in (23). If k > 0 then α ∈ (0, 1) and we deal with true processes and random walks.

If k ∈ R * -\ {-1} then α ∈ R \ (0, 1
) and we deal with pseudo processes and pseudo random walks. Our computations englobe both cases, but the results are new for the case α ∈ R \ (0, 1) (for the case α ∈ (0, 1) see e.g. [START_REF] Harrison | On skew Brownian motion[END_REF]).

Before getting deeper into our subject, let us give here a brief account of the several notations that we are going to use in the remaining of this section. * Regarding processes :

-W stands for a standard Brownian motion (classical symmetric).

-M stands for a standard symmetric random walk on Z.

-S stands for the pseudo random walk on Z that is α-skewed at 0. -X n stands for the properly normalized pseudo random walk α-skewed at 0. -X denotes the pseudo-process under study.

-The superscript † denotes killing at zero: this path operation will be only performed for classical processes (either the symmetric random walk M or the Brownian motion W ).

For example, it should be clear that W † denotes a standard Brownian motion killed when it hits zero for the first time. * Regarding functions :

-The subscript 0 (or superscript 0 ) will always refer to a quantity concerning the classical standard symmetric Z-valued random walk M . -The hat superscript 'ˆ' will always refer to quantities that are related to a reflected process at 0 (the standard symmetric random walk on Z or the standard Brownian motion). -The check superscript 'ˇ' will always refer to quantities that are related to a process with killing at 0 (the standard symmetric random walk on Z or the standard Brownian motion). -The letters ψ, Ψ denote characteristic functions.

-The order of variables in a characteristic function will always be : time / starting point / argument.

For example, it should be plain that ψ0(j; m0, v) refers to the characteristic function of the classical standard symmetric random walk killed at zero and starting from m0, taken at time j and evaluated at v, namely E m 0 [exp(iv M † j )]. These notations will be recalled when needed below. Let m0 ∈ Z denote an arbitrary integer. Constructed on some pseudo probability space (Ω, F, P m 0 ) we consider (Sn) n≥0 the pseudo skewed random walk on the integers starting P m 0 -a.s. from S0 = m0 with pseudo transition probabilities given by

P m 0 (Sn+1 = Sn + 1|S0, . . . , Sn) = α if Sn = 0 1/2 otherwise. ( 30
)
P m 0 (Sn+1 = Sn -1|S0, . . . , Sn) = 1 -α if Sn = 0 1/2 otherwise. ( 31
)
We attach to (Sn) n≥0 its natural filtration F S n n≥0 defined by F S n := σ (S k : k ≤ n) for n ≥ 0. The pseudo random sequence (Sn) n≥0 constructed likewise is a pseudo markovian random walk. Remark 3.4. Note that -contrary to the case of pseudo processes in continuous time -any F S n adapted functional of the pseudo skewed random walk (S k ) k≥0 can be expressed as a functional of (S0, S1, . . . , Sn), hence a finite dimensional pseudo distribution functional. Consequently, there is no ambiguïty in what is meant by the pseudo law under P m 0 of a pseudo skewed random walk on the whole path set Z N .

The following reflection principle is the key to get to the main result of this section. Lemma 3.5. For all j ∈ N * and m ∈ Z * ,

P 0 (Sj = m) =    αP 0 (|Sj| = m) if m > 0 (1 -α) P 0 (|Sj| = |m|) if m < 0. ( 32
)
Moreover, for any integer m0 > 0

P m 0 (Sj = m) =      αP m 0 (|Sj| = m) + (1 -α) P m 0 (|Sj| = m ; ∀n ∈ [[1, j]], |Sn| ̸ = 0) if m > 0 (1 -α) P m 0 (|Sj| = |m|) -P m 0 (|Sj| = m ; ∀n ∈ [[1, j]], |Sn| ̸ = 0) if m < 0
and for any integer m0 < 0,

P m 0 (Sj = m) =      (1 -α)P |m 0 | (|Sj| = |m|) + αP |m 0 | (|Sj| = |m| ; ∀n ∈ [[1, j]], |Sn| ̸ = 0) if m < 0 α P |m 0 | (|Sj| = m) -P |m 0 | (|Sj| = |m| ; ∀n ∈ [[1, j]], |Sn| ̸ = 0) if m > 0.
Proof. We only treat the case where the pseudo random walk starts from 0. The other cases, although tedious, can be analysed in the same fashion.

Let n ∈ N * and m ∈ Z * fixed.

Let us introduce Gn = sup (j ≤ n : Sj = 0) = sup (j ≤ n : |Sj| = 0) that is an F S n -measurable functional on the path space. It is easy to check that the pseudo-path (SG n +j ) j∈{0,...,n-Gn} has the same pseudo-law under the pseudo conditional probability P 0 (.|SG n +1 = 1) as the pseudo-path (|SG n +j |) j∈{0,...,n-Gn} under the original pseudo probability P 0 (this assertion makes sense remembering Remark 3.4). Observe also that P 0 (SG n +1 = 1) = α.

If m > 0, then a.s. on the set {Sn = m}, the value of Gn cannot be equal to n, so that Gn ≤ n -1 a.s. conditionally on this set. Hence, for m > 0 we have

P 0 (Sn = m) = P 0 ({Sn = m} ∩ {SG n +1 = 1}) = P 0 S Gn+(n-Gn) = m | SG n +1 = 1 P (SG n +1 = 1) = P 0 |S Gn+(n-Gn) | = m P 0 (SG n +1 = 1) = αP 0 (|Sn| = m) .
The case m < 0 is proved in a similar way. We have P 0 (SG n +1 = 1) = k/(1 + k) and it is easily checked that the pseudo-path (SG n +j ) j∈{0,...,n-Gn} has the same pseudo-law under the conditionnal probability P 0 (.|SG n+1 = -1) as the pseudo-path (-|SG n +j |) j∈{0,...,n-Gn} under P 0 (again, this assertion makes sense because of Remark 3.4).

For an arbitrary m0 ∈ Z, we denote by (M k ) k≥0 a standard symmetric random walk on Z constructed on (Ω, F, P m 0 ) such that M0 = m0 P m 0 -a.s.

Let bE (Z, C) denote the set of bounded complex valued functions defined on Z. We introduce the following family of operators T 0 j j∈N acting from bE (Z, C) to bE (Z, C) and defined by

T 0 j f : m →    αE m f (|Mj|) + (1 -α)E m f (-|Mj|) + (1 -α) E m f (|M † j |) -E m f (-|M † j |) if m ≥ 0 αE |m| f (|Mj|) + (1 -α)E |m| f (-|Mj|) + α E |m| f (-|M † j |) -E |m| f (|M † j |) if m < 0. ( 33
)
We have the following lemma.

Lemma 3.6. For any function f ∈ bE (Z, C),

E m (f (Sj)) = T 0 j f (m).
Proof. Let m0 ∈ N arbitrary. It is plain from the definition of the transitions of the pseudo skewed random walk (Sn) n≥0 (given in (30)) that (|Sn|) n≥0 and (|Mn|) n≥0 share the same pseudo law on the path space N N under P m 0 : this assertion makes sense recalling Remark 3.4. In particular, we are allowed to rewrite the equalities of Lemma 3.5 in the following manner

P m 0 (Sj = m) =    αP m 0 (|Mj| = m) + (1 -α) P m 0 (|Mj| = m ; ∀n ∈ [[1, j]], |Mn| ̸ = 0) if m > 0 (1 -α) (P m 0 (|Mj| = |m|) -P m 0 (|Mj| = |m| ; ∀n ∈ [[1, j]], |Mn| ̸ = 0)) if m < 0. ( 34
)
For the same reasons, we also have P m 0 (Sj = 0) = P m 0 (|Sj| = 0) = P m 0 (|Mj| = 0). The announced equality follows then easily by taking the expectation for f ∈ bE (Z, C). Indeed, we have

E m 0 (f (Sj)) = m∈Z f (m)P m 0 (Sj = m) = f (0)P m 0 (Sj = 0) + m∈N * f (m)P m 0 (Sj = m) + m∈N * f (-m)P m 0 (Sj = -m) .
Now using (34) (and the definition of killing), we find

E m 0 (f (Sj)) = f (0)P m 0 (|Mj| = 0) + α m∈N * f (m)P m 0 (|Mj| = m) + (1 -α) m∈N * f (-m)P m 0 (|Mj| = m) + (1 -α) m∈N * f (m)P m 0 |M † j | = m - m∈N * f (-m)P m 0 |M † j | = m = T 0 j f (m0).
The same kind of arguments may be invoqued for negative m0 and the result of the lemma follows.

Analogously to the above, let us introduce (Wt) t≥0 a standard Brownian motion constructed on (Ω, F, P x ) starting from x (i.e. W0 = x, P x -a.s).

Let bE (R, C) stand for the set of Borel bounded complex valued functions defined on R. Similarly as before, we introduce the following family of operators (Tt) t≥0 acting from bE (R, C) to bE (R, C) and defined by

Ttf : x →    αE x f (|Wt|) + (1 -α)E x f (-|Wt|) + (1 -α) E x f (|W † t |) -E x f (-|W † t |) if x ≥ 0 αE |x| f (|Wt|) + (1 -α)E |x| f (-|Wt|) + α E |x| f (-|W † t |) -E |x| f (|W † t |) if x < 0. ( 35
) Lemma 3.7. For any f ∈ bE (R, C), E x (f (Xt)) = Ttf (x).
Proof. Suppose x > 0. We use ( 22) and check at once that

E x (f (Xt)) = +∞ -∞ f (y)p(t, x, y)dy = α +∞ 0 f (y)ĝ(t, x, y)dy + (1 -α) +∞ 0 f (y)ǧ(t, x, y)dy + (1 -α) 0 -∞ f (y)ĝ(t, x, -y)dy - 0 -∞ f (y)ǧ(t, x, -y)dy = α +∞ 0 f (y)ĝ(t, x, y)dy + (1 -α) +∞ 0 f (y)ǧ(t, x, y)dy + (1 -α) +∞ 0 f (-y)ĝ(t, x, y)dy - +∞ 0 f (-y)ǧ(t, x, y)dy = Ttf (x). ( 36
)
The same type of arguments may be invoqued for negative x and the result of the lemma follows.

The main result of this section is the following statement.

Proposition 3.8. The rescaled asymmetric random walk converges in the sense of Definition 3.2 to the pseudo SBM :

X n := n -1/2 S ⌊nt⌋ t≥0 pseudo-w -→ (Xt) t≥0 as n → +∞.
Proof. Our main concern is to ensure that all arguments for the convergence in the classical case still hold true in our pseudo probability context. Since we deal with pseudo stochastic processes, the difficulties are two fold:

-we are not allowed to use a Skorokhod embedding, a tool that is often used for the convergence of random walks;

-we are not allowed to perform transformations directly on the path space such as reflection or killing for pseudo processes.

Note also that, except in the case α = 1/2, the skewed random walk and the skewed Brownian motion (classical or pseudo) do not have independent increments. Following [START_REF] Harrison | On skew Brownian motion[END_REF] in the classical case where α ∈ (0, 1), the idea is to observe the similarity of structure shared by the family of operators (T 0 j ) and (Tt) (whose definition involve only classical processes) and take advantage of known results concerning the convergence of random walks.

We only give the main arguments of the proof and leave the computational details to the reader.

Let us fix x ≥ 0, t ≥ 0, u ∈ R and set mn := ⌊ √ nx⌋. As mentionned above, in the whole proof ψ0(j; mn, u) := E mn (exp (iu |Mj|)) , resp. ψ0(j; mn, u) := E mn (exp(iu|M † j |))

stands for the characteristic function of the classical standard symmetric random walk on Z reflected at 0 (resp. killed at 0) starting from ⌊ √ nx⌋ and taken at time j and evaluated at u. Similarly ψ(t; x, u)

:= E x (e iu|W t | ), resp. ψ(t; x, u) := E x (e iu|W † t | ) (37)
stands for the characteristic function of a standard Brownian motion (Ws) reflected at 0 on [0, ∞) (resp. killed at 0) at time t and evaluated at u. Let us denote hu ∈ bE (R, C) the function defined by hu : z → exp(iu z).

From the result of Lemma 3.7 we infer that

E x (exp (iu Xt)) = Tthu(x) = α ψ(t; x, u) + (1 -α) ψ(t; x, -u) + (1 -α) ψ(t; x, u) -ψ(t; x, -u) (38)
Similarly, from the result of Lemma 3.6, we have

E mn exp iu n -1/2 S ⌊nt⌋ = T 0 ⌊nt⌋ h n -1/2 u (mn)
= α ψ0(⌊nt⌋; mn, n -1/2 u) + (1 -α) ψ0(⌊nt⌋; mn, -n -1/2 u)

+ (1 -α) ψ0(⌊nt⌋; mn, n -1/2 u) -ψ0(⌊nt⌋; mn, -n -1/2 u) . (39)

Classical results regarding the convergence of characteristic functions for normalized symmetrized random walks reflected at 0 (resp. killed at 0) towards the characteristic functions of the standard Brownian motion reflected at 0 (resp. killed at 0) ensure that ψ0 ⌊nt⌋; mn, n -1/2 u (resp. ψ0 ⌊nt⌋; mn, n -1/2 u ) converges simply to ψ(t; x, u) (resp. to ψ(t; x, u)) as n tends to infinity.

By comparison of (39) with (38) and by the linearity of the convergence, we get

E mn exp iu n -1/2 S ⌊nt⌋ -----→ n→+∞ E x [exp (iuXt)] .
Similar arguments show that this convergence also holds for x < 0 and x = 0. Recasting this convergence using our operators writes: for any x ∈ R,

T 0 ⌊nt⌋ h n -1/2 u (mn) -----→ n→+∞ Tthu(x). (40)
We leave it to the reader to check using very close arguments that we have also

T 0 ⌊nt⌋+1 h n -1/2 u (mn) -----→ n→+∞ Tthu(x). (41)
Let us now turn to the convergence of the characteristic function for joint pseudo distributions. As above, let x ≥ 0, 0 ≤ t1 < t2, u1, u2 ∈ R and set mn = ⌊ √ nx⌋. Applying the Markov property of the pseudo random walk and the result of Lemma 3.6 gives Ψ0(⌊nt1⌋, ⌊nt2⌋; mn, u1, u2)

:= E mn exp iu1 n -1/2 S ⌊nt 1 ⌋ + iu2 n -1/2 S ⌊nt 2 ⌋ = E mn exp iu1 n -1/2 S ⌊nt 1 ⌋ E S ⌊nt 1 ⌋ exp iu2 n -1/2 S ⌊nt 2 ⌋-⌊nt 1 ⌋ = T 0 ⌊nt 1 ⌋ j → h n -1/2 u 1 (j)T 0 ⌊nt 2 ⌋-⌊nt 1 ⌋ h n -1/2 u 2 (j) (mn). ( 42
)
Note that (⌊nt2⌋ -⌊nt1⌋) ∈ {⌊n(t2 -t1)⌋, ⌊n(t2 -t1)⌋ + 1}. As before, we are in position to use standard convergence results regarding the characteristic functions of normalized classical symmetric random walks.

Using (40) and (41), we prove that Ψ0(⌊nt1⌋, ⌊nt2⌋; mn, u1, u2) converges as n tends to infinity to Ψ(t1, t2; x, u1, u2) := Tt 1 [y → hu 1 (y)Tt 2 -t 1 hu 2 (y)] (x).

Note that from the definitions of the operators (T 0 j ) j∈N and (Tt) t≥0 this convergence results from the properties of the weak convergence of classical symmetric random walks (reflected or killed) towards their corresponding standard Brownian motion (reflected or killed).

It remains to check that

Ψ(t1, t2; x, u1, u2) = E x [exp (i(u1Xt 1 + u2Xt 2 )] (43)
in order to conclude our proof. But, exactly as before this last equality is easily verified by the use of the Markov property for the pseudo-process (Xt) at time t1 (which itself derives directly from the Chapman-Kolmogorov identity). Indeed, by application of Lemma 3.7 and the definition of hu, we have

Ψ(t1, t2; x, u1, u2) := Tt 1 [y → hu 1 (y)Tt 2 -t 1 hu 2 (y)] (x) = E x exp (iu1 Xt 1 ) E X t 1 (exp (iu2 Xt 2 -t 1 )) = E x [exp (i(u1Xt 1 + u2Xt 2 )]
which gives (43). Hence,

lim n→+∞ E mn exp i u1 n -1/2 S ⌊nt 1 ⌋ + u2 n -1/2 S ⌊nt 2 ⌋ = E x [exp (i (u1 Xt 1 + u2 Xt 2 ))] .
Finally, we may extend the previous limit result by induction to prove that for any ℓ ∈ N * , u1, . . . , u ℓ ∈ R and times 0

≤ t1 < • • • < t ℓ , lim n→+∞ E mn exp i ℓ j=1 uj n -1/2 S ⌊nt j ⌋ = E x exp i ℓ j=1 ujXt j .
The proof is completed.

The evolution equation on a finite interval

In this section we assume 0 < a < ∞. Again, we first address the case k > 0, in which we compute the transition function of the skew Brownian motion killed at the end-points of I = (-a, a), then we treat the general case k ∈ R * \{-1}.

4.1

The transition function of the SBM killed at -a or a (case k > 0)

In the same manner as for [START_REF] Raviart | Introduction à l'analyse numérique des équations aux dérivées partielles[END_REF], it can be shown that the solution to the parabolic problem ( 2) on (-a, a) × R can be represented as

u(t, x) = E x [u0( Xt)] ( 44
)
where X is the SBM of parameter β = 1-k 1+k , killed at -a or a. This process behaves like the SBM X as long as it does not exit from (-a, a). When it touches -a or a it is sent at a cemetery point ∂. By convention, for any function f one has f (∂) = 0, which ensures that the homogeneous Dirichlet boundary condition in ( 4) is satisfied.

Let us assume that x ∈ (-a, a) and let Ť(-a,a) = inf{t ≥ 0 : Xt / ∈ (-a, a)}. We compute

E x [u0( Xt)] = E x [u0( Xt); Ť(-a,a) ≤ t] + E x [u0( Xt); Ť(-a,a) > t] = E x [u0(∂); Ť(-a,a) ≤ t] + E x [u0( Xt); Ť(-a,a) > t] = E x [u0(Xt); T (-a,a) > t], ( 45
)
where X is the SBM considered in Section 3.2 and T (-a,a) = inf{t ≥ 0 : Xt / ∈ (-a, a)}. We first derive the expression of the kernel p(t, x, y) that satisfies P x (Xt ∈ dy; T (-a,a) > t) = p(t, x, y)dy, i.e. the fundamental solution to (2) : Proposition 4.1. Let β ∈ (-1, 1). Let X be the solution of [START_REF] Revuz | Continuous martingales and Brownian motion[END_REF], i.e. the SBM of parameter β. Let a > 0.

For any x, y ∈ (-a, a), any t > 0, one has P x (Xt ∈ dy; T (-a,a) > t) = p(t, x, y)dy with

p(t, x, y) =                      (1 -β)p (-a,a) W (t, x, y) if x ≥ 0, y < 0 -βp (0,a) W (t, x, y) + (1 + β)p (-a,a) W (t, x, y) if x > 0, y > 0 (1 + β)p (-a,a) W (t, x, y) if x ≤ 0, y > 0 βp (-a,0) W (t, x, y) + (1 -β)p (-a,a) W (t, x, y) if x < 0, y < 0. ( 46
)
Proof of Proposition 4.1. We use Lemma 3.1 and take advantage of the symmetry of the space interval (-a, a): then it is possible to adapt the arguments in [START_REF] Walsh | A diffusion with a discontinuous local time[END_REF], using for example in particular the transition function of the Brownian motion killed at -a or a instead of the transition function of the Brownian motion. We detail all the steps for the sake of completeness.

We set T0 = inf{t ≥ 0 : Xt = 0} and remark that T (-a,a) = inf{t ≥ 0 : Xt / ∈ (-a, a)} = inf{t ≥ 0 : |Xt| = a}. We denote τ0 = inf{t ≥ 0 : Wt = 0} and τ (-a,a) = inf{t ≥ 0 : |Wt| = a}. In the computations below, by a slight abuse of notation we may denote by P x either P(•|X0 = x) or P(•|W0 = x). This will be clear from the context. Let t > 0. We first treat the case x > 0, y < 0. Using the fact that P x (Xt ∈ dy; T0 > t) = 0, and using the strong Markov property of X, it follows that P x (Xt ∈ dy; T (-a,a) > t) = P x (Xt ∈ dy; T0 ≤ t; T (-a,a) > t) + P x (Xt ∈ dy; T0 > t; T (-a,a) > t)

= E x 1 T 0 ≤t P x (Xt ∈ dy; T (-a,a) > t|FT 0 ) = t 0 P 0 ( Xt-u ∈ dy; T(-a,a) > t -u)f x T 0 (u)du,
where X is another SBM of parameter β, starting from zero under P 0 , T(-a,a) = inf{t ≥ 0 : | Xt| = a} and f x T 0 (u)du is the law of T0 under P x . From Lemma 3.1 we have

P 0 ( Xt-u ∈ dy; T(-a,a) > t -u) = P 0 (| Xt-u| ∈ -dy; T(-a,a) > t -u; Xt-u < 0) = 1-β 2 P 0 (| Wt-u| ∈ -dy; τ(-a,a) > t -u) = 1-β 2 P 0 ( Wt-u ∈ -dy; τ(-a,a) > t -u) + P 0 (-Wt-u ∈ -dy; τ(-a,a) > t -u) = (1 -β)P 0 ( Wt-u ∈ dy; τ(-a,a) > t -u)
where W is a Brownian motion starting from zero under P 0 and τ(-a,a) = inf{t ≥ 0 : | Wt| = a}. Thus denoting by f x τ 0 (u)du the law of τ0 under P x , noticing that f x T 0 (u) = f x τ 0 (u), and using this time the strong Markov property of W , we calculate

P x (Xt ∈ dy; T (-a,a) > t) = (1 -β) t 0 P 0 ( Wt-u ∈ dy; τ(-a,a) > t -u)f x τ 0 (u)du, = (1 -β)P x (Wt ∈ dy; τ0 ≤ t; τ (-a,a) > t) = (1 -β)P x (Wt ∈ dy; τ (-a,a) > t) = (1 -β)p (-a,a) W (t, x, y)dy,
and obtain the first line of (46). Next, we treat the case x > 0 and y > 0. We have P x (Xt ∈ dy; T (-a,a) > t)

= P x (Xt ∈ dy; T0 ≤ t; T (-a,a) > t) + P x (Xt ∈ dy; T0 > t; T (-a,a) > t)

As 0 < x, y < a, the term P x (Xt ∈ dy; T0 > t; T (-a,a) > t) corresponds to the transition of a Brownian motion killed at 0 or a, or in other words,

P x (Xt ∈ dy; T0 > t; T (-a,a) > t) = p (0,a) W (t, x, y)dy.
As for the term P x (Xt ∈ dy; T0 ≤ t; T (-a,a) > t), we use similar computations as in the case x > 0, y < 0. Since Xt is positive this time, we obtain P x (Xt ∈ dy; T0 ≤ t; T (-a,a) > t) = (1 + β)P x (Wt ∈ dy; τ0 ≤ t; τ (-a,a) > t).

Notice also that P x (Wt ∈ dy; τ0 ≤ t; τ (-a,a) > t)

= P x (Wt ∈ dy; τ (-a,a) > t) -P x (Wt ∈ dy; τ0 > t; τ (-a,a) > t)

= p (-a,a) W (t, x, y)dy -p (0,a) W (t, x, y)dy.

Putting all the pieces together yields

P x (Xt ∈ dy; T (-a,a) > t) = p (0,a) W (t, x, y) + (1 + β){p (-a,a) W (t, x, y) -p (0,a) W (t, x, y)} dy = -βp (0,a) W (t, x, y) + (1 + β)p (-a,a) W (t, x, y) dy,
and thus, the second line of (46). The remaining cases can be treated in a similar manner.

Remark 4.2. Note the consistence of ( 21) with (46), as a → ∞ in (46).

Remark 4.3. Note that as a -a p(t, x, y)dy = P x ( T (-a,a) > t), p(t, x, •) does not necessarily integrate to 1. The definition of a genuine family of transition probability measures from p(t, x, y)dy would require combining this density measure with a Dirac measure that charges ∂ with probability P x ( T (-a,a) ≤ t), see [START_REF] Revuz | Continuous martingales and Brownian motion[END_REF], p84.

Still by a slight abuse of language we call this kernel the transition function of X. Note that in the hereafter examined case of a negative coefficient, the integrals of p(t, x, •) do not even necessarily correspond to true probabilities.

The fundamental solution to

(2) in the general case k ∈ R * \ {-1}
When k > 0, the transition function yields a fundamental solution to the evolution equation. When k < 0, k ̸ = -1, as in the case when I = R, we check that its expression also yields a fundamental solution.

Lemma 4.4. Let k ∈ R * \ {-1} and define β = 1-k 1+k . The kernel p(t, x, y) defined by (46) is the fundamental solution to (3)-( 4), or to [START_REF]T-coercivity for the Maxwell problem with sign-changing coefficients[END_REF]. Equivalently it satisfies for any y ∈ (-a, a),

∂t p(t, x, y) = 1 2 ∂ 2 xx p(t, x, y), ∀x ∈ (-a, 0) ∪ (0, a), (47) and p(t, 0-, y) = p(t, 0+, y), (48)
k∂x p(t, 0-, y) = ∂x p(t, 0+, y), Proof. If a kernel q(t, x, y) is a fundamental solution to (3)-( 4), then by definition, the function

(t, x) → u(t, x) = a -a u0(y)q(t, x, y)dy (51) solves (3)
-(4) for any u0 ∈ C ∞ c (-a, a). As ∂tu(t, x) = a -a u0 (y)∂tq(t, x, y)dy, ∂xu(t, x) = a -a u0(y)∂xq(t, x, y)dy,
and

∂ 2 xx u(t, x) = a -a u0(y)∂ 2 xx q(t, x, y)dy, letting u0 vary in C ∞ c (-a, a)
, it is easy to see that q(t, x, y) satisfies ( 47)-(50). Conversely, if (47)-(50) hold for q(t, x, y), one may check that the function defined by (51) solves ( 3)-( 4). Thus the equivalence.

Consider first k > 0 and β = 1-k 1+k . From (44), (45) and Proposition 4.1, it is clear that the kernel p(t, x, y) defined by ( 46) is the fundamental solution to (3)-( 4), and therefore satisfies (47)-(50).

Considering (47) and for example the first line of (46) it is clear that for x < 0 and y > 0 one has

∂tp (-a,a) W (t, x, y) = 1 2 ∂ 2 xx p (-a,a) W (t, x, y) (52)
(of course this fact can formally be seen from [START_REF] Martinez | One-dimensional parabolic diffraction equations: pointwise estimates and discretization of related stochastic differential equations with weighted local times[END_REF]; see also Remark 4.5 below).

Consider now k ∈ R * -\ {-1}, set β = 1-k 1+k
, and define the kernel p(t, x, y) by ( 46). From the first line of (46), and (52) it is clear that (47) holds for x < 0 and y > 0. The other cases can be treated similarly, so that p(t, x, y) solves (47) for any y ∈ (-a, a), x ∈ (-a, 0) ∪ (0, a).

Concerning the transmission conditions, we first note that when y > 0, Eq. ( 49) can be rewritten

(1 + β)∂x p(t, 0+, y) = (1 -β)∂x p(t, 0-, y) (53)
or equivalently, using the definition of p(t, x, y) in (46) (with y < 0)

(1 -β 2 )∂xp (-a,a) W (t, 0-, y) = (-β -β 2 )∂xp (0,a) W (t, 0+, y) + (1 + 2β + β 2 )∂xp (-a,a) W (t, 0+, y),
which in turn, due to the continuity of ∂xp (-a,a) W (t, 0, y) at x = 0, is equivalent to

∂xp (0,a) W (t, 0+, y) = 2∂xp (-a,a) W (t, 0, y). (54)
Note that condition (54) only bears on the properties of p (-a,a) W and p (0,a) W

, and is equivalent to (49) whatever the sign of k, provided k ̸ = -1. Thus we know from the case k > 0 that (54) holds, and thus (49) must hold for all k ̸ = -1. The same argument applies to y < 0. Equation ( 48) is a consequence of ( 17) and of the continuity at x = 0 of the kernel p (-a,a) W (t, x, y). Finally, the Dirichlet condition (50) follows from [START_REF] Nguyen | Limiting absorption principle and well-posedness for the Helmholtz equation with sign changing coefficients[END_REF] and the form of p(t, x, y). We conclude that the latter is indeed the fundamental solution of ( 3)-( 4) when k < 0, k ̸ = -1.

Remark 4.5. In fact it is possible to start from (46), and to use [START_REF] Martinez | One-dimensional parabolic diffraction equations: pointwise estimates and discretization of related stochastic differential equations with weighted local times[END_REF] in order to formally check (47)-(49), for any k ∈ R * \ {-1}. We especially want to explain how one can derive the transmission condition in (49).

As we have already seen in the proof of Lemma 4.4, it is enough to check (54). On one hand one has ∂xp

(0,a) W (t, 0+, y) = ∞ n=-∞ y -2na t g(t, 0, y -2na) - -y -2na t g(t, 0, -y -2na) = ∞ n=-∞ y -2na t g(t, 0, y -2na) + y + 2na t g(t, 0, y + 2na) = 2 ∞ n=-∞ y + 2na t g(t, 0, y + 2na)
On the other hand one has

∂xp (-a,a) W (t, 0, y) = ∞ n=-∞ y -4na t g(t, 0, y -4na) - 2a -y -4na t g(t, 0, 2a -y -4na) = ∞ n=-∞ y -4na t g(t, 0, y -4na) + y + 4na -2a t g(t, 0, y + 4na -2a) = ∞ n=-∞ y -2(2n)a t g(t, 0, y -2(2n)a) + y + 2(2n -1)a t g(t, 0, y + 2(2n -1)a) = ∞ n=-∞ y + 2na t g(t, 0, y + 2na).
Therefore (54).

Note that in the case k ∈ R * -\ {-1} the kernel p(t, x, y) defined by ( 46) is not positive: indeed in that case β > 1 and it suffices to examine the first line of (46).

So even if we complement the density measure p(t, x, y)dy with Dirac measures in order to obtain a family of transition pseudo probability measures (see our Remark 4.3), the latter will not define a Markov process. One could however introduce the pseudo-process associated to this family, in the spirit of Section 3.2. Remark 4.6. We end up this section by noticing that as expected, the form of the kernel (46) coincides with [START_REF] Lachal | From pseudorandom walk to pseudo-brownian motion: First exit time from a one-sided or a two-sided interval[END_REF]. Indeed, the transition density of the Brownian motion on an interval (a, b), killed at a or b has the following spectral representation (see Appendix 1 in [START_REF] Nikolaevitch | Handbook of brownian motion : facts and formulae, Probability and its applications[END_REF])

p (a,b) W (t, x, y) = 2 b -a n≥1 exp - n 2 π 2 2(b -a) 2 t sin nπ b -a (x -a) sin nπ b -a (y -a) .
Thus one has

p (-a,a) W (t, x, y) = 1 a n≥1 exp - n 2 π 2 8a 2 t sin nπ 2a (x -a) sin nπ 2a (y -a) .
Regrouping the terms with odd and even indices, the above expression rewrites, when x ≥ 0 and y < 0

p (-a,a) W (t, x, y) = + 1 a q≥1 exp - q 2 π 2 2a 2 t sin qπ a (x -a) sin qπ a (y -a) + 1 a q≥1 exp - (2q -1) 2 π 2 8a 2 t cos (2q -1)π 2a (x -a) cos (2q -1)π 2a (y -a) = 1 a q≥1 1 k e -µ 2
q t/2 gq(x)gq(y) + e -λ 2 q t/2 fq(x)fq(y)

Multiplying by 1 -β = 2k k+1 to retrieve the first line of (46), one easily recovers the expression (13) when x ≥ 0 and y < 0. A similar calculation applies for the other possible choices of x and y.

5 Various numerical schemes for the PDE [START_REF]T-coercivity for the Maxwell problem with sign-changing coefficients[END_REF] In this section we want to construct several numerical schemes for the approximation of the solution u of (2), inspired by the theoretical results of the previous sections.

In Section 5.1 we construct a scheme ūN spec inspired by the spectral representation of the semigroup of Section 2.

In Section 5.2 we explain how we can infer a finite difference type scheme ūn RW from the scaled pseudo asymmetric random walk X n of Section 3.3.

In Section 5.3 we construct a scheme ūh,N f und , which is inspired by the fact that the fundamental solution p(t, x, y) computed in Section 4.1 involves transition functions of killed Brownian motions, which can be seen as the fundamental solutions of simple heat equations, with homogeneous Dirichlet boundary conditions.

An initial condition u0 is given (in Section 5.1 it is of class L 2 (I), in Sections 5.2 and 5.3 we can imagine it is continuous and bounded).

We have k ∈ R * \ {-1} and set β = 1-k 1+k and α = 1+β 2 . In the case k ∈ (0, ∞) several numerical schemes are available (one can for example first perform a finite element discretization w.r.t. the space variable, and then a Crank-Nicholson scheme, e.g. [START_REF] Raviart | Introduction à l'analyse numérique des équations aux dérivées partielles[END_REF]), including probabilistic ones (e.g. [START_REF] Étoré | On random walk simulation of one-dimensional diffusion processes with discontinuous coefficients[END_REF], [START_REF] Martinez | One-dimensional parabolic diffraction equations: pointwise estimates and discretization of related stochastic differential equations with weighted local times[END_REF], [START_REF] Étoré | Exact simulation of one-dimensional stochastic differential equations involving the local time at zero of the unknown process[END_REF]). Thus, for this well explored case the schemes presented hereafter may provide additional methods (up to our knowledge they have never been proposed in this form).

Their main interest is that they allow to handle the case k ∈ R * -\ {-1}, for which in particular classical probabilistic methods (e.g. [START_REF] Étoré | On random walk simulation of one-dimensional diffusion processes with discontinuous coefficients[END_REF], [START_REF] Martinez | One-dimensional parabolic diffraction equations: pointwise estimates and discretization of related stochastic differential equations with weighted local times[END_REF], [START_REF] Étoré | Exact simulation of one-dimensional stochastic differential equations involving the local time at zero of the unknown process[END_REF]) cannot be applied. In the case k > 0 the latter deeply rely on stochastic simulations of the trajectories of the SBM of parameter β ∈ (-1, 1). When k ∈ R * -\ {-1} one can define a pseudo-SBM, but it proves difficult to define associated trajectories. Note also that the techniques for proving the convergence of deterministic schemes do not apply in the case k ∈ R * -\ {-1}.

Scheme inspired by the spectral representation of the semigroup

We assume that u0 ∈ L 2 (I) and fix a truncation order N ∈ N * . Then ūN spec is defined by

ūN spec (t, x) = N n=1
ane -λ 2 n t/2 fn(x) + bne -µ 2 n t/2 gn(x)

where the λn's, µn's, gn's and fn's are as in Section 2, and the an's and bn's are defined as in [START_REF] Bonnetier | Spectral bases for graphs mixing dielectric and metamaterials[END_REF]. That is to say ūN spec (t, x) is obtained by keeping the first N terms in the spectral representation (11) of u(t, x).

Of course in order to use this method we have to be able to compute the an's and bn's. In Section 6 we consider two examples where these computations can be done explicitly. If this is not the case we can resort to numerical integration. We sum up hereafter the proposed algorithm (Algorithm 1).

ALGORITHM 1: Computation of ūN spec (t, x). Parameters of the method: A time 0 ≤ t < ∞ at which we want to compute the approached solution. A truncation order N ∈ N * . Remember that the λn's, µn's, gn's and fn's are as in Section 2. (either exactly or by numerical integration).

2) Return

ūN spec (t, x) = N n=1
ane -λ 2 n t/2 fn(x) + bne -µ 2 n t/2 gn(x).

As the series in ( 11) is normally convergent we get immediately the following convergence result.

Proposition 5.1. Let us consider u the solution of (2). Let t ≥ 0 and let us consider for any N ∈ N * the function ūN spec (t, •) defined by Algorithm 1. We have

||u(t, •) -ūN spec (t, •)|| L 2 (I) ----→ N →∞ 0.

Scheme inspired by the scaled pseudo asymmetric random walk

Let u0 : R → R be continuous and bounded and consider first [START_REF] Gall | One-dimensional stochastic differential equations involving the local times of the unknown process, Stochastic analysis and applications[END_REF]. Let us fix n ∈ N * a discretization order. From Proposition 3.8, u(T, x) = E x [u0(XT )] is approached by E x [u0( X n T )] = E x [u0(n -1/2 S ⌊nT ⌋ )]. We recall that the expectation symbols have to be understood as pseudo expectations, and that X is the pseudo SBM of Section 3.2.

Let T > 0 and let us assume that nT = N , an integer. Denoting u n 0 the function defined by

u n 0 (z) = u0(n -1/2 z), z ∈ Z, one has E x [u0( X n T )] ≈ E ⌊ √ nx⌋ [u n 0 (SN )]
where the expectation in the right hand side is computed under

P ⌊ √ nx⌋ s.t. P ⌊ √ nx⌋ (S0 = ⌊ √ nx⌋) = 1. Let us denote vN = u n 0 .
Using the (pseudo) Markov property of the (pseudo) random walk S we get

E ⌊ √ nx⌋ [u n 0 (SN )] = E ⌊ √ nx⌋ E[ vN (SN ) | S0, . . . , SN-1 ] = E ⌊ √ nx⌋ vN-1(SN-1) (55)
where we have denoted vN-1

(z) = E[ vN (SN ) | SN-1 = z ], z ∈ Z. In fact, defining more generally, vm-1(z) = E[ vm(Sm) | Sm-1 = z ], z ∈ Z, 1 ≤ m ≤ N,
and proceeding to further conditionings in (55) we get

E ⌊ √ nx⌋ [u n 0 (SN )] = E ⌊ √ nx⌋ vN-2(SN-2) = . . . = E ⌊ √ nx⌋ v0(S0) = v0(⌊ √ nx⌋).
Note that from the (possibly pseudo) transition probabilities of the random walk S we have for any

1 ≤ m ≤ N , vm-1(z) = 1 2 [vm(z + 1) + vm(z -1)]1 z̸ =0 + [αvm(z + 1) + (1 -α)vm(z -1)]1z=0, z ∈ Z.
ALGORITHM 2: Computation of ūn RW (t, x). Parameters of the method: A time horizon 0 < T < ∞ and a discretization order n ∈ N * . We set N = nT and Na = √ n a and assume this quantities are integers. We set xj = j/ √ n for -Na ≤ j ≤ Na.

Algorithm: 1) Set U 0 j = u0(xj) for any -Na + 1 ≤ j ≤ Na -1. 2) For 0 ≤ m ≤ N -1, compute U m+1 j = 1 2 U m j+1 + 1 2 U m j-1 for -Na + 1 ≤ j ≤ Na -1 U m+1 0 = αU m+1 1 + (1 -α)U m+1 -1
with the convention that U m -Na = U m Na = 0.

3) Return a piecewise constant function ūn

RW (t, x) satisfying ūn RW ( m n , xj) = U m j , ∀ -Na + 1 ≤ j ≤ Na -1, ∀0 ≤ m ≤ M, and ūn RW ( m n , ±a) = 0, ∀1 ≤ m ≤ M.
It should be possible to adapt the results of Proposition 3.8 to prove convergence of the above scheme. This would require considering the trajectories of killed scaled pseudo asymmetric random walks, which presents difficulties we have decided not to address in the present paper.

Nevertheless, we suspect that the function ūn RW : [0, T ] × [-a, a] → R defined by Algorithm 2 should converge towards u the solution of (2): sup

(t,x)∈[0,T ]×[-a,a] |u(t, x) -ūn RW (t, x)| ----→ n→∞ 0. (62)
This will be illustrated by the numerical experiments of Section 6.

Scheme inspired by the expression of the fundamental solution involving the transition function of the killed Brownian motion

Let u0 ∈ C(I; R). We denote u + 0 = u01 R + and u - 0 = u01 R * -. Let x ∈ (0, a), from Lemma 4.4 and Eq. (46) we have

u(t, x) = a -a u0(y)p(t, x, y)dy = (1 -β) a -a u - 0 (y)p (-a,a) W (t, x, y)dy + (1 + β) a -a u + 0 (y)p (-a,a) W (t, x, y)dy -β a 0 u + 0 (y)p (0,a) W (t, x, y)dy = (1 -β)u1(t, x) + (1 + β)u2(t, x) -βu3,+(t, x), (63)
where the functions u1, u2 and u3,+ are respectively solution of the following heat equations:

         ∂tu1(t, x) = 1 2 ∂ 2 xx u1(t, x), x ∈ I, t > 0, u1(0, x) = u - 0 (x) x ∈ I, u1(t, ±a) = 0, t > 0,          ∂tu2(t, x) = 1 2 ∂ 2 xx u2(t, x), x ∈ I, t > 0, u2(0, x) = u + 0 (x) x ∈ I, u2(t, ±a) = 0, t > 0, and
         ∂tu3,+(t, x) = 1 2 ∂ 2 xx u3,+(t, x), x ∈ (0, a), t > 0, u3,+(0, x) = u + 0 (x) x ∈ (0, a), u3,+(t, 0) = u3,+(t, a) = 0, t > 0.
Indeed p (-a,a) W (t, x, y) (resp. p (0,a) W (t, x, y)) may be viewed as the fundamental solution of the heat equation with half Laplacian on the domain (-a, a) (resp. (0, a)), with homogeneous Dirichlet boundary conditions ( [START_REF] Nikolaevitch | Handbook of brownian motion : facts and formulae, Probability and its applications[END_REF], Appendix I, N� 6).

In the same manner, for x < 0 we have

u(t, x) = (1 -β)u1(t, x) + (1 + β)u2(t, x) + βu3,-(t, x), (64)
with u1 and u2 as before and u3,-the solution of

         ∂tu3,-(t, x) = 1 2 ∂ 2 xx u3,-(t, x), x ∈ (-a, 0), t > 0, u3,-(0, x) = u - 0 (x) x ∈ (-a, 0), u3,-(t, -a) = u3,-(t, 0) = 0, t > 0.
Our idea is to perform finite different schemes for u1, u2, u3,± and to combine them through (63)(64) in order to get a scheme for the approximation of u. We sum up the procedure in Algorithm 3 where we use implicit finite different schemes, which are known to be unconditionably stable ( [START_REF] Allaire | Analyse numérique et optimisation: une introduction à la modélisation mathématique et à la simulation numérique[END_REF]). In the present case they are also consistent and therefore convergent by Lax principle.

ALGORITHM 3: Computation of ūh,N f und (t, x). Parameters of the method: A time horizon 0 < T < ∞. A time dicretization order N is given and we set δt = T /N . A space step h is given and we set Na = a/h (we assume this is an integer). We set xj = jh for -Na ≤ j ≤ Na.

Algorithm: 1) Set U 0 1,j = u - 0 (xj) for any -Na + 1 ≤ j ≤ Na -1. Set U 0 2,j = u + 0 (xj) for any -Na + 1 ≤ j ≤ Na -1. Set U 0 3+,j = u + 0 (xj) for any 1 ≤ j ≤ Na -1. Set U 0 3-,j = u - 0 (xj) for any -Na + 1 ≤ j ≤ -1. 2) For 0 ≤ m ≤ N -1, compute the vectors U m+1 From the convergence of the finite difference schemes we immediately get the following convergence result. Remark 5.3. In fact we could infer several other numerical schemes from Eq. (46). For example, consider x > 0 and t > 0. Eq. ( 63) also implies that u(t, x) = (1 -β)E x [u - 0 (Wt); t < τ (-a,a) ] + (1 + β)E x [u + 0 (Wt); t < τ (-a,a) ] -βE x [u + 0 (Wt); t < τ (0,a) ].

So we could consider approaching each of the above expectations by Monte Carlo sums involving samples of independent Brownian motions and the corresponding stopping times τ (-a,a) or τ (0,a) (see e.g. [START_REF] Karoui | Les outils stochastiques des marchés financiers : une visite guidée de einstein à black-scholes[END_REF]). However, as the space dimension is one, we know that this Monte Carlo method would be slower than the finite difference approach described in Algorithm 3. Nevertheless such an approach could be interesting if we would address the problem in a space of higher dimension.

Numerical experiments

Example 1. We take I = (-1, 1) (i.e. a = 1) and k = -0.5. We choose the following initial condition u0(x) = 10x 3 -3x 2 -9x + 4 2 , ∀x ∈ I.

Indeed in order to use Algorithm 1 we have to compute the an's and bn's, through Eq. [START_REF] Bonnetier | Spectral bases for graphs mixing dielectric and metamaterials[END_REF]. By the polynomial nature of the initial condition u0 these coefficients will be made explicit, providing thus a benchmark for the finite difference scheme inspired algorithms (Algo. 2 and 3). Remember that the quantites I A(x)|fn(x)| 2 dx and I A(x)|gn(x)| 2 dx in Eq. ( 9) are given by Eq. [START_REF] Karatzas | Brownian motion and stochastic calculus[END_REF]. Besides one can compute A(x)gn(x)u0(x) dx = -k (-1) n π 2 n 2 -60 π 3 n 3 .

Recalling that here λn = (2n-1)π 2 and µn = nπ, n ≥ 1, we have everything at hand to perform Algorithm 1.

For performing Algorithms 2 and 3 no previous computation is needed. We see a very good concordance between the three methods (which actually can be observed for coarser discretizations).

In particular we can numerically check the convergence of the Algorithm 2 announced in Eq. (62). To that aim we consider ū200 spec (T, •) as the reference solution and report in (t, •), at times t = 0.012, t = 0.12 and t = 0.4, for the initial condition u 0 (x) = 1 x<0 -0.5.

Algorithm: 1 )

 1 Compute an= I A(x)u0(x)fn(x) dx I A(x)|fn(x)| 2 dx , bn = I A(x)u0(x)gn(x) dx I A(x)|gn(x)| 2 dx

2 2 3+ 2 3 )

 2223 for -Na + 1 ≤ j ≤ Nafor -Na + 1 ≤ j ≤ Na -1 U m+1 3+,j -U m 3+,j δt = U m+1 3+,j+1 -2U m+1 3+,j + U m+1 for -Na + 1 ≤ j ≤ -1 with the conventions that U m+1 1,±Na = U m+1 2,±Na = U m+1 3+,Na = U m+1 3+,0 = U m+1 3-,-Na = U m+1 3-,0 = 0. Return a piecewise constant function ūh,N f und (t, x) satisfying ūh,N f und (m δt, xj) = (1 -β)U m 1,j + (1 + β)U m 2,j -βU m 3+,j , ∀1 ≤ j ≤ Na -1, ∀0 ≤ m ≤ M, ūh,N f und (m δt, xj) = (1 -β)U m 1,j + (1 + β)U m 2,j + βU m 3-,j , ∀ -Na + 1 ≤ j ≤ -1, ∀0 ≤ m ≤ M, ūh,N f und (m δt, 0) = (1 -β)U m 1,0 + (1 + β)U m2,0 , ∀0 ≤ m ≤ M, and ūh,N f und (m δt, ±a) = 0, ∀1 ≤ m ≤ M.

Proposition 5 . 2 .

 52 Let us consider u the solution of (2). Let 0 < T < 0 and let us consider for any N ∈ N * and any h ∈ (0, a) the function ūh,N f und : [0, T ] × [-a, a] → R defined by Algorithm 3. We have sup(t,x)∈[0,T ]×[-a,a] |u(t, x) -ūh,N f und (t, x)| -------→ h↓0, N →∞ 0.

1 - 1 A 2 + 1 ( 2 n 2 and 1 - 1

 11212211 (x)fn(x)u0(x) dx = k (2 n -1) 4 π 4 -528 (n -1/2) π (-1) n -480 -72 (n -1/2) 2 π -1) 4 π 4 -528 π 1/33 (n -1/2) 2 π 2 -9 11(n -1/2) (-1) n + 480 + 72 (n -1/2) 2 π

Figure 1

 1 Figure 1 represents the graphs of ū200 spec (T, •), ū2.5×10 5 RW (T, •) and ū2×10 -3 , 500 f und (T, •) at T = 0.4 with the following choices of parameters: N = 200 for Algorithm 1, n = 2.5 × 10 5 for Algorithm 2, a space step h = 2 × 10 -3 and a time discretization order N = 500 for Algorithm 3.We see a very good concordance between the three methods (which actually can be observed for coarser discretizations).In particular we can numerically check the convergence of the Algorithm 2 announced in Eq. (62). To that aim we consider ū200 spec (T, •) as the reference solution and report in Table1the value of

  ū200 spec (T, x) -ūn RW (T, x)for increasing values of n, up to n = 2.5 × 10 5 . Convergence is indeed observed.

Figure 1 : 5 RW

 15 Figure 1: Plot of an approximation of the function u(T = 0.4, •), by ū200 spec (T, •), ū2.5×10 5 RW

Figure 2 :

 2 Figure 2: Plot of the initial condition u 0 and of an approximation of the function u(t, •), by ū2×10 -3 , 500 f und(t, •), at times t = 8 × 10 -4 , t = 0.12 and t = 0.4, for the initial condition u 0 (x) = (10x 3 -3x 2 -9x + 4)/2.

Figure 3 : 5 RW

 35 Figure 3: Plot of an approximation of the function u(T = 0.4, •), by ū200 spec (T, •), ū2.5×10 5 RW

Figure 4 :

 4 Figure 4: Plot of the initial condition u 0 and of an approximation of the function u(t, •), by ū2×10 -3 , 500 f und

Table 1

 1

	the value of

Fundings

This research has been supported by the IRS projects "FKSC" and "Spectral properties of graphs with negative index materials" (for E. Bonnetier and P. Étoré) and by the Bézout Labex, funded by ANR, reference ANR-10-LABX-58 (for Miguel Martinez).

To sum up, one may approach E x [u0( X n T )] by computing v0(⌊ √ nx⌋) through the dynamical programming procedure

[vm(z + 1) + vm(z -1)]1 z̸ =0 ∀z ∈ Z, ∀1 ≤ m ≤ N.

The algorithm (56) is written in a recursive form. We now rewrite it in an iterative form -using also a new set of notations, in order to stress the fact that it is very similar to an explicit finite difference scheme, with space step h = n -1/2 and time step δt = n -1 . Let us consider the space grid {xj} j∈Z defined by xj = j/ √ n for any j ∈ Z, and the scheme {U m j }, for j ∈ Z, 0 ≤ m ≤ M , defined by

and, for 0 ≤ m ≤ N -1,

It is obvious that (57)-( 59) is equivalent to (56), in other words U M j = v0(j) for any j ∈ Z. Let us explain briefly why we may interpret (57)-(59) as an explicit finite difference scheme. The simplest way is to examine the case k = 1. Then α = 1 2 and (58)-(59) becomes

Besides ( 14) becomes simply the heat equation

Performing an explicit finite difference scheme with space step h and time step δt for Eq. ( 61) amounts to considering a space grid {x h j } j∈Z defined by x h j = jh for any j ∈ Z, and to compute

2h 2 for any j ∈ Z for any 0 ≤ m ≤ N -1. Taking h = n -1/2 and δt = n -1 (note that this corresponds to touching the bound giving the CFL condition, see e.g. [START_REF] Allaire | Analyse numérique et optimisation: une introduction à la modélisation mathématique et à la simulation numérique[END_REF]) we get (57) and (60). Therefore the interpretation.

In fact it seems that the transition (pseudo) probabilities of the random walk S suggests how to take into account the transmission condition in (15) in a finite difference scheme for ( 14), leading to condition (59).

Note that by applying the scheme (57)-(59) we get for any j ∈ Z, and any 0

For computational purposes we have to consider a PDE problem with a bounded space domain, and this is our PDE of interest [START_REF]T-coercivity for the Maxwell problem with sign-changing coefficients[END_REF].

Firstly, the domain [-a, a] is discretized with a grid {xj} Na j=-Na with Na = √ n a (we assume this quantity is an integer) and xj = j/ √ n for -Na ≤ j ≤ Na. spec (T, •) is considered as the benchmark).

On Figure 2 we check another interesting phenomenon. We plot the initial condition u0 and an approximation of u(t, •) by ūh,N f und (t, •) (we keep the same parameters h = 2 × 10 -3 and N = 500 as previously) at times t = 8 × 10 -4 , t = 0.12 and t = 0.4. Note that t = 8 × 10 -4 = T /N corresponds to the first time step in the finite difference scheme.

Observe that the transmission condition destroys the regularity of the initial datum. Also observe that the slope of the graph of u(t, •) is of negative sign at 0-, but twice bigger in absolute value than the positive slope at 0+. This is what we expect as k = -0.5.

Example 2. In this second example we keep I = (-1, 1) and k = -0.5, but choose the following initial condition u0(x) = 1x<0 -0.5, ∀x ∈ I.

(65)

Indeed we want to test the robustness of our numerical schemes to a non smooth initial condition, especially if this initial condition presents a discontinuity at the interface point x = 0. Of course for Algorithm 1 the initial condition can be taken in L 2 (I) so one knows that the algorithm converges if we take u0 defined by (65).

The analysis of probabilistic schemes such as Algorithm 2 usually relies on arguments of convergence in pseudo law, which are usually valid only for smooth functions. We suspect that the smoothing properties of the operator 1 2A ∇ • (A∇ ) could be used to prove convergence of Algorithm 2, even for non smooth initial conditions.

In order to use Algorithm 1 we have to compute the an's and bn's again, the µn's and λn's remaining unchanged. Easy computations show that an = -(-1) n 2(k -1) (2n -1)π(k + 1) and bn = 4 nπ(k + 1)

We take again an order of truncation N = 200 for Algorithm 1.

Again for the Algorithm 2 we take a discretization order n = 2.5 × 10 5 , and for the Algorithm 3 we take a space step h = 2 × 10 -3 and a time discretization order N = 500. (T, •) for T = 0.4 and again shows very good agreement between the results of the three schemes.

On Figure 4 we plot the initial condition u0 and an approximation of u(t, •) by ūh,N f und (t, •) (we keep the same parameters h = 2×10 -3 and N = 500 as previously) at times t = 0.012, t = 0.12 and t = 0.4. In particular, the plot illustrates the fact that the solution u(t, x) does not satisfy a maximum principle, one of the main arguments classically used in the numerical analysis of deterministic finite difference schemes for parabolic equations.

Data availability and conflict of interest statements

No datasets were generated or analysed during the current study. The authors declare no competing interests.