One-Dimensional Parabolic Diﬀraction Equations: Pointwise Estimates and Discretization of Related Stochastic Diﬀerential Equations With Weighted Local Times

In this paper we consider one-dimensional partial diﬀerential equations of parabolic type involving a divergence form operator with a discontinuous coeﬃcient and a transmission compatibility condition. We prove existence and uniqueness result by stochastic methods which also allow us to develop a low complexity Monte Carlo numerical resolution method. We get accurate pointwise estimates for the derivatives of the solution from which we get sharp convergence rate estimates for our stochastic numerical method. approximation method of this solution and to get sharp convergence rate estimates.


Introduction
Given a finite time horizon T and a positive matrix-valued function a(x) which is smooth except at the interface surfaces between subdomains of R d , consider the parabolic diffraction problem div(a(x)∇)u(t, x) = 0 for all (t, x) ∈ (0, T ] × R d , u(0, x) = f (x) for all x ∈ R d , Compatibility transmission conditions along the interfaces surfaces. (1) Suppose that 1 2 div(a(x)∇) is a strongly elliptic operator. Existence and uniqueness of continuous solutions with possibly discontinuous derivatives along the surfaces hold true: see, e.g. Ladyzenskaya et al. [12,chap.III,sec.13]. Our first objective is to provide a probabilistic interpretation of the solutions which allows us to get pointwise estimates for partial derivatives of the solution u(t, x). These estimates, which are interesting in their own, allow us to complete our second objective, that is, to develop an efficient stochastic numerical approximation method of this solution and to get sharp convergence rate estimates.
For reasons which we will describe soon, in this paper we complete our program when d = 1 only. Thus our results are first steps to address various applications where divergence form operators with a discontinuous coefficient arise and stochastic simulations are used: for example, the numerical resolution of solute transport equation in Geophysics (see, e.g., Salomon et al. [25] and references therein), the numerical resolution of the Poisson-Boltzmann equation in Molecular Dynamics (see, e.g., Bossy et al. [5] and references therein); another motivation comes from Neurosciences, more precisely from an algorithm of identification of the magnetic permittivity around the brain (see [6,7]). This algorithm actually solves an inverse problem: it consists in an iterative procedure aimed to compute the permittivity such that the solution of a Maxwell equation parametered by this permittivity fits with a good accuracy measurements obtained by sensors located on the patient's brain; this Maxwell equation depends on the values taken at the locations of the sensors by the solution of a Poisson equation involving a divergence form operator with discontinuous coefficients. Monte Carlo methods allow one to obtain this small set of values without solving the Poisson equation in its whole domain, which may significantly reduce the CPU time at each step of the iterative procedure.
Whatever is the dimension d, the theory of Dirichlet forms allows one to construct Markov processes whose generators in suitable Sobolev spaces are 1 2 div(a(x)∇) (see, e.g., the monography by Fukushima et al. [11]). However such a process constructed this way is expressed as the sum of a martingale and an abstract additive functional with finite quadratic variation; equivalently, it satisfies a Lyons-Zheng decomposition which involves its natural time reverse filtration and the logarithmic derivative of the (unknown) fundamental solution of (1) (see, e.g., Roskosz [23] for details). It thus seems difficult both to derive from these Markov processes, either poinwise estimates on partial derivatives of the function u(t, x), or to develop an efficient stochastic numerical resolution method for (1).
In the particular case of piecewise constant functions a(x), stochastic representations of u(t, x) can be obtained by transformations of Brownian motions and of the geometry of the transmission surfaces: see Bossy et al. [5].
For more general discontinuous functions a but in the one dimensional case d = 1, one can prove that 1 2 ∂ x (a(x)∂ x ) is the generator of the stochastic process solution of a stochastic differential equation (SDE) involving its own local time: see, e.g., Bass and Chen [2],Étoré [8], Lejay [14], Martinez [16]. This new description is the starting point for recent numerical studies: Lejay and Martinez [15] andÉtoré [8,9] proposed simulation methods for this solution based on approximations of a(x) and random walks simulations, and they analyzed the convergence rates of these methods. Here we propose a simpler numerical method and we interpret the strong solution to (1) in terms of the exact process.
To simplify the presentation, we now suppose that the function a(x) is discontinuous at point 0 only. See the section 8 for the case where a(x) has a finite number of discontinuities.
Thus, let a(x) = (σ(x)) 2 be a real function on R which is right continuous at point 0 and differentiable on R − {0} with a bounded derivative. Consider the one-dimensional stochastic differential equation with weighted local time Here L 0 t (X) is the right-sided local time corresponding to the sign function defined as sgn(x) := 1 for x > 0 and sgn(x) := −1 for x ≤ 0 (for properties of local times, see, e.g., Meyer [18]) and σ − is the left derivative of σ. Under conditions weaker than those of theorem 3.3 below, the equation (2) has a unique weak solution which is a strong Markov process : see Le Gall [13]. For all real number x 0 we thus may consider a probability space (Ω, F, P x 0 ), a one-dimensional standard Brownian motion (B t , t ≥ 0) on this space, and a solution X := (X t ) to (2) satisfying X 0 = x 0 , P x 0 − a.s. However, even simpler than Lyons-Zheng decompositions, this Markov process is not easy to simulate because of the difficulty to numerically approximate the local time process (L 0 t (X)). This leads us to apply a transformation which removes the local time of X (as Le Gall [13] did it to construct a solution to (2); Lejay [14] also used this transformation). We thus get a new stochastic differential equation without local time which can easily be discretized by the standard Euler scheme. As the transformation is one-to-one and its inverse is explicit, one then readily deduces an approximation X of X. Choosing X 0 = X 0 we then approximate u(t, x 0 ) by E x 0 f (X t ), the latter being computed by Monte Carlo simulations of X.
Our results are two-fold. First, we use probabilistic techniques to show that, for a wide class of functions f , the solution of the PDE (1) with d = 1 can be represented as and to get pointwise estimates for partial derivatives of this solution. Second, owing to these estimates, we prove a sharp convergence rate estimate for E x 0 f (X t ) to u(t, x 0 ). This convergence rate is unknown in the literature because the SDE obtained by removing the local time has discontinuous coefficients: whereas the convergence rates of discretizations of SDEs are well established when the coefficients are smooth (see a review in [27]) our estimates open the understanding of the discretization of SDEs with discontinuous coefficients. To our knowledge, the only results in that direction are due to Yan [28] who proves weak convergence of the Euler scheme for general SDEs with discontinuous coefficients but does not precise convergence rates. The paper is organized as follows. In section 2 we construct our transformed Euler discretization scheme for the SDE (2). In section 3, we state our main results. Our first results concern our stochastic representation of u(t, x) and pointwise estimates for its partial derivatives. They are respectively proven in sections (4) and (5). Our next results describe the convergence rate of our transformed Euler scheme: we distinguish the case where the initial function is flat around the discontinuity point 0 and the general case where this assumption is no longer true. The corresponding proofs are in sections (6) and (7). We discuss possible extensions of our results in section (8). In Appendix we remind technical results that we use in our proofs, namely, a representation of the density of the first passage time at 0 of an elliptic diffusion, and a recent estimate from [4] for the expected number of visits in small balls of Itô processes observed at discrete times.

