On the time evolution of Bernstein processes associated with a class of parabolic equations

In this article dedicated to the memory of Igor D. Chueshov, I first summarize in a few words the joint results that we obtained over a period of six years regarding the long-time behavior of solutions to a class of semilinear stochastic parabolic partial differential equations. Then, as the beautiful interplay between partial differential equations and probability theory always was close to Igor's heart, I present some new results concerning the time evolution of certain Markovian Bernstein processes naturally associated with a class of deterministic linear parabolic partial differential equations. Particular instances of such processes are certain conditioned Ornstein-Uhlenbeck processes, generalizations of Bernstein bridges and Bernstein loops, whose laws may evolve in space in a non trivial way. Specifically, I examine in detail the time development of the probability of finding such processes within two-dimensional geometric shapes exhibiting spherical symmetry. I also define a Faedo-Galerkin scheme whose ultimate goal is to allow approximate computations with controlled error terms of the various probability distributions involved.


Introduction and outline
This article is a tribute to some of the works and achievements of our friend and colleague Igor D. Chueshov, who unfortunately and unexpectedly passed away on April 23rd, 2016. The qualitative analysis of the behavior of solutions to various stochastic partial differential equations, henceforth SPDEs, was one of Igor's strong points. I have therefore deemed it appropriate to briefly summarize here the results that he and I obtained in that area over a period stretching from 1998 to 2004. As far as the presentation of the many other facets of his activities is concerned, I am thus referring the reader to the other contributions in this volume.
When Igor and I first met in 1994 on the occasion of an international conference on SPDEs in Luminy, we set out to investigate the behavior of solutions to those stochastic parabolic equations which specifically occur in population dynamics, population genetics, nerve pulse propagation and related topics, given the fact that there were already a substantial number of works in those areas concerning the deterministic case (see, e.g., [3], [20] and the many references therein). But instead of starting up front with partial differential equations driven by some kind of noise, we first considered a class of random parabolic initial-boundary value problems mainly for the sake of simplification. Assuming then various statistical and dynamical properties such as those of the central limit theorem and the Ornstein-Uhlenbeck process for the lower-order coefficients of the equations, we eventually elucidated the ultimate behavior of the corresponding solution random fields in [6]. In particular, we established the existence of a global attractor, determined its detailed structure and were able to compute the Lyapunov exponents explicitly in some cases. We then extended these results to the case of parabolic SPDEs driven by a homogeneous multiplicative white noise defined in Stratonovitch's sense in [7], investigated there various stability properties of the non-random global attractor and established the existence of a recurrent motion of sorts among its components. Furthermore, in [8] we analyzed the same type of equations as in [7] but with the noise defined in Itô's sense. In this way we were able to establish the existence and many properties of a random global attractor and excluded in particular the existence of any kind of recurrence phenomena, thereby obtaining radically different results than in [7]. The analysis carried out in [8] was further deepened in [2], where it was shown that the stabilization of the solution random fields toward the global attractor is entirely controlled by their spatial average, thereby obtaining exchange of stability results particularly relevant to the description of certain migration phenomena in population dynamics. Finally, in [9] we proved the existence of invariant sets under the flow generated by certain systems of SPDEs including those of Lotka-Volterra and Landau-Ginzburg.
But Igor's interests did not limit themselves to investigations of solutions to SPDEs as he was also genuinely interested in the many possible connections that exist between systems of differential equations on the one hand, and the theory of random dynamical systems and stochastic processes on the other hand (see, e.g., [5]). This prompted me to present here some very recent and preliminary results concerning the time evolution of certain Bernstein processes naturally associated with a class of deterministic linear partial differential equations. Accordingly, the remaining part of this article is organized as follows: In Section 2 I recall what a Bernstein process is, and state there a theorem that shows how to associate such a process with the two adjoint parabolic Cauchy problems and where T > 0 is arbitrary and where △ x stands for Laplace's operator with respect to the spatial variable. In these equations N > 0 is a normalization factor whose significance I explain below. Moreover, V is real-valued while ϕ 0 and ψ T are positive data which are assumed to be either Gaussian functions of the form where σ 0,T > 0 and a 0,T ∈ R d are arbitrary vectors with |.| the usual Euclidean norm, or In (5)-(6), x j and a 0,T,j denote the j th component of x and a 0,T , respectively. Furthermore these initial-final conditions have localization properties which are more clear-cut than those of (3)-(4) in that they vanish identically outside hypercubes in R d . The cases where with δ 0 the Dirac measure concentrated at the origin and ψ T given by (4) or (7) are also considered. An important observation here is that (3)-(4) and (5)-(6) are not normalized as standard probability distributions, for the only normalization condition needed below involves ϕ 0 , ψ T and N in a rather unexpected way which is inherently tied up with the construction of Bernstein processes. Finally, the following hypothesis is imposed regarding the potential function in (1)-(2): (H) The function V : R d → R is continuous, bounded from below and satisfies V (x) → +∞ as |x| → +∞.
An immediate consequence of this hypothesis is that the resolvent of the usual self-adjoint realization of the elliptic operator on the right-hand side of (1)-(2) is compact in L 2 C R d , the usual Lebesgue space of all square integrable, complex-valued functions on R d . This means that the operator in question has an entirely discrete spectrum (E n ) n∈N d , and that there exists an orthonormal basis (f n ) n∈N d ⊂ L 2 C R d consisting entirely of its eigenfunctions (see, e.g., Section XIII.14 in [16]). In the context of this article the convergence of the series for every t ∈ (0, T ] is also required. Then, under the above conditions the construction of a Markovian Bernstein process rests on two essential ingredients, namely, Green's function (or heat kernel) associated with (1)- (2), which satisfies the symmetry and positivity conditions g(x, t, y) = g(y, t, x) > 0 for all x, y ∈ R d and every t ∈ (0, T ], and the probability measure µ on R d × R d whose density is given by which satisfies the normalization condition Notice that (11) may be considered as the definition of N , and that the inequality in (9) is a consequence of two-sided Gaussian bounds for g whose existence follows from the general theory developed in [1] and further refined in Chapter 3 of [10]. Moreover, as a consequence of (H) and (8), Green's function admits an expansion of the form which converges strongly in L 2 C R d × R d for every t ∈ (0, T ] (unless more detailed information about the f n 's or ultracontractive bounds become available, in which case the convergence can be substantially improved, see, e.g., Chapter 2 in [10]). Thus, in Section 2 the knowledge of g and µ is used to state a theorem about the existence of a probability space which supports a Markovian Bernstein process Z τ ∈[0,T ] whose state space is the entire Euclidean space R d , and which is characterized by its finite-dimensional distributions, the joint distribution of Z 0 and Z T and the probability of finding Z t at any time t ∈ [0, T ] in a given region of space. In that section a very simple result regarding the time evolution of Z τ ∈[0,T ] is also proved when considering (1)-(2) with (5)- (6). Section 3 is devoted to the analysis of the function that determines the time evolution of the probability of finding Z τ ∈[0,T ] in particular two-dimensional geometric shapes that exhibit spherical symmetry in the case of the so-called harmonic potential and for various combinations of the initial-final data given above. Finally, a simple Faedo-Galerkin scheme is proposed whose ultimate goal is to allow approximate computations of all the probability distributions involved.
2 An existence result for a class of Bernstein processes in R d As a stochastic process a Bernstein process may be defined independently of any reference to a system of partial differential equations, and there are several equivalent ways to do so (see, e.g., [13]). I shall restrict myself to the following: Definition. Let d ∈ N + and T ∈ (0, +∞) be arbitrary. An R d -valued process Z τ ∈[0,T ] defined on the complete probability space (Ω, F , P) is called a Bernstein process if for every bounded Borel measurable function f : R d → R and for all r, s, t satisfying r ∈ (s, t) ⊂ [0, T ]. In (14), the σ-algebras are and where B d stands for the Borel σ-algebra on R d . Moreover, E denotes the (conditional) expectation functional on (Ω, F , P).
The dynamics of such a process are, therefore, solely determined by the properties of the process at times s and t, irrespective of its behavior prior to instant s and after instant t. Of course, it is plain that this fact generalizes the usual Markov property.
In what follows an important rôle is played by the positive solution to (1) and the positive solution to (2), namely, and respectively. Taken together, (1) and (2) may thus be looked upon as defining a decoupled forward-backward system of linear deterministic partial differential equations, with (17) wandering off to the future and (18) evolving into the past. The functions and with F ∈ B d , both being well defined and positive for all x, y, z ∈ R d and all r, s, t satisfying r ∈ (s, t) ⊂ [0, T ], are equally important as is the probability measure µ whose density is (10), namely, where G ∈ B d × B d , which satisfies the normalization condition (11). The corresponding initial and final marginal distributions then read respectively, as a consequence of (17) and (18). It is the knowledge of (20) and (21) that makes it possible to associate with (1) and (2) a Bernstein process in the following sense: Theorem. Assume that V satisfies Hypothesis (H), that condition (8) holds and that P and µ are given by (20) and (21), respectively. Then there exists a probability space (Ω, F , P µ ) supporting an R d -valued Bernstein process Z τ ∈[0,T ] such that the following properties are valid: (a) The process Z τ ∈[0,T ] is Markovian, and the function P is its transition function in the sense that for each F ∈ B d and all r, s, t satisfying r ∈ (s, t) ⊂ [0, T ]. Moreover, The finite-dimensional probability distributions of the process are given by for every integer n ≥ 2, all F 1 , ..., F n ∈ B d and all t 0 = 0 < t 1 < ... < t n < T , where u and v are given by (17) and (18), respectively.
(c) The probability of finding the process in a given region F ⊂ R d at time t is given by for each F ∈ B d and every t ∈ [0, T ] .
(d) P µ is the only probability measure leading to the above properties.
I omit the proof of this theorem, which can be adapted either from the abstract arguments in [13] or from the more analytical approach in [21], and will rather focus on its consequences regarding the time evolution of Z τ ∈[0,T ] . Prior to that some comments are in order: Remarks. (1) Hypothesis (H) and condition (8) are sufficient but not necessary for the theorem to hold. However, the advantage of having (12) is that such an expansion greatly simplifies some calculations and also has the virtue of making theoretical results amenable to approximations and computations. I will dwell a bit more on this point in the next section.
(2) Bernstein processes may be Markovian but in general they are not. Independently of that they have played an increasingly important rôle in various areas of mathematics and physics over the years. It is not possible to give a complete bibliography here, but I will refer instead the interested reader to [13], [17] and [21] which contain many references describing the history and earlier works on the subject, tracing things back to the pioneering works [4] and [18]. Moreover, Bernstein processes have also lurked in various forms in more recent applications of Optimal Transport Theory, as testified by the monographs [12] and [19]. In this regard it is worth mentioning that they are also referred to as Schrödinger processes or reciprocal processes in the literature.
(3) The probability measure µ of a non-Markovian Bernstein process does not have as simple a structure as that given by (21). A case in point is the so-called periodic Ornstein-Uhlenbeck process, which is one of the simplest stationary Gaussian non-Markovian processes that can be viewed as a particular Bernstein process, as was recently proved in [22] (see also, e.g., [17] and the references therein for other analyses of the periodic Ornstein-Uhlenbeck process). In this case the construction of the measure µ is much more complicated than in the Markovian case, as it involves a weighted average of a sequence of suitably constructed signed measures naturally associated with an infinite hierarchy of forward-backward linear parabolic equations.
Coming back to the main theme of this article, it is interesting to note that the probability of finding the process at any given time t ∈ [0, T ] in an arbitrary region of space is expressed as an integral of the product of u and v through the simple formula (24). This is a manifestation of the fact that the process Z τ ∈[0,T ] is actually reversible and exhibits a perfect symmetry between past and future, a property already built to some extent into the definition given at the beginning of this section. It is of course difficult to say more about the time evolution of Z τ ∈[0,T ] unless we know more about the potential function V .
However, at the very least the following result holds, which in effect describes a recurrence property of the process in a particular case: Proposition 1. Let Z τ ∈[0,T ] be the Bernstein process associated with (1)- (2) in the sense of the above theorem, where ϕ 0 and ψ T are given by (5) and (6), respectively, and let C a0,σ0 = x ∈ R d : |x j − a 0,j | < σ 0 , j = 1, ..., d be the hypercube outside which ϕ 0 vanishes identically, that is, ϕ 0 = 0 on F a0 , σ 0 = R d \ C a0,σ0 . Let C aT ,σT be defined in a similar way. Then Proof. This is an immediate consequence of (24), for Thus, in this case the process certainly starts its journey within C a0,σ0 and ends it within C aT ,σT . Since this is true no matter how small σ 0,T are, that constitutes a generalization of the class of Bernstein bridges constructed in [22]. In particular, if a 0 = a T and if σ T ≤ σ 0 the inclusion C aT ,σT ⊆ C a0,σ0 holds, so that the process goes back to the region where it started from with probability one, independently of its unknown whereabouts at intermediary times t ∈ (0, T ). These properties and Proposition 1 remain true for all choices of ϕ 0 , ψ T that vanish identically outside of a given Borel set, for instance for the isotropic version of (5)-(6), namely, provided the sets C a0,T ,σ0,T are replaced by the d-dimensional open balls B a0,T ,σ0,T = x ∈ R d : |x − a 0,T | < σ 0,T of radius σ 0,T centered at a 0,T . It would be interesting to carry out a numerical simulation in real time of the behavior of the processes generated in this way. The preceding result fails to hold if the initial-final data are not of the above form. In the next section I investigate this issue more closely in case the potential function is given by (13).

Some new results for the harmonic case
The starting point is thus the forward-backward system and Green's function associated with (25)-(26) is known to be Mehler's multidimensional kernel where (., .) R d denotes the usual inner product in R d (see, e.g., the Appendix in [22]). Then if ϕ 0 , ψ T are given by (3)-(4), the solutions (17)- (18) and the integral on the left-hand side of (11) can all be computed explicitly since the integrals are Gaussian. For instance, the forward solution reads for every t ∈ (0, T ], while the backward solution is obtained from (28) by replacing σ 0 by σ T , a 0 by a T and t by T − t, respectively. The downside is that these expressions are complicated, cumbersome and in any case unsuited to extract valuable information out of (24) unless particular choices are made for these parameters. For example, if σ 0 = σ T = 1 and a 0 = a T = 0, the forward solution (28) and the related backward solution reduce to respectively, while an explicit computation from (11) gives for the corresponding normalization factor. Therefore, the substitution of these expressions into (24) leads to for each t ∈ [0, T ] and every F ∈ B d , so that the probability of finding the process in any region of space is here independent of time. The reason for this independence can easily be understood by means of the substitution of (27) and (29)-(30) into (23), which first leads to the Gaussian law of (Z t1 , ..., Z tn ) ∈ R nd and from there eventually to the covariance for all s, t ∈ [0, T ] and all i, j ∈ {1, ..., d}, where E µ denotes the expectation functional on the probability space of the theorem. Therefore, the Bernstein process thus constructed identifies in law with the standard d-dimensional Ornstein-Uhlenbeck velocity process, so that the choice of (3)-(4) as initial-final data corresponds in a sense to an equilibrium situation whereby the law remains stationary (see, e.g., [14] for general properties of this and related processes). For instance, if is the two-dimensional annulus centered at the origin with R 1 ≥ 0 and R 2 > 0, then The situation is quite different if the system (25)-(26) is considered with ϕ 0 given by (7) and ψ T given by (4) where σ T = 1 and a T = 0. In this case for the forward and backward solutions, respectively, and furthermore the value of N can again be determined directly from (11). Indeed the relevant integral is by virtue of (27), so that Therefore, one obtains in particular so that the process is conditioned to start at the origin since Moreover, for positive times an explicit evaluation from (24) leads to where the width parameter is identified as It is then instructive to consider again the case of Z τ ∈ [0,T ] wandering in the two-dimensional annulus A R1,R2 , and to investigate the way that as soon as t > 0. Moreover, if T is sufficiently large there exists a t * ∈ (0, T ) such that the function t → P µ (Z t ∈ A R1,R2 ) is monotone decreasing for every t ∈ [t * , T ] .
Proof. Statement (a) follows immediately from (33) and (35) for R 1 = 0, as does the very first part of (b) since then o / ∈ A R1,R2 . Now where and for any fixed t ∈ (0, T ] this function is monotone increasing for R < 2ρ(t) and monotone decreasing for R > 2ρ(t). Furthermore, (34) and t → 2ρ(t) are monotone increasing and concave with 2ρ(t) < 1 uniformly in t. Therefore, if 0 < R 1 < R 2 < 1 and if T is large enough, there exists a t * ∈ (0, T ) such that A natural interpretation of Statement (a) is that the process leaves the origin as soon as t > 0, and tends to quickly "leak out" of the disk A 0,R2 when R 2 is sufficiently small. Moreover, Statement (b) means that the probability of finding the process in the annulus increases for small times, then reaches a maximal value and eventually decreases for large times when R 1 and R 2 are sufficiently small, in sharp contrast to Statement (c) where the probability in question is monotone increasing for all times if R 1 and R 2 are sufficiently large. Finally, the substitution of (27) and (31)-(32) into (23) again determines the projection of the law onto R nd and, after long algebraic manipulations, the covariance for all s, t ∈ [0, T ] and all i, j ∈ {1, ..., d}. Therefore, the Bernstein process thus constructed is identical in law with the Ornstein-Uhlenbeck process conditioned to start at the origin of R d . A last example can be provided by choosing ϕ 0 and ψ T both of the form (7) in (25)-(26). In this case one gets  Arguing as in the preceding example one then obtains so that the process is conditioned to start and end at the origin, thereby representing a random loop in R d . Moreover, for positive times one still gets from (24) and in particular in the case of the two-dimensional annulus, but with a width parameter now given by for every t ∈ [0, T ]. This function is quite different from (34), and the following result is valid: The following statements hold: Moreover, the function t → P µ (Z t ∈ A 0,R2 ) is monotone decreasing on 0, T 2 and monotone increasing on T 2 , T , thereby taking the minimal value Moreover, the function t → P µ (Z t ∈ A R1,R2 ) is monotone increasing on 0, T 2 and monotone decreasing on T 2 , T , thereby taking the maximal value Proof. While (41) follows from (38), Relation (39) with R 1 = 0 leads to where ρ ′ (t) ≥ 0 for t ∈ 0, T 2 and ρ ′ (t) ≤ 0 for t ∈ T 2 , T according to (40), which implies Statement (a). Statement (b) follows from these properties of ρ ′ and an analysis similar to that of Statement (c) in Proposition 2. Indeed, we remark that the curve ρ : [0, T ] → [0, +∞) given by (40) is concave aside from satisfying ρ(0) = ρ(T ) = 0, and that it takes on the maximal value at the mid-point of the time interval. Therefore, the inequalities hold for every t ∈ [0, T ], which implies that (37) is monotone decreasing throughout the time interval as a function of R, a consequence of the hypothesis regarding the radii.
The above properties of (40) thus show that the Bernstein process of Proposition 3 constitutes a generalization of a Brownian loop, that is, of a particular case of a Brownian bridge (see, e.g., [14]). This renders the preceding result quite natural, in that the probability of finding the process in the disk A 0,R2 is minimal at the mid-point of the time interval where there is maximal randomness. At the same time, the situation is reversed if the annulus is relatively far away from the origin.
As long as the regions of interest are spherically symmetric, the preceding calculations may be performed in any dimension and not merely for d = 2. However, I shall refrain from doing that and rather focus briefly on what to do when the values of the parameters σ 0,T and a 0,T are arbitrary, or when other combinations of the above initial-final data are chosen. It is here that an expansion of the form (12) is essential, and I will now show what (12) reduces to in the case of (27). First, the spectral decomposition of the elliptic operator on the right-hand side of (25)-(26) is known explicitly (the operator identifies up to a sign with the Hamiltonian of an isotropic system of quantum harmonic oscillators, see, e.g., [15]). Indeed, let (h n ) n∈N be the usual Hermite functions where the H n 's stand for the Hermite polynomials Then, it is easily verified that the tensor products ⊗ d j=1 h nj where the n j 's run independently over N provide an orthonormal basis of eigenfunctions in L 2 C R d which satisfy the eigenvalue equation for each n ∈N and every x ∈ R d , where n = (n 1 , ..., n d ) ∈ N d and The immediate consequences are that (8) holds, and that expansion (12) for (27) takes the form where the series is now absolutely convergent for each t ∈ (0, T ] uniformly in all x, y ∈ R d . This very last statement follows from Cramér-Charlier's inequality valid uniformly in n, x and y, where k ≤ 1.086435 (see, e.g., Section 10.18 in [11] and the references therein). The advantage of having (46) is that the forward solution (17) may now be rewritten in terms of the Fourier coefficients of ϕ and ψ along the basis (h n ) n∈N d , namely, which in case of Gaussian initial-final data provides a nice representation of (28). In a similar way the backward solution (18) is so that the normalization condition (11) now reads This way of formulating things, in turn, leads to the possibility of constructing a sequence of Faedo-Galerkin approximations to the problem at hand. Thus for any positive integer N ≥ 1, let E N R d be the N d -dimensional subspace of L 2 C R d generated by the h n 's where n j ∈ {0, ..., N − 1} for each component of n. Green's function (46) may then be approximated by and to (48) and (50), respectively. Consequently, various numerical computations and controlled approximations of the probability distributions of interest now become possible. I complete this short article by a simple illustration of this fact stated in Proposition 4 below, whose proof is based on the following result which provides an approximate value for N : Lemma. Let (52) be written as n exp [−T E n ]β n = 1 whereα n = N −1 α n , β n = N −1 β n for every n ∈ N d . Then for all σ 0,T > 0, a 0,T ∈ R d , the unique positive solution to is of the form where c > 0 is a constant depending only on σ 0,T and a 0,T . Moreover, with the value (57) in (49) and (51) one gets for T sufficiently large.
Proof. It is clear that (57) holds because of (56) sinceα 0 > 0,β 0 > 0 by virtue of the fact that the eigenfunction h 0 associated with the bottom of the spectrum is strictly positive in R d . Then, the proof that the remaining term satisfies follows from the fact that theα n 's and theβ n 's are uniformly bounded in n, and from the summation of the underlying geometric series which is made possible thanks to the explicit form (44).
In particular, (b) If σ 0 = σ T : = σ and a 0 = a T := a and if the process Z τ ∈[0,T ] is stationary, the preceding relations reduce to for T large enough, each t ∈ [0, T ] and every F ∈ B d , where ρ = σ 1+σ .
Proof. From (21) and (22) one has where g is given by (46), that is, according to (44) and (45) for n = 0. One then obtains and replacing N by N 0,T together with the explicit evaluation of these Gaussian integrals gives the leading term in (58). It remains to show that the contribution to (58) coming from the second term on the right-hand side of (61) is exponentially small. Writing momentarilỹ g(x, T, y) = because of (57), as desired. Finally (58) implies (59), and also (60) under the hypothesis in (b) since the function t → P µ (Z t ∈ F ) given by (24) is then independent of t.
Using once more (44) and (45) with n = 0, the corresponding approximation (54) for t = T then reads so that replacing N by (57) and arguing as in the above proofs one eventually gets for T sufficiently large. A similar approximation procedure applies to the backward solution, so that in the end one obtains yet another algorithm to compute (59) sinceα 0 andβ 0 can be determined explicitly in case of Gaussian initial-final data. It would have been difficult to evaluate (64) directly from (24) given the complicated form (28). As a matter of fact, the technique used also works if the data are of the form (5)-(6) sinceα 0 andβ 0 are then easily determined by numerical calculations. More generally, there is an important computational issue about (54) and (55), namely, that of knowing how large one has to choose N as a function of the desired degree of precision to reconstruct u and v. As long as error terms of the form O (exp [−T ]) are considered satisfactory, the above considerations show that the choice N = 1 is sufficient. If not, larger values of N will do.
Finally, thanks to an expansion of the form (12), similar Faedo-Galerkin approximation methods may be applied to the forward-backward solutions of (1)-(2) when the potential function satisfies Hypothesis (H) and (8), or even more general conditions, provided that precise information be available about the spectrum (E n ) n∈N d and the corresponding sequence of eigenfunctions (f n ) n∈N d . The detailed results will be published elsewhere.