Notation.
For all left continuous function g we denote either by g − (x) or by g(x−) the left limit of g at point x. When g is right continous, we denote either by g + (x) or by g(x+) the right limit of g at point x.
We denote by C b (R) the set of all bounded continuous functions with bounded continuous derivatives up to order .
In all the paper, for all integers 0 ≤ < ∞ and 1 ≤ p ≤ ∞ we denote the L p (R) norm of the function g by g p and we set where ∂ i x g is the i-th derivative of g. The positive real numbers denoted by C may vary from line to line; they only depend on the functions f and σ, the point x 0 , and the time horizon T . This means that, in particular, C does not depend on the discretization step h n of the Euler scheme and the smoothing parameter δ introduced in section 7. In addition, the quantifiers will make it clear when C does not depend on the function f .
The expectation E x 0 refers to the probability measure P x 0 under which X 0 = X 0 = x 0 a.s.

Our transformed Euler scheme
Suppose that a(0+) − a(0−) is strictly positive. Using the symmetric local timeL as in [13] the equation (2) writes so that the hypotheses of theorem 2.3 in [13] are well satisfied since Therefore Girsanov's theorem implies that the stochastic differential equation (2) has a unique weak solution. To construct a practical discretization scheme for this SDE we use a transformation which removes the local time. Set Denote by β the piecewise linear function β with slope β + on R + and slope β − on R − , and by β −1 is inverse map: Set alsoσ andb Adapting in an obvious way the calculation in Le Gall [13, p.60] we apply Itô-Tanaka's formula (see, e.g., Revuz and Yor [22, Chap.VI]) to β(X t ). The process Y := β(X) satisfies the SDE with discontinuous coefficients: Remark 2.1. The above function β is not the single possible choice to get a stochastic differential equation without local time. One can as well choose any linear by parts function β such that For additional comments in this direction, seeÉtoré [8].
Now denote by h n the step-size of the discretization, that is, h n := T n . For all 0 ≤ k ≤ n set t n k := k h n . Let (Y n t ) be the Euler approximation of (Y t ) defined by Y n 0 = β(X 0 ) and, for all t n k ≤ t ≤ t n k+1 , We then define our approximation of (X t ) by the transformed Euler scheme 3 Main Results

A probabilistic interpretation of the one-dimensional PDE (1)
Suppose that a(x) is smooth everywhere except along smooth discontinuity hypersurfaces S i .  [12] for the definition of this Banach space); this solution is continuous, twice continuously differentiable in space and once continuously differentiable in time on For the sake of completeness and because of its importance in our analysis, we will prove this existence and uniqueness theorem in the one-dimensional case by using stochastic arguments essentially; this approach allows us to get the precise pointwise estimates on partial derivatives of u(t, x) which are necessary to get sharp convergence rate estimates for our transformed Euler scheme.
From now on we limit ourselves to the case d = 1 and we restrict the set of discontinuity points of a(x) = (σ(x)) 2 to {0} (extensions are discussed in section 8). We rewrite (1) and its transmission condition as Suppose also that the function σ is of class C 3 b (R − {0}) and is left and right continuous at point 0. Suppose finally that the first derivative of the function σ has finite left and right limits at 0. Let (X t ) be the solution to (2). Let the bounded function f be in the set Then the function ) and continuous on [0, T ]×R which satisfies (12).
To prove the compatibility transmission condition ( ) and to get uniqueness of the solution to (12), it seems easy to adapt the standard proof of the stochastic representations of solutions v(t, x) of parabolic equations with smooth coefficients (see, e.g., Friedman [10]) which relies on the application of Itô's formula to v(t, X t ). Here, as the first space derivative of u(t, x) is discontinuous for all t, one would rather need to apply a formula of Itô-Tanaka type. However the classical Itô-Tanaka's formula cannot be extended to functions which depend on time and space: see, e.g., Protter and San Martin [24]. To circumvent this difficulty we will use a trick due to Kurtz used in Peskir [21].

Smoothness properties in
The next theorem, which will be proven in section 5, is interesting in its own right from a PDE point of view since it provides accurate pointwise estimates on the derivatives of the solution to (12) when a(x) is discontinuous. To the best of our knowledge, these estimates are new because they are expressed in terms of f γ,1 norms for reasons which will be clear when we derive the theorem 3.5 below from the theorem 3.4.
Theorem 3.3. (i) Under the hypotheses on the function σ made in theorem 3.1 the probability distribution of X t under P x has a density q X (x, t, y) which satisfies: and (ii) Suppose in addition that the function σ is of class C 4 b (R − {0}) and that its three first derivatives have finite left and right limits at 0. Set Then, for all j = 0, 1, 2 and i = 1, . . . , 4 such that 2j + i ≤ 4, where γ = 1 if 2j + i = 1 or 2, and γ = 3 if 2j + i = 3 or 4, and · γ,1 is defined as in (4).

Convergence rate of our transformed Euler scheme
Our next theorem states that the discretization error of the transformed Euler scheme is of order 1/n 1/2− for all 0 < < 1 2 when the function f belongs to W 4 . It improves the results announced in [17]. The precise description in the right-hand side of (20), and the use of the L 1 (R) norms of the derivatives of f , is necessary to prove the theorem 3.5 below.
We now relax the condition that the functions f and Lf satisfy the transmission conditions in the definition (17) of W 4 .
Under the hypotheses on the function σ made in theorem 3.3-(ii), there exists a positive number C (depending on f ) such that, for all 0 < < 1 2 , all n large enough, and all x 0 in R, Remark 3.6. When the coefficient a(x) is smooth, the convergence rate of the classical Euler scheme is of order 1/n and the discretization error can even be expanded in terms of powers of 1/n: for a survey, see, e.g., Talay [27]. Here the coefficientsb andσ are discontinuous; this explains that we are not able to prove better convergence rates as 1/n 1/2− , Notice also that our Euler scheme (X n t ) converges weakly to (X t ) since (Y n t ) converges weakly to (Y t ): see Yan [28]. Remark 3.7. One cannot let tend to 0 in (20) and (22) in spite of the fact that the constants C do not depend on . A more precise statement would be that the absolute value of the error is bounded by C √ hn φ(n) , where φ is a function which, as n tends to infinity, tends to infinity more slowly than any power of n.

Proof of Theorem 3.1
In the calculations below we will use several times the two following observations. First, for all function g of class C 2 b (R − {0}) having a second derivative in the sense of the distributions which is a Radon measure and satisfying the transmission condition the Itô-Tanaka formula applied to g(X t ) and the definition (18) Second, let σ + (x) be an arbitrary C 3 b (R) extension of the function σ(x)I x> which satisfies, for a + (x) : Denote by (X + t ) the unique strong solution to Let τ 0 (X) be the first passage time of the process (X t ) at point 0: Of course, a similar representation holds true for all x < 0 provided the introduction of a diffusion process X − obtained by smoothly extending σ(x)I x< .
First step: smoothness and boundedness. In this paragraph we prove that the function In the rest of this paragraph, w.l.g. we limit ourselves to the case x > 0.
In view of the representation (24) with φ ≡ f and the theorem A.1 in Appendix, we easily deduce the continuity of u(t, x) w.r.t. t and x. Notice that, in particular, the second and third equalities in (12) are satisfied.
Next, to study the boundedness of the function ∂ x u(t, x), we differentiate the flow of (X + t ): Integrate by parts the stochastic integral in the right-hand side; there exists a bounded continuous function G such that: Therefore We then consider the two last terms of the right-hand side of (24). In view of (23) we are in a position to use the lemma A.6 in Appendix with where L + is the infinitesimal generator of the process (X + t ), that is, Therefore We proceed similarly to prove that noticing that, from (25), and that we here can apply the lemma A.6 in Appendix with α = 0. We finally justify that ∂ x u(t, x) has left and right limits when x tend to 0. Indeed, let (x n ) n≥0 a sequence of positive real number tending to 0. We deduce from (26) that (∂ x u(t, x n )) n≥0 is a Cauchy sequence. Denote by M its limit. Let (x n ) n≥0 be another sequence of positive real numbers tending to 0. As the sequence (∂ x u(t,x n )) n≥0 also tends to M , which shows that ∂ x u(t, 0+) is well defined. We similarly obtain that ∂ x u(t, 0−) is also well defined.
Second step: u(t, x) satisfies the first equality in (12). In view of (23) we have, for all 0 < t < T , 0 < < T − t and x in R, Changing φ into Lf in (24) shows In addition, we have already reminded that in [13] the process (X t ) is shown to be strong Markov. Therefore, Itô's formula leads to Divide by the left and right-hand sides and observe that, for all x = 0, Applying Lebesgue's Dominated Convergence theorem we deduce In view of the representation (67) of the density r x 0 (s) in Appendix we have r x 0 (0) = 0 and thus, again applying Lebesgue's Dominated Convergence theorem, we have, for all x = 0, Third step: u(t, x) satisfies the transmission condition ( ). In view of of the preceding first step, for all fixed t the second partial derivative w.r.t. x of u(t, x) is a Radon measure. Thus we may apply the Itô-Tanaka formula to u(t, X s ) for 0 ≤ s ≤ and fixed time t. Our first step also ensures that the resulting Brownian integrals are martingales. Therefore (30) Observe that the equality (28) holds true for x = 0 since it only results from the Markov property of (X t ) and that, combined with (27) it leads to Therefore we deduce from (30) that Since Lf and Lu(t, ·) are bounded functions, the compatibility transmission condition ( ) will be proved if we show that To this end, set Φ(x) := x 0 1 a(y) dy. Observe that the condition (13) implies that Φ is one-to-one.
Similarly to what we did to get (23) we apply Itô-Tanaka's formula to Φ(X t ) and get where (β t ) is the DDS Brownian Motion of the martingale M t := The desired result (31) follows.
Last step: uniqueness. We finally prove that u(t, x) := E x f (X t ) is the unique solution to (12) in the sense of theorem 3.1. We adapt a trick due to Kurtz used in Peskir [21,Sec.3].
As, for all real number x, x ∨ 0 = 1 2 (x + |x|) and x ∧ 0 = 1 2 (x − |x|), Itô-Tanaka's formula implies (We remind that we use the non-symmetric local time corresponding to sgn(x) = I x> − I x≤ .) Now, let U (t, x) be an arbitrary solution to (12). For all fixed t in [0, T ] the function and its partial derivatives have left and right limits when x tends to 0. Thus we may apply the classical Itô's formula to this function and the semimartingales (X s ∨ 0) and (X s ∧ 0). As the resulting Brownian integrals are martingales we obtain: Similarly, We finally use that . In view of the first equality in (12) it comes: It now remains to use that, by hypothesis, U (t, x) satisfies the transmission condition ( ). That ends the proof. In this subsection we closely follow a part of the proof of Aronson's estimate (see, e.g., Bass [3, chap.7, sec.4] and Stroock [26]). We detail the modifications of the classical calculations for the sake of completeness.
We start with observing that, owing to the condition transmision satisfied by fonctions in W 2 , integrating by parts leads to similarly, in view of theorem 3.1, P t f (x) satisfies the transmission condition ( ), two successive integrations by parts lead to Next, setting P t f (x) := E x f (X t ) we have An obvious approximation argument shows that all function g such that g 1 ≤ 1 can be approximated in L 1 (R) norm by a sequence of functions in Therefore In view of the theorem 3.1, P s f and P s g satisfy the condition transmission ( ) for all 0 < s < T . Therefore (33) implies that, for all 0 < s < t < T , It thus remains to prove: since a density argument and the linearity of the operator P t would then also lead to We observe that, in view of (32) and (13), As in the proof of Aronson's estimates we apply Nash's inequality (see, e.g., Stroock [26]) where H 1 (R) is the Sobolev space of functions in L 2 (R) with derivative in the sense of the distributions also in L 2 (R). It comes: from which (34) follows. We thus have proven (16). Notice that, by choosing f as a smooth approximation of the indicatrix function of an open interval not including 0, the preceding inequality implies that the probability distribution of X t under P x has a density and that this density, denoted by q X (x, t, y), satisfies (15). We thus have proven the part (i) of theorem 3.3.

Part (ii): Estimates for time partial derivatives of u(t, x)
In all this subsection, all the constants C do not depend on the function f in W 4 . The calculation is directed by the need to get bounds in terms of · ,1 norms of f rather than in · ∞ norms of its derivatives.
Proof. As above, w.l.g. we may and do consider x > 0. We start from (24) and write where We have and thus Successively using the inequality (16) and the lemma A.5 in Appendix we obtain We now use the following well known estimate (see, e.g., Friedman [10]): for all t > 0, the probability density q X + (x, t, y) of X + t under P x satisfies From the Itô formula and the preceding inequality we have In view of (36) we thus are in a position to obtain (35).
As already noticed, the representation (67) of r x 0 (s) shows that r x 0 (0) = 0. Thus, from the equality (23) with g ≡ Lf (remember that f belongs to W 4 ) and the above calculations, we may deduce We have shown the following proposition: Proposition 5.2. There exists C > 0 such that, for all t ∈ (0, T ],

Part (ii) (cont.): Estimates for space partial derivatives of u(t, x)
The objective of this subsection is to prove estimates for the four first spatial derivatives of u(t, x). As in the preceding subsection, all the constants C do not depend on the function f in W 4 .
Proposition 5.3. There exists C > 0 such that, for all t ∈ (0, T ], Proof. In this proof the various constants C do not depend on the function f . As above w.l.g.
In view of (25) and the Gaussian estimate (40) we have Therefore it suffices to prove In view of (38) this inequality results from lemma A.6 applied to the function noticing that, in view of (16), we may choose Corollary 5.4. There exists C > 0 such that, for all t in (0, T ], Proof. It suffices to use ∂ t u(t, x) = Lu(t, x) for all t in (0, T ] and all x = 0, and to use the estimates in propositions 5.1 and 5.3.
Proposition 5.5. There exists C > 0 such that, for all t ∈ (0, T ], Proof. Since As in the proof of proposition 5.3 we fix x > 0 and start from equality (36): Proceeding as in the proof of (42) we first get Second, to estimate ∂ x ∂ t v(t, x) we use (39) and proceed as in the proof of proposition 5.3. In particular, we apply lemma A.6 to the functioñ , noticing that, in view of (16) and (23) so that we may choose CH := C f 3,1 and α = 1 2 .
Corollary 5.6. There exists C > 0 such that, for all t ∈ (0, T ], 6 Convergence rate of our Euler scheme (I): Proof of theorem 3.4

Error decomposition
For all k ≤ n set θ n k := T − t n k . The proof of theorem 3.4 proceeds as follows. Since u(0, x) = f (x) and u(T, x) = E x f (X T ) for all x, the discretization error at time T can be decomposed as follows: and thus Let us check that the last term in the right-hand side can be satisfyingly bounded from above. As u(0, x) = f (x) for all x, we have Since f is in L 1 (R), f is bounded and thus f • β −1 is Lipschitz. In view of the inequality (35) we deduce The rest of this section is devoted to the analysis of where the time increment T k is defined as and the space increment is defined as In all the calculation below, we use the following notations: given some real number r(n) depending on n, and two positive numbers µ and ν, and We briefly sketch our methodology to study the convergence rate of our Euler scheme. We then distinguish two cases. On the one hand, when Y is small, we explicit the expansion of u(t n k+1 , ·) around 0 and use the theorem 3.1; these two calculations allow us to cancel the lower order term in the expansion. We emphasize that using the transmission condition ( ) is natural: it results from the construction of the approximation scheme by means of the function β −1 whose derivatives are discontinuous at 0.
We again emphasize that the estimate (20) is made necessary to prepare the proof of theorem 3.5 which relies on approximations of functions f in W by sequences of functions in W 4 .
In all the sequel x 0 is arbitrarily fixed.

A preliminary estimate on our Euler scheme
Lemma 6.1. Under P x 0 , for all k ≥ 1, the random variable Y n t k has a density p n t n k w.r.t. Lebesgue's measure. The function p n t n k belongs to C ∞ (R). In addition, there exist C k (n) > 0 and λ k (n) > 0 such that Proof. To prove existence and smoothness of the density p n t n k , we aim to apply the classical lemma 2.1.5 in Nualart [19]. Denote by µ n t n k (dy) the law of Y n t k . Conditionnally to the past up to time t k−1 , the law of Y n t k is Gaussian. Therefore, for all integer α, there exists C α (n) > 0 such that for all test function φ in C ∞ (R) with compact support. Inequality (53) is deduced by induction.

Estimate for the time increment T k
Remember the definition (49) of T k and that θ n k = T − t n k . We have Similarly, In view of theorem 3.3 we have From the preceding we deduce (55)

Expansion of the space increment S k
Let S k be defined as in (50). Set We emphasize that, due to the dissymmetry of the definition β −1 , k+1 X n does not coincide with X n t n k+1 − X n t n k when X n t n k+1 and X n t n k have opposite signs, which explains the two notations and . However the definitions (8) and (10) imply We need to introduce the four following events: In view of the definition of the function β −1 in section 2 we have Similarly, We now use that Ω ++ belongs to the σ-field generated by (B t ) up to time t n k+1 . In view of (56) we get Proceeding similarly and expliciting the conditional expectation of ( k+1 X n ) 2 w.r.t. the past of (B t ) up to time t n k , we obtain In addition, in view of theorem 3.3 we have To summarize the calculations of this subsection, we have obtained We now estimate the remaining term E x 0 R k .
6.5 Estimate for E x 0 R k : localization around 0 Arbitrarily fix 0 < < 1 2 . We aim to show To get this precise estimate we need to use the transmission condition ( ) in equation (12). This explains that we localize on the event where Y n t n k is close to 0. We start with checking that we may neglect the complementary event.
Define Γ(y) by Observe that , and, similarly, In addition,in view of (56) the sum of the three first terms in the right-hand side reduces to so that now are in a position to use the transmission condition ( ) in equation (12). Remembering the definition (5) of β + and β − we deduce that the preceding expression is null. We may proceed similarly as above on the event Ω −+ * k . We thus have proven (60), which ends the proof of (59).

Summing up
Gather the expansions (55) and (58). Use the equality (54) and the inequalities (59), (48). It comes: To deduce (20) it now remains to apply the theorem A.9 in Appendix to the Itô process (Y We thus approximate f in W by functions f δ in W 4 . This approximation will be studied in the L 1 (R)-norm for a reason explained in remark 7.1 below. Let f be in W. For all 0 < δ < 1 define the approximating function f δ in W 4 as follows: where the polynomial functions (p j ) 0≤j≤4 are defined on [0, 1] and are solutions of the following interpolation problem: p As, for all i = 1, . . . , 4, where the constant C here depends on f .
Remark 7.1. Our final error estimates highly depend on the fact that the family (f δ ) approximates f at a good rate in L 1 (R) norm. This explains that f 1 is involved in these error estimates. When the error estimates are expressed in terms of f ∞ or f L 2 (R) the convergence rates are lower than those obtained here: see the section 8 below for additional comments.

Error analysis
The function f now belonging to W, we consider the function f δ as in subsection 7.1. We have Estimate for I 1 (δ). In view of the estimate (15) we have Estimate for I 2 (δ). Using theorem 3.4 and the inequalities (61) we have Notice that the constant C here depends on f because of the smoothing procedure.
the theorem 3.4 and the inequality (15) lead to Thus Global estimate. Gathering the preceding estimates we obtain Choose δ of the type h α n . The optimal value for α is α = 1 4 . The desired result follows.

Extensions and conclusion
Some extensions of our results can readily be obtained.
In the case where a(x) has a finite number of discontinuities, one can split the real line into intervals whose boundary points are the discontinuity points of a(x) and introduce transmission conditions at each of these points. One can also construct an explicit transformation β removing the local time of (X t ) at these discontinuity points. Thus one can readily extend our transformed Euler Scheme. All the results in section 3 still hold true provided straightforward modifications in the calculations made in section 6. Now consider the equation Compatibility transmission conditions at the discontinuity points of a(x).
If the bounded function b is smooth enough (e.g. b is in C 6 b (R)), one can represent the solution of (62) by means of a SDE similar to (2) except that the drift term involves b(X t ), a new modified Euler scheme can easily be constructed, and all our results remain true.
Another easy extension concerns the norm of f used in the convergence rate estimates. Under the hypothesis of theorem 3.4 (that is, the initial function f belongs to W 4 ) one can prove that there exists a positive number C (depending on f ) such that When one wants to use directly (63) to study the case where f is not in W 4 (f belongs to W, say), it appears that the approximation procedure in section 7.1 leads to bounds which explode more rapidly than ours in terms of the smoothing parameter δ: where C now depends on the L ∞ (R) norms f (i) ∞ (0 ≤ i ≤ 4) of f . Choose δ of the type h α n . For α = 1/4 the right-hand side becomes Of course, one also gets error estimates in terms of L ∞ (R − I) norms and L 1 (I) norms of f , where I is a suitably chosen interval around 0.
Other issues to address in the future are: first, to extend to the multi-dimensional case our probabilistic interpretation of parabolic diffraction problems in terms of Markov processes which can easily be simulated, and to prove that the transition densities of these processes satisfy Aronson's estimates; second, to construct a Euler type scheme whose simulation has a weak complexity and to extend our error analysis to a multidimensional setting, that is, when the divergence form operator writes and the matrix valued function a(x) is discontinuous along hypersurfaces. This program is easy to realize when the discontinuity hypersurfaces reduce to hyperplanes. In more general situations, one possible direction is to extend the work done in Bossy et al. Theorem A.1 (Pauwels [20]). Suppose that there exists k ∈ N such that µ is of class C k+1 b (R) and γ is of class C k+2 b (R). Suppose also Let T > 0 and x < 0. Under P x the first passage time of (Z t ) at point 0 before time T , τ 0 (Z)∧T , has a smooth density r x 0 (s) which is of class C k ((0, T ] × (−∞, 0)) and satisfies r x 0 (s) = Ψ(s, S(x))ρ(T − s, S(x))κ(s, S(x)). (67) The formula analogous to (67) is obtained by everywhere changing S(x) into −S(x) and Q(z) into −Q(−z).
Below we need the following (crude) estimates on the function Ψ and its derivatives: Under the hypotheses of theorem A.1 with k ≥ 2, there exists C > 0 such that, for all 0 < s ≤ T and all x < 0, Proof. For the sake of completeness we sketch here the easy proof. The first inequality results from the boundedness of g.
The second inequality results from the boundedness of g and g and the following estimates derived from (66): To obtain the third inequality we use the change of variable γ = sθ s−θ to get It then remains to differentiate w.r.t. s and to observe that ER −x γ ≤ C(|x| + √ γ).
Similar arguments lead to the last inequality. The shortest way seems to differentiate w.r.t. x the expression just obtained ∂ s Ψ(s, x), keeping in mind that |∂ x R x t | ≤ 1.
Concerning the function ρ, we have the following result: Integrating by parts two times the last integral w.r.t. ξ leads to |∂ s ρ(T − s, x)| ≤ C exp(C|x|).
We are now in a position to prove the lemmas A.5 and A.6 that we used in section 5.3.
From the definition of κ, we see that Consequently, in order to prove (69), it suffices to obtain Notice that In addition, as κ(s, x) is the density of the first passage time of x by a Brownian motion starting at 0 at time 0, We thus have proven (69).
Lemma A.6. There existsC > 0 such that, for all 0 ≤ α < 1, and all function H bounded on Proof. W.l.g. we again suppose x < 0. We use the theorem A.1 to represent the density r x 0 and observe that  We deduce from inequality (69) that, for some possibly new constantC, Similarly, As for some positive constantsC and , exp(C|x| − x 2 4(t−s) ) ≤C exp(− x 2 4(t−s) ) for all x ∈ R and 0 ≤ s < t ≤ T, we deduce that for some possibly new constantC, We now turn to the bound for ∂ 2 xx t 0 r x 0 (t − s)H(s)ds . From the above we have Remark A. 10. In [4] the statement of the preceding theorem claims that h 0 depends on ε. In fact the proof shows that one only needs that exp(− 1 h ε ) ≤ Ch 3/2−ε .