Stochastic geometric mechanics in nonequilibrium thermodynamics: Schrödinger meets Onsager^*

Instead to repeat what has been done and elaborated in the last 35 years [16], we summarize briefly here the main results relevant to our present stochastic dynamical purpose.

Schrödinger considers the simplest system of classical Brownian particles whose dynamics is described by the self-adjoint Hamiltonian operator

$$\begin{equation} \mathcal H = \mathcal H^* = -\frac{1}{2} \Delta + U(x), \end{equation} \tag{ 1 }$$

on $L^2({\mathbf{R}}^n)$, where $U:{\mathbf{R}}^n\to{\mathbf{R}}$ is a bounded below potential. Observing that the associated quantum unitary group $\exp (-it \mathcal H)$ can be analytically continued, via $t \to -it$, into the self-adjoint semigroup in $L^2(\mathbf{R}^n)$:

$$\begin{equation*} T(t) = e^{-t \mathcal H}, \end{equation*}$$

so that, for any positive ψ in a dense domain of $\mathcal H$, $w^*_\psi(x,t) = (e^{-(t-t_0) \mathcal H} \psi) (x)$ solves the initial value problem, $t\geqslant t_0$

$$\begin{equation} -\frac{\partial w^*}{\partial t} = \mathcal H^* w^*, \quad w^*(t_0,\cdot) = \psi. \end{equation} \tag{ 2 }$$

Schrödinger referred to this first Euclidean step as a banal 'superficial analogy' with his wave equation. Strangely enough, for decades, his few readers will not go beyond this point.

But Schrödinger takes also note of the critical fact that, if $t\leqslant t_1$, we could define as well, on appropriate domain, a final value problem for $w_\phi(x,t) = (e^{-(t_1-t) \mathcal H} \phi) (x)$, namely,

$$\begin{equation} \frac{\partial w}{\partial t} = \mathcal Hw, \quad w(t_1,\cdot) = \phi. \end{equation} \tag{ 3 }$$

Let us stress that our $*$ notation does not mean that $w^*$ and w are trivially related as for complex L² functions. Conditions on ψ and φ will be given later.

Now, considering the same system, for $t\in I = [t_0,t_1]$, Schrödinger's original thought experiment and problem [17] was:

Imagine a Brownian particle with Hamiltonian $\mathcal H$ evolving from time t₀ with an initial distribution $\rho_{t_0}$. Suppose we observed the state of the particle at time $t_1\gt t_0$ with another distribution $\rho_{t_1}$. What is the 'most likely' evolution of the probability distribution ρ_t of this particle at time $t\in [t_0,t_1]$?

Schrödinger used an entropic argument and a discretization technique, familiar in statistical physics, but he knew very well what he wanted to achieve and suggest. His justification has been made rigorous many times since the sixties. As a sample see [3, 4, 18–20]. His result is that this most probable diffusion, say X(t), satisfies a Born-type formula, that is, for any Borelian U,

$$\begin{equation} \mathbf{P}(X(t)\in U) = \int_U w w^* (t,x) \mathrm dx, \end{equation} \tag{ 4 }$$

or, equivalent, the probability distribution ρ_t has a specific (integrable) product form

$$\begin{equation} \rho_t(\mathrm dx) = w w^* (t,x) \mathrm dx, \end{equation} \tag{ 5 }$$

where $w^*$ and w are positive solutions of (2) and (3) respectively.

First, a remark about the physical interpretation of Schrödinger's idea. He is clearly describing relatively rare events. For $\mathcal H$ we expect, on the basis of our macroscopic experience, irreversible dynamics given by equation (2). But, in fact, any good physics experiment is organized in such a way that mostly interesting events are collected. And without a careful organization of the set-up, they could be pretty rare. In this sense, any physics experiment is very conditioned. This plays a key role in the interpretation of the consequences of Schrödinger 'gedanken experiment'.

With the benefits of knowing 90 years of scientific history since Schrödinger's observation, we can affirm now what he could not guess, back then. His suggestion was a paradigm shift about stochastic dynamics, whose impact was accidentally delayed for historical reasons.

The unorthodox data of ρ at the boundary $\partial I$ of a time interval has indeed remarkable consequences. Formally, w solving (3) is a kind of time reversal of $w^*$. Let us see why this is more than a formal relation. Since $\{\rho_{t_0}(\mathrm dx), \rho_{t_1}(\mathrm dz)\}$ are arbitrary, this diffusion should not be time-homogeneous in general. It has been proved (e.g. [2, 20], or [21] for the forward filtration case) that its forward and backward transition probabilities for $s\leqslant t\in I = [t_0,t_1]$ are, respectively,

$$\begin{align} p(s,x;t,\mathrm dy) &= h(x, t-s, y) \frac{w(t,y)}{w(s,x)} \mathrm dy, \end{align} \tag{ 6 }$$

$$\begin{align} p^*(s,\mathrm dx;t,y) &= \frac{w^*(s,x)}{w^*(t,y)} h(x,t-s,y) \mathrm dx, \end{align} \tag{ 7 }$$

where h denotes the heat kernel $h(x, t-s, y) = (e^{-(t-s) \mathcal H} \delta_x)(y)$ associated with equation (2). Given the product form in (4) $\rho(x,t) = ww^* (t,x)$ of Schrödinger's most probable process, his nonequilibrium generalization of detailed balance condition (see also [2, 16]) becomes

$$\begin{equation} \rho(x,s) p(s,x,t,\mathrm dy) \mathrm dx = p^*(s,\mathrm dx,t,y) \rho(y,t) \mathrm dy, \quad s\leqslant t \in I, \end{equation} \tag{ 8 }$$

which also means that the joint distribution of the Markovian time-reversible diffusion (in Schrödinger sense) of two instants $s\leqslant t\in I$ is of the specific form

$$\begin{equation} \mu_{s,t}(\mathrm dx,\mathrm dy) = w^*(s,x) h(x, t-s, y) w(t,y) \mathrm dx\mathrm dy \end{equation} \tag{ 9 }$$

for $\{w^*(t_0), w(t_1)\} = \{\psi,\phi\}$ unspecified (but positive) boundary conditions of the two heat equations (2) and (3). Since h is positive, these boundary conditions must be as well, for a meaningful joint probability µ.

By construction, the marginals of µ should therefore satisfy the system for $\{\psi,\phi\}$,

$$\begin{equation}\left\{ \begin{aligned} \psi(x) \int h(x,t_1-t_0,y) \phi(y) \mathrm dy = \rho_{t_0}(x), \\ \phi(y) \int h(x,t_1-t_0,y) \psi(x) \mathrm dx = \rho_{t_1}(y), \end{aligned} \right. \end{equation} \tag{ 10 }$$

as found originally by Schrödinger. Also notice that the joint distribution (9) is the Euclidean version of what Feynman will call, years after, a (complex) 'transition element' in [22]. It is crucial to stress that $\psi, \phi$ have no other relation between them than to solve the system (10), so that equation (4) represents, in fact, an Euclidean scalar product.

For a fixed Hamiltonian $\mathcal H$, Bernstein sketched a method [23] to construct systematically processes like the ones of Schrödinger, from the data of an 'arbitrary' joint distribution $\mu_{t_0,t_1}(\mathrm dx,\mathrm dy)$ at the frontiers of a time interval $[t_0,t_1]$. He suggested that such stochastic processes could be called reciprocal processes, as they are defined in the following intrinsically time-reversible fashion:

$$\begin{equation*} \mathbf{P}(X_t\in A | X_u, X_v, X_{s_1}, \ldots, X_{s_n}) = \mathbf{P}(X_t\in A | X_u, X_v), \end{equation*}$$

for any $u\lt t\lt v\in[t_0,t_1]$, $\{s_1, \ldots, s_n\} \subset [t_0,u]\cup[v,t_1]$ and Borelian A.

This strategy was followed by Jamison [21, 24] using Beurling's crucial proof of existence and uniqueness of solutions $(\psi,\phi)$ for system (10) (see [18] and section 5). He observed that, in general, such processes have no reason to be Markovian. They are Markovian when and only when the end-point joint distribution $\mu_{t_0,t_1}$ has the form (9) advocated by Schrödinger, i.e.

$$\begin{equation} \mu_{t_0,t_1}(\mathrm dx,\mathrm dy) = \psi(x) h(x, t_1-t_0, y) \phi(y) \mathrm dx\mathrm dy. \end{equation} \tag{ 11 }$$

In the latter case, in particular, if the Hamiltonian $\mathcal H$ is of the following diffusion-type, with constant diffusive coefficients $(D^{\alpha\beta})$,

$$\begin{equation*} \mathcal H = - D^{\alpha\beta} \frac{\partial^2}{\partial x^\alpha\partial x^\beta}, \end{equation*}$$

then resulting reciprocal processes are again diffusions on I, with infinitesimal (forward) generators,

$$\begin{equation} D^{\alpha\beta} \left[ \frac{\partial^2}{\partial x^\alpha\partial x^\beta} + \frac{\partial (\log w)}{\partial x^\beta} \frac{\partial}{\partial x^\alpha} \right], \end{equation} \tag{ 12 }$$

where w is the solution of equation (3).

Föllmer used large deviation theory to reformulate Schrödinger's variational problem as an entropy minimization problem [19]. Notice that, originally, Schrödinger applied the Monte Carlo method, that is to repeat independently his thought experiment a large number of times, say, N. From a modern viewpoint, this will produce a sequence of i.i.d. processes $X^i(t)$, $i = 1,\ldots, N$, with common law R where R is the law of the Brownian particle with Hamiltonian $\mathcal H$ at time interval $[t_0,t_1]$ and initial distribution $\rho_{t_0}$. The empirical law of the system is

$$\begin{equation*} \mathcal L^N = \frac{1}{N} \sum_{i=1}^{N} \delta_{X^i}. \end{equation*}$$

Schrödinger's problem can be translated into finding the following limit of conditional laws:

$$\begin{equation} \lim_{\epsilon\to0^+} \lim_{N\to\infty} \mathbf{P}\left( \mathcal L^N\in \cdot \mid \mathcal L^N(b)\in B(\rho_{t_1},\epsilon) \right), \end{equation} \tag{ 13 }$$

where B means the ball in the space of all probability laws of processes on $[t_0,t_1]$ with initial marginal distribution $\rho_{t_0}$, equipped with some compatible distance.

A straightforward application of Sanov's theorem implies that for any 'regular' subset A of that space,

$$\begin{equation*} \mathbf{P}\left( \mathcal L^N\in A \right) \stackrel{n\to\infty}{\asymp} \exp\left[ - n \left( \inf_{P\in A} D_{\mathrm{KL}}(P|R) \right) \right], \end{equation*}$$

where $D_{\mathrm{KL}}$ is the Kullback–Leibler divergence (also called relative entropy) defined by

$$\begin{equation*} D_{\mathrm{KL}}(P|R) = \begin{cases} \mathbf{E}_P \left[ \log\left( \frac{\mathrm dP}{\mathrm dR} \right) \right], & P \ll R, \\ +\infty, & \mathrm{otherwise}. \end{cases} \end{equation*}$$

So that the limit in equation (13) is

$$\begin{equation*} \lim_{\epsilon\to0^+} \lim_{N\to\infty} \mathbf{P}\left( \mathcal L^N\in \cdot \mid \mathcal L^N(t_1)\in B(\rho_{t_1},\epsilon) \right) = \delta_{P_*}, \end{equation*}$$

where $P_*$ is the minimizer of $D_{\mathrm{KL}}(P|R)$ over all probability laws P whose two endpoint marginal distributions are given by $\rho_{t_0}$ and $\rho_{t_1}$. The time t marginal of $P_*$ reduces to the solution ρ_t of Schrödinger's problem in equation (5).

As an analogue of Benamou–Brenier formula in classical Optimal Transport theory, the previous relative entropy minimization problem, with Hamiltonian $\mathcal H$ given by (1) for example, is equivalent to minimize the following average action [3],

$$\begin{equation} \int_{t_0}^{t_1} \int_{\mathbf{R}^n} \frac{|v(t,x)|^2}{2} \rho(t,\mathrm dx) \mathrm dt \end{equation} \tag{ 14 }$$

among all pairs $(\rho,v)$ solving the Fokker–Planck (FP) equation,

$$\begin{equation*} \left\{ \begin{aligned} &\frac{\partial}{\partial t} \rho + \mathrm{div}\, \left(\rho v \right) - \frac{1}{2} \Delta \rho = 0, \\ &\rho(t_0) = \rho_{t_0}, \ \rho(t_1) = \rho_{t_1}. \end{aligned} \right. \end{equation*}$$

The minimizer of (14) is the pair $(\rho,\nabla S)$ where S solves the following Hamilton–Jacobi–Bellman (HJB) equation,

$$\begin{equation}\left\{ \begin{aligned} & \frac{\partial S}{\partial t} + \frac{1}{2} |\nabla S|^2 + \frac{1}{2} \Delta S = 0, \\ & S(t_1) = \log\phi. \end{aligned} \right. \end{equation} \tag{ 15 }$$

In fact, this HJB equation is related to the heat equation (3) via the Cole–Hopf transformation (also called logarithmic transformation in stochastic optimal control [25]):

$$\begin{equation} S = \log w. \end{equation} \tag{ 16 }$$

Suppose that we study a system that requires n + 1 thermodynamic coordinates (i.e. its thermodynamic space is n + 1-dimensional), consisting of n work coordinates $(X^i)$ and the internal energy $E = :X^0$. The variables $(X^i)$ are extensive quantities that we can observe and control mechanically (and electromagnetically) in a macroscopic way. Each Xⁱ can be associated to a intensive quantity Y_i conjugate to it (with respect to the internal energy). These constitute a system of n conjugate pairs $(Y_i, X^i)$, $i = 1,\ldots, n$. The most commonly considered conjugate pairs are: pressure/volume pair $(-P, V)$, chemical potential/particle number pair $(\mu, N)$, magnetic field/moment pair $(\mathbf{B}, \mathbf{M})$, etc. The temperature/entropy pair (T, S) is not included here, as the entropy S is not considered as a work coordinate.

The first law of thermodynamics is essentially the conservation of energy of the system. The change of the internal energy E cannot be explained solely in terms of the total work W supplied to the system, and the discrepancy Q is understood as the energy transferred in the form of heat to the system:

$$\begin{equation*} \mathrm dE= \delta Q + W. \end{equation*}$$

If the system is closed and isolated, the first law of thermodynamics yields the conservation of internal energy E of the system.

The fundamental thermodynamic relation of Gibbs expresses E as a function of S and work coordinates $(X^i)$. Its differential form is an infinitesimal version of the first law of thermodynamics for quasistatic processes:

$$\begin{equation*} \mathrm d E=T \mathrm d S + Y_i \mathrm dX^i. \end{equation*}$$

Here and after, we adopt Einstein's summation convention. The entropic form of Gibbs relation reads

$$\begin{equation*} \mathrm d S=\frac{1}{T} \mathrm d E - \frac{Y_i}{T} \mathrm dX^i. \end{equation*}$$

We define the conjugate variables F₀ of E and F_i of Xⁱ , with respect to entropy, by

$$\begin{equation*} F_0 = \frac{1}{T}, \quad F_i = -\frac{Y_i}{T}. \end{equation*}$$

Then Gibbs relation reduces to

$$\begin{equation*} \mathrm d S = F_\alpha \mathrm d X^\alpha, \end{equation*}$$

where we use the Greek alphabet $(\alpha,\beta,\ldots)$ to indicate the summation from 0 to n, to differ from the Latin alphabet $(i,j,\ldots)$ which ranges from 1 to n. From the chain rule it follows that

$$\begin{equation*} F_\alpha = \frac{\partial S}{\partial X^\alpha}. \end{equation*}$$

The entropy S will be a maximum when the system is in an equilibrium state $\bar X = (\bar X^\alpha)$. Any fluctuation around the equilibrium state must cause the entropy to decrease. We let $(\hat X^\alpha)$ denote the fluctuations

$$\begin{equation*} \hat X^\alpha = X^\alpha - \bar X^\alpha. \end{equation*}$$

Assuming the system is close to its equilibrium, the fluctuations are small and we can expand the entropy around its equilibrium value to obtain

$$\begin{equation} S(X) \approx S(\bar X) + \frac{1}{2} A_{\alpha\beta} \hat X^\alpha \hat X^\beta \end{equation} \tag{ 17 }$$

where

$$\begin{equation*} A_{\alpha\beta} = \frac{\partial^2 S}{\partial X^\alpha \partial X^\beta}\bigg|_{X= \bar X} = \frac{\partial^2 S}{\partial \hat X^\alpha \partial \hat X^\beta}\bigg|_{\hat X = 0}. \end{equation*}$$

The matrix $(A_{\alpha\beta})$ is symmetric due to Maxwell relations, and it is negative semi-definite since the entropy is a concave function of all the thermodynamic coordinates. Equation (17) contains no first-order terms in X to ensure that spontaneous fluctuations about the equilibrium do not cause an increase in the entropy.

From chain rule it follows that

$$\begin{equation} F_\alpha = \frac{\partial S}{\partial X^\alpha} = \frac{\partial S}{\partial \hat X^\alpha} = A_{\alpha\beta} \hat X^\beta. \end{equation} \tag{ 18 }$$

We call $(F_\alpha)$ the conjugate variables of the fluctuation variables $(\hat X^\alpha)$.

The generalized flux J is defined as the time derivative of conditional average

$$\begin{equation*} J^\alpha:= \frac{\mathrm d \langle \hat X^\alpha \rangle_0}{\mathrm dt}, \end{equation*}$$

where $\langle \cdot \rangle_0$ stands for the average given the initial value $\hat X^\alpha(0)$. Then the rate of entropy production due to fluctuations is

$$\begin{equation} \frac{\mathrm d\langle S \rangle_0}{\mathrm dt} \stackrel{\times}{=} \left\langle \frac{\partial S}{\partial \hat X^\alpha} \right\rangle_0 \frac{\mathrm d \langle \hat X^\alpha \rangle_0}{\mathrm dt} = \langle F_\alpha \rangle_0 J^\alpha. \end{equation} \tag{ 19 }$$

Remark 3.1. From a modern point of view, Onsager's original derivations are partially incorrect, as indicated by the symbol × in equation (19). Indeed, stochastic calculus was still unknown in Onsager's time. Since $\hat X$ describe fluctuations, it is a stochastic process in general, modelled by a stochastic differential equation (SDE), see section 3.3. So in the calculus involving $\hat X$, Itô's formula must be applied. For example, in the derivation (19) of entropy production rate, we need first apply Itô' formula to get the stochastic differential of S,

$$\begin{equation*} \mathrm d S = \frac{\partial S}{\partial \hat X^\alpha} \mathrm d\hat X^\alpha + \frac{1}{2}\frac{\partial^2 S}{\partial \hat X^\alpha\partial \hat X^\beta} \mathrm d\hat X^\alpha \cdot \mathrm d\hat X^\beta, \end{equation*}$$

where the term $\frac{1}{2}\frac{\partial^2 S}{\partial \hat X^\alpha\partial \hat X^\beta} \mathrm d\hat X^\alpha \cdot \mathrm d\hat X^\beta$, additional to (19), is called Itô's correction.

Boltzmann's entropy formula provides

$$\begin{equation} S(X) = k_\mathrm B \log \Omega(X), \end{equation} \tag{ 20 }$$

where $\Omega(X)$ is the number of microstates compatible with the macrostate $X = (X^\alpha)$ in equilibrium. It is clear that the probability for X to have a value in the volume element $\mathrm dX$ of the thermodynamical space is proportional the number function $\Omega(X)$, i.e. $P(X) \propto \Omega (X) = e^{S(X)/K_\mathrm B}$. Indeed, we can now substitute (20) into (17) and obtain the following expression for the probability distribution of fluctuations around an equilibrium state,

$$\begin{equation} \begin{aligned} P(\hat X) &= \frac{\Omega(X)}{\int\Omega(X) \mathrm dX} = \frac{1}{\int e^{S(X)/k_\mathrm B} \mathrm dX} \exp\left( \frac{S(X)}{k_\mathrm B} \right) \\ &= \sqrt{\frac{\det(A)}{(2\pi k_\mathrm B)^n}} \exp\left( \frac{1}{2k_\mathrm B} A_{\alpha\beta} \hat X^\alpha \hat X^\beta \right). \end{aligned} \end{equation} \tag{ 21 }$$

By Laplace transform, one finds the covariance of $\hat X$ as follows,

$$\begin{equation} \langle \hat X^\alpha \hat X^\beta \rangle = -k_\mathrm B A^{\alpha\beta}, \end{equation} \tag{ 22 }$$

where $(A^{\alpha\beta})$ is the inverse matrix of $(A_{\alpha\beta})$.

In Onsager's seminal papers [13, 14], he assumed that,

'$\ldots$ the average regression of fluctuations will obey the same laws as the corresponding macroscopic irreversible process.'

Nowadays this is known as Onsager's regression hypothesis. It means that if the system is in a state fairly close to its equilibrium, i.e. if $\hat X$ may be regarded as small, the conditional average fluctuation $\langle \hat X \rangle_0$ may be expanded as the linear sums

$$\begin{equation} \frac{\mathrm d\langle \hat X^\alpha \rangle_0}{\mathrm dt} \approx -\lambda_\beta^\alpha \langle \hat X^\beta \rangle_0, \end{equation} \tag{ 23 }$$

where $\lambda = (\lambda_\beta^\alpha)$ is an $(n+1)\times(n+1)$ constant matrix. Its solution has the following short-time approximation:

$$\begin{equation} \langle \hat X^\alpha(t) \rangle_0 = \hat X^\alpha(0) - t \lambda_\beta^\alpha \hat X^\beta(0) + \mathrm{O}(t^2). \end{equation} \tag{ 24 }$$

The principle of microscopic reversibility states that the microscopic dynamics of particles and fields is time-reversible, like Newton's equation. According to Boltzmann, this microscopic reversibility implies a principle of detailed balance for collisions. As a consequence, the correlation functions for macroscopic fluctuations $\hat X$ obey the relations

$$\begin{equation} \langle \hat X^\alpha(0) \hat X^\beta(t) \rangle = \langle \hat X^\alpha(t) \hat X^\beta(0) \rangle. \end{equation} \tag{ 25 }$$

Substitute (24) into (25), and notice that $\langle \hat X^\alpha(0) \hat X^\beta(t) \rangle = \langle \hat X^\alpha(0) \langle \hat X^\beta(t) \rangle_0 \rangle$. We obtain

$$\begin{equation*} \lambda_\gamma^\beta \langle \hat X^\alpha(0) \hat X^\gamma(0) \rangle = \lambda_\gamma^\alpha \langle \hat X^\gamma(0) \hat X^\beta(0) \rangle. \end{equation*}$$

Then we apply equation (22) to get

$$\begin{equation*} \lambda_\gamma^\beta A^{\alpha\gamma} = \lambda_\gamma^\alpha A^{\beta\gamma}. \end{equation*}$$

Now we can define a new matrix $L^{\alpha\beta} = -\lambda_\gamma^\beta A^{\alpha\gamma}$ and obtain the following Onsager's reciprocal relations:

$$\begin{equation*} L^{\alpha\beta} = L^{\beta\alpha}. \end{equation*}$$

The matrix $(L^{\alpha \beta})$ is called the Onsager matrix.

If we make use of the generalized force F in equation (18), the time rate of change of the fluctuation can be written

$$\begin{equation} J^\alpha = \frac{\mathrm d\langle \hat X^\alpha \rangle_0}{\mathrm dt} = -\lambda_\beta^\alpha \langle \hat X^\beta \rangle_0 = L^{\alpha\gamma} A_{\gamma\beta} \langle \hat X^\beta \rangle_0 = L^{\alpha\gamma} \langle F_\gamma \rangle_0 = L^{\alpha\gamma} \left\langle \frac{\partial S}{\partial \hat X^\gamma} \right\rangle_0. \end{equation} \tag{ 26 }$$

Substituting (26) in (19), we obtain the following expression for the entropy production,

$$\begin{equation} \frac{\mathrm d\langle S \rangle_0}{\mathrm dt} \stackrel{\times}{=} \langle F_\alpha \rangle_0 J^\alpha = L^{\alpha\beta} \langle F_\alpha \rangle_0 \langle F_\beta \rangle_0. \end{equation} \tag{ 27 }$$

In out-of-equilibrium conditions the right-hand side of this equation has to be a positive quantity, which vanishes at equilibrium. So the Onsager matrix $(L^{\alpha \beta})$ is positive definite.

In view of equation (26), let us write the following Langevin-type equation

$$\begin{equation} \frac{\mathrm d \hat X^\alpha}{\mathrm dt} = L^{\alpha\beta} \frac{\partial S}{\partial \hat X^\beta} + \xi, \end{equation} \tag{ 28 }$$

where $\xi = (\xi^\alpha)$ is a white noise with correlation function

$$\begin{equation} \langle\xi^\alpha(t) \xi^\beta(s)\rangle = 2 D^{\alpha\beta} \delta(t-s), \end{equation} \tag{ 29 }$$

and $D = (D^{\alpha\beta})$ is a constant matrix describing the magnitude of the fluctuation.

We know, since Einstein's work on Brownian motion, that D cannot be chosen arbitrarily for a given system, if it is in equilibrium. The noise amplitude D and the dissipation rate L are related.

The (time-independent) FP equation for equation (28) in equilibrium reads (the drift is $L\nabla S$ and the diffusion matrix is 2D),

$$\begin{equation} D^{\alpha\beta} \frac{\partial^2 \rho}{\partial x^\alpha \partial x^\beta} - L^{\alpha\beta} \frac{\partial}{\partial x^\alpha} \left( \rho \frac{\partial S}{\partial x^\beta} \right) = 0. \end{equation} \tag{ 30 }$$

Since (21) is a solution of equation (30), we deduce

$$\begin{equation} D^{\alpha\beta} = k_\mathrm B L^{\alpha\beta}. \end{equation} \tag{ 31 }$$

This relation is a generalization of the fluctuation-dissipation relation in classical Brownian motion theory.

Remark 3.2. For the same reason as remark 3.1, the calculation of equation (27) is not completely correct. On the other hand, the close-to-equilibrium approximations used in (17) and (23) are important to achieve equation (28). All of these show that Onsager's original derivations are a bit defective. What we are going to do in the upcoming section 4 is to reformulate Onsager's idea from Schrödinger's viewpoint and fix those problems.

From a modern point of view, what Onsager is telling us is that mesoscopic observables obey the large deviation principle around the macroscopic law due to the law of large numbers:

$$\begin{equation*} \mathbf{P}\left[ \frac{1}{\delta t} \int_0^{\delta t} \left(\frac{\mathrm d X(s)}{\mathrm d s}- J(t) \right) \mathrm ds \in \mathrm dx \right] \approx 4d \exp\left(- \frac{x^T D^{-1} x}{4d} \delta t \right), \end{equation*}$$

where d is the spatial dimension.

From (28), (29) and (31), as well as Boltzmann's formula (20), we know that the infinitesimal generator of the fluctuation process $\hat X$ is

$$\begin{equation*} D^{\alpha\beta} \left[ \frac{\partial^2}{\partial x^\alpha\partial x^\beta} + \frac{\partial (\log\Omega)}{\partial x^\beta} \frac{\partial}{\partial x^\alpha} \right], \end{equation*}$$

where Ω is the number function in equation (20). This generator is exactly of the form (12). This suggests that Onsager's way to approach nonequilibrium statistical physics can be compatible with Schrödinger's thought experiment.

Indeed, we can replace the Brownian particles of Schrödinger's problem by a thermodynamic system close to its equilibrium state, and the position variables of the particle by the deviation $(\hat X^\alpha)$ of thermodynamic coordinates from the equilibrium. From time t₀, the observables $(\hat X^\alpha)$ start to spontaneously evolve, with an initial distribution near to the equilibrium distribution (21). At a sufficiently large time t₁, the system reaches the equilibrium state and the observables $(\hat X^\alpha)$ end up with the equilibrium distribution (21).

Since the time scale of t₁ is very large, mathematically we may treat it as infinity. The heat equation (3) and its associated HJB equation (15) becomes approximately time-independent, so the Cole–Hopf transformation (16) has the same form as Boltzmann's entropy formula (20). This suggests that the entropy S, as a function of fluctuation variables $(\hat X^\alpha)$, should satisfy the HJB equation.

Moreover, the fact that the Onsager matrix $(L^{\alpha \beta})$ is symmetric and positive definite suggests that one can treat $(L^{\alpha \beta})$ as a Riemannian metric on the thermodynamic space of fluctuation variables $(\hat X^\alpha)$. The conjugation between $(\hat X^\alpha)$ and $(F_\alpha)$ and their relation with S in (18) indicate that Legendre transform should underlie the framework.

In view of these observations and remarks 3.1 and 3.2, it is possible to connect Schrödinger's problem with Onsager's theory in systematical way, compatible with Itô's stochastic calculus. In our paper [30], we have developed a stochastic framework of geometric mechanics, based on the theory of SO differential geometry due to Schwartz and Meyer [31, 32]. These geometric theories are indeed motivated by and compatible with the stochastic nature of Brownian randomness and Itô's calculus.

Consider a general m-dimensional manifold M, which can be the flat Euclidean spaces, or the thermodynamic space. The coordinates on M are denoted by $(x^\alpha)$.

The SO tangent bundle $\mathcal T^S M$ is a fiber bundle such that its sections are SO operators of the following form

$$\begin{equation} A = A^\alpha \frac{\partial}{\partial x^\alpha} + A^{\alpha\beta} \frac{\partial^2}{\partial x^\alpha\partial x^\beta}. \end{equation} \tag{ 32 }$$

A natural frame of it is the following set of elementary differential operators of first-order and SO:

$$\begin{equation} \left\{\frac{\partial}{\partial{x^\alpha}}, \frac{\partial^2}{\partial x^\alpha \partial x^\beta}: 1\leqslant \alpha\leqslant \beta \leqslant m \right\}. \end{equation} \tag{ 33 }$$

The dual bundle is the SO cotangent bundle, denoted by $\mathcal T^{S*} M$, of which the frame dual to (33) is

$$\begin{equation} \left\{\mathrm d^2 x^\alpha, \textstyle{\frac{1}{2}} \mathrm dx^\alpha\cdot \mathrm dx^\alpha, \mathrm dx^\alpha\cdot \mathrm dx^\beta: 1\leqslant \alpha< \beta \leqslant m \right\}, \end{equation} \tag{ 34 }$$

where d is the classical differential operator, $\mathrm d^2$ is called SO differential operator and · is called symmetric product. They acts on smooth functions $f, g$ on M as follows,

$$\begin{equation*} \langle \mathrm d^2 f, A \rangle = \frac{\partial f}{\partial x^\alpha} A^\alpha + \frac{\partial^2 f}{\partial x^\alpha \partial x^\beta} A^{\alpha\beta}, \quad \langle \mathrm df\cdot \mathrm dg, A \rangle = \frac{\partial f}{\partial x^\alpha} \frac{\partial g}{\partial x^\beta} A^{\alpha\beta}, \end{equation*}$$

where A is the operator (32). From equation (34) and after, we adopt the convention that $\mathrm dx^\beta\cdot \mathrm dx^\alpha = \mathrm dx^\alpha\cdot \mathrm dx^\beta$ for all $1\leqslant \alpha\lt \beta \leqslant m$. The factor $\frac{1}{2}$ in front of in $\mathrm dx^\alpha\cdot \mathrm dx^\alpha$ in (34) is a normalized constant since the off-diagonal elements $\mathrm dx^\alpha\cdot \mathrm dx^\beta$ are folded to the lower triangular ones with $\alpha\lt \beta$. The coordinates $(x^\alpha)$ on M induce a canonical coordinate system on $\mathcal T^{S*} M$, denoted by $(x^\alpha, p_\alpha, o_{\alpha\beta})$, where extra variables $(o_{\alpha\beta})$ describe SO effects. Sections of $\mathcal T^{S*} M$ are called SO forms.

As indicated by the notation, the expression $\mathrm d^2 f$ can indeed be understood as $\mathrm d(\mathrm d f)$, where the differential $\mathrm d$ inside the parentheses is de Rham's exterior differential on M, while the one outside is the exterior differential on TM by regarding the first differential $\mathrm df$ as a function on TM [34]. The bundle $\mathcal T^{S*} M$ is not an unacquainted object, since $\mathcal T^{S*} M\times\mathbf{R}$ is bundle diffeomorphic to the classical SO jet bundle of the trivial bundle $(M\times \mathbf{R}, \pi, M)$. The bundle $\mathcal T^S M$ can also be endowed with a jet-like structure by introducing the notion of 'stochastic jets' for diffusion processes, as opposed to the higher-order tangent bundle of higher-order jets for deterministic curves introduced in [35]. The reason why we make use of the two odd-looking bundles is that, according to Schwartz and Meyer's heuristic principle, any geometric statement for such SO (co)tangent vectors will have a probabilistic content. See [30] for more details.

Let be given a probability space $(\Omega,\mathcal F,\mathbf{P})$ and a nondecreasing filtration $\mathcal P_t\subset \mathcal F$. Roughly speaking, $\mathcal F$ contains all the probabilistic information while $\mathcal P_t$ only contains all past information before the present time t. For a diffusion process $(x^\alpha(t))$ valued on M, the coefficients of its generator can be characterized by

$$\begin{equation} (D^\alpha x)(t) := \lim_{\epsilon\to0^+} \mathbf{E} \textstyle{{\left[ \frac{x^\alpha(t+\epsilon)-x^\alpha(t) }{\epsilon} \Big| \mathcal P_t \right]}}, \end{equation} \tag{ 35 }$$

$$\begin{equation*} (Q^{\alpha\beta} x)(t) := \lim_{\epsilon\to0^+} \mathbf{E} \textstyle{{\left[ \frac{(x^\alpha(t+\epsilon)-x^\alpha(t)) (x^\beta(t+\epsilon)-x^\beta(t))}{\epsilon} \Big| \mathcal P_t \right]}}. \end{equation*}$$

The pair $(D^\alpha x(t),Q^{\alpha\beta} x(t))$ is a process taking values in $\mathcal T^S M$, and called the mean derivatives of $(x^\alpha(t))$. The mean derivative of a function f acting on x(t) can be defined similarly and expressed by

$$\begin{equation} \mathbf{D}_t \,f := \lim_{\epsilon\to0^+} \mathbf{E} \textstyle{{\left[ \frac{f(x(t+\epsilon))-f(x(t)) }{\epsilon} \Big| \mathcal P_t \right]}} = \frac{\partial f}{\partial t} + D^i x \frac{\partial f}{\partial x^i} + \frac{1}{2} Q^{\,jk} x \frac{\partial^2 f}{\partial x^j \partial x^k}. \end{equation} \tag{ 36 }$$

In general, $(D^\alpha x)$ does not transform as a vector field, which can be verified by applying Itô's formula for changes of coordinates. In order to overcome this problem, we equip M with a linear connection $\nabla$, and use it to compensate the correction term resulting from Itô's formula. That is, we define the following $\nabla$-dependent mean derivative, in forms of Christoffel's symbols,

$$\begin{equation*} D_\nabla^\alpha x = D^\alpha x + \textstyle{{\frac{1}{2}}} \Gamma^\alpha_{\beta\gamma} Q^{\beta\gamma}x. \end{equation*}$$

Then $D_\nabla X$ does transform as a vector field.

The SO Poincaré form θ is a SO form on $\mathcal T^{S*} M$, given by

$$\begin{equation*} \theta = p_\alpha \mathrm d^2 x^\alpha + \textstyle{\frac{1}{2}} o_{\alpha\beta} \mathrm dx^\alpha\cdot \mathrm dx^\beta. \end{equation*}$$

And the canonical SO symplectic form becomes

$$\begin{equation*} \omega = \mathrm d^2 \theta = \mathrm d^2 x^\alpha \wedge \mathrm d^2 p_\alpha + \textstyle{\frac{1}{2}} \mathrm dx^\alpha\cdot \mathrm dx^\beta \wedge \mathrm d^2 o_{\alpha\beta}. \end{equation*}$$

See [30, section 6.2] for the definition of the action of the above SO form ω on SO operators of the form (32). It can be verified that ω is nondegenerate in the sense that $\omega(A,B) = 0$ for all SO operators B implies A = 0. We call the pair $(\mathcal T^{S*} M, \omega)$ a SO symplectic manifold. A SO Hamiltonian is a smooth function $H: \mathcal T^{S*} M\times\mathbf{R} \to \mathbf{R}$. Its SO Hamiltonian vector field is a SO vector field on $\mathcal T^{S*} M$ given by

$$\begin{equation*} A_H \approx \frac{\partial H}{\partial p_\alpha} \frac{\partial}{\partial{x^\alpha}} - \frac{\partial H}{\partial x^\alpha} \frac{\partial}{\partial{p_\alpha}} + \frac{\partial H}{\partial o_{\alpha\beta}} \frac{\partial^2}{\partial x^\alpha \partial x^\beta} - \frac{\partial^2 H}{\partial x^\alpha \partial x^\beta} \frac{\partial}{\partial{o_{\alpha\beta}}}, \end{equation*}$$

where ≈ stands for principal parts. The stochastic Hamilton's equations are the following SDEs on $\mathcal T^{S*} M$,

$$\begin{equation} \left\{ \begin{aligned} D^\alpha x &= \frac{\partial H}{\partial p_\alpha}, \\ Q^{\alpha\beta} x &= 2\frac{\partial H}{\partial o_{\alpha\beta}}, \\ D_\alpha p &= - \frac{\partial H}{\partial x^\alpha}, \\ D_{\alpha\beta} o &\approx - \frac{\partial^2 H}{\partial x^\alpha \partial x^\beta}. \\ \end{aligned} \right. \end{equation} \tag{ 37 }$$

The above equations are not necessarily solvable. But if we set $p_\alpha = p_\alpha(t,x)$ and $o_{\alpha\beta} = o_{\alpha\beta}(t,x)$, they can be simplified. Indeed, the first two equations give the generator of $(x^\alpha(t))$, which can be used to apply Itô's formula to the last two equations. It follows that

$$\begin{equation*} o_{\alpha\beta} = \frac{\partial p_\alpha}{\partial x^\beta} = \frac{\partial p_\beta}{\partial x^\alpha}. \end{equation*}$$

These are known as Maxwell's relations in thermodynamics (also called Onsager's relations somewhat inaccurately in [15]).

A change of coordinates leaving the form of the stochastic Hamilton's equations (37) unchanged give rise to the following SO HJB equation,

$$\begin{equation} \frac{\partial S}{\partial t} + H(\mathrm d^2 S,t) = \frac{\partial S}{\partial t} + H\left( x^\alpha, \frac{\partial S}{\partial x^\alpha}, \frac{\partial^2 S}{\partial x^\alpha \partial x^\beta}, t \right) = 0. \end{equation} \tag{ 38 }$$

This HJB equation can also be related to the stochastic Hamilton's equations (37) via the following statement: for any process x(t) satisfying the first two equations of (37), the process

$$\begin{equation*} (x^\alpha(t), p_\alpha(t), o_{\alpha\beta}(t) ) = \left( x^\alpha(t), \frac{\partial S}{\partial x^\alpha}(t, x(t)), \frac{\partial^2 S}{\partial x^\alpha\partial x^\beta}(t, x(t)) \right) \end{equation*}$$

solves (37) if and only if S solves (38).

There is no unique way to obtain a SO Hamiltonian from a given classical Hamiltonian $H_0: T^* M\times\mathbf{R} \to \mathbf{R}$. But if we specify a Riemannian metric g, then there is a canonical way, that is,

$$\begin{equation} H(x,p,o,t) = H_0(x,p,t) + \textstyle{\frac{1}{2}} g^{\alpha\beta}(x) \left(o_{\alpha\beta} - \Gamma^\gamma_{\alpha\beta}(x)p_\gamma \right). \end{equation} \tag{ 39 }$$

In this case, the stochastic Hamilton's equations can be rewritten to the following global form on $T^* M$,

$$\begin{equation} \left\{ \begin{aligned} D^{}_\nabla x = \nabla_p H_0, \\ \frac{\overline{\mathbf{D}}}{\mathrm dt} p = -\mathrm d_x H_0, \end{aligned} \right. \end{equation} \tag{ 40 }$$

subject to $Q^{\alpha\beta}x(t) = g^{\alpha\beta}(x(t))$, where $\frac{\overline{\mathbf{D}}}{\mathrm dt} = \frac{\partial}{\partial t} + \nabla_{D_\nabla X} + \frac{1}{2} \Delta_{\mathrm{LD}}$. The corresponding HJB equation for this case is

$$\begin{equation} \frac{\partial S}{\partial t} + H_0(\mathrm d S, t) + \frac{1}{2}\Delta S = 0, \end{equation} \tag{ 41 }$$

where the Laplacian term express clearly the stochastic deformation of classical Hamilton–Jacobi equation.

The stochastic Lagrangian mechanics is established by considering the following stochastic variational problem: to minimize the action functional

$$\begin{equation} \mathcal S[x] = \mathbf{E} \int_{t_0}^{t_1} L_0\left(t, x(t), D_\nabla^{} x(t) \right) \mathrm dt \end{equation} \tag{ 42 }$$

over all diffusion processes $(x(t):t_0\leqslant t\leqslant t_1)$, with initial and final marginal distributions µ₀ and µ₁, and satisfying $Q^{\alpha\beta}x(t) = g^{\alpha\beta}(x(t))$. Solving the stochastic variational problem via the Cameron–Martin type variation leads to the following stochastic Euler–Lagrange equation

$$\begin{equation} \frac{\overline{\mathbf{D}}}{\mathrm dt} \left( d_{\dot x} L_0 \right) = d_x L_0. \end{equation} \tag{ 43 }$$

It turns out that the stochastic Euler–Lagrange equation is equivalent, as it should, to the global version of stochastic Hamilton's equation (40) via the Legendre transform

$$\begin{equation} D_\nabla^\alpha x = \frac{\partial H_0}{\partial p_\alpha}, \quad L_0(t,x,D_\nabla^{} x) = p_\alpha D_\nabla^\alpha x - H_0(x,p,t). \end{equation} \tag{ 44 }$$

For the canonical SO Hamiltonian induced by H₀, as in (39), we have

$$\begin{equation*} H(x,p,o,t) = p_\alpha D^\alpha x + \frac{1}{2} o_{\alpha\beta} Q^{\alpha\beta} x - L_0(t,x,D_\nabla^{} x). \end{equation*}$$

Then we apply the mean derivative (36) to the solution S of the special HJB equation (41),

$$\begin{equation} \mathbf{D}_t S = \frac{\partial S}{\partial t} + D^\alpha x \frac{\partial S}{\partial x^\alpha} + \frac{1}{2} Q^{\alpha\beta} x \frac{\partial^2 S}{\partial x^\alpha \partial x^\beta} = - H + p_\alpha D^\alpha x + \frac{1}{2} o_{\alpha\beta} Q^{\alpha\beta} x = L_0. \end{equation} \tag{ 45 }$$

Hence, the action functional (42) can be expressed as

$$\begin{equation} \mathcal S[x] = \mathbf{E} \left[ S(t_1,x(t_1)) - S(t_0,x(t_0)) \right]. \end{equation} \tag{ 46 }$$

Now we take M to be the thermodynamic space of fluctuation variables, and equip it with a Riemannian metric $(L^{\alpha\beta})$. Consider the following SO Hamiltonian

$$\begin{equation} H(x,p,o) = \frac{1}{2} L^{\alpha\beta}(x) p_\alpha p_\beta + k_\mathrm B L^{\alpha\beta}(x) \left( o_{\alpha\beta} - \Gamma_{\alpha\beta}^\gamma(x) p_\gamma \right), \end{equation} \tag{ 47 }$$

which is canonically induced by $H_0 = \frac{1}{2} L^{\alpha\beta}(x) p_\alpha p_\beta$. We first solve the last two partial differential equations of (37) by considering the following backward heat equation:

$$\begin{equation*} \frac{\partial w}{\partial t} + k_\mathrm B \Delta w = 0, \end{equation*}$$

where the Laplacian Δ is under the Riemannian metric $(L^{\alpha\beta})$. Then the Cole–Hopf transformation

$$\begin{equation} S(t,x) = 2k_\mathrm B \log w(t,x) \end{equation} \tag{ 48 }$$

solves the following HJB equation

$$\begin{equation*} \frac{\partial S}{\partial t} + \frac{1}{2} |\nabla S|^2 + k_\mathrm B \Delta S = 0, \end{equation*}$$

where the gradient $\nabla$ is also with respect to $(L^{\alpha\beta})$. It is easy to verify that

$$\begin{equation} p_\alpha = \frac{\partial S}{\partial x^\alpha}, \quad o_{\alpha\beta} = \frac{\partial^2 S}{\partial x^\alpha\partial x^\beta} \end{equation} \tag{ 49 }$$

solve the last two equations. Therefore, the first two mean differential equations reduce to

$$\begin{equation} \frac{\mathrm d x^\alpha}{\mathrm dt} = L^{\alpha\beta} \frac{\partial S}{\partial x^\beta} + \xi, \end{equation} \tag{ 50 }$$

where the white noise ξ has correlation

$$\begin{equation*} \langle\xi^\alpha(t) \xi^\beta(s)\rangle = 2k_\mathrm B L^{\alpha\beta} \delta(t-s). \end{equation*}$$

Comparing (50) with (28), we see that the Riemannian metric $(L^{\alpha\beta})$ is a generalization of the Onsager matrix to inhomogeneous transport processes (as $(L^{\alpha\beta})$ can vary for different thermodynamic variables), and the function S in (48) provides a definition of entropy for nonequilibrium systems. The relation (46) between action functional $\mathcal S$ and entropy S suggests that the quantity $S(t,x(t))$ can be regarded as a trajectory-dependent entropy [36] and $\mathcal S$ is the totality of change of entropy. The relations (49) and (18) indicate that the generalized force $(F_\alpha)$ is just the conjugate momentum $(p_\alpha)$.

It follows from the Legendre transform (44) that $D_\nabla^\alpha x = L_{\alpha\beta} p_\alpha$ and the Lagrangian associated to H of (47) is given by

$$\begin{equation*} L_0(x, D^{}_\nabla x) = \frac{1}{2} L_{\alpha\beta}(x) D^\alpha_\nabla x D^\beta_\nabla x = \frac{1}{2} L^{\alpha\beta}(x) p_\alpha p_\beta. \end{equation*}$$

This lead to the following entropy production formula, by (45):

$$\begin{equation*} \frac{\mathrm d\langle S\rangle}{\mathrm dt} = \langle L_0\rangle = \frac{1}{2} \langle L^{\alpha\beta} p_\alpha p_\beta \rangle, \end{equation*}$$

which recovers Onsager's formula (27).

Let us consider a Langevin equation of a form more general than (28), but now interpreted in a phase space,

$$\begin{equation} \left\{\begin{aligned} \mathrm dx &= v \mathrm dt, \\ \mathrm dv &= -\nabla U(x)\mathrm dt - \gamma v \mathrm dt + \sqrt{2\gamma k_\mathrm B T} \mathrm dW, \end{aligned} \right. \end{equation} \tag{ 51 }$$

where $U:\mathbf{R}^n\to \mathbf{R}$ is a smooth potential such that $e^{-U(x)} \in L^1(\mathbf{R}^n)$, γ a positive constant describing the damping strength, W is a standard n-dimensional Brownian motion. Written formally as an ordinary differential equation, it is generally regarded as a Newton equation with a dissipative term proportional to the velocity, and a source of noise

$$\begin{equation*} \ddot x = -\nabla U(x) - \gamma \dot x + \sqrt{2\gamma k_\mathrm B T} \dot W, \end{equation*}$$

whose 'acceleration', given definition (51), is particularly problematic (the last term of r.h.s. being the white noise). The process (x, v) is Markovian and hypoelliptic (there is no noise in the first equation of (51)). Moreover, this process is ergodic, with a (Gibbs) product invariant measure

$$\begin{equation*} \mu(\mathrm dx,\mathrm dv) = Z^{-1} \exp\left\{-\frac{H_0(x,v)}{k_\mathrm B T}\right\} \mathrm dx\mathrm dv, \end{equation*}$$

where H₀ is the Hamiltonian $\frac{1}{2}|v|^2 + U(x)$, and whose normalization constant Z is the partition function on $\mathbf{R}^{2n}$. The process (x, v) is non-reversible in Kolmogorov equilibrium sense [37] (which we shall recall below), i.e. its generator $\mathcal L$ is not self-adjoint in $L^2(M)$. The FP equation of the Langevin equation (51) takes the form $\frac{\partial p}{\partial t} = \mathcal L^* p$, for $\rho = \rho_t(x,v)$, with

$$\begin{equation*} \mathcal L^* \rho = -v \cdot \nabla_x \rho + \nabla_x U \cdot \nabla_v \rho + \gamma \nabla_v (v \rho) + \gamma k_\mathrm B T \Delta_v \rho, \end{equation*}$$

the existence and uniqueness of solution of this equation are known, as well as its rate of convergence to equilibrium, but the proofs are difficult [38]. Moreover, as can be expected, those equations cannot be solved explicitly, except for quadratic potentials U.

It is relevant here, to recall that this analytical notion of reversibility [37] was expressed originally for a single stationary diffusion process with invariant measure $\mu(\mathrm dx) = \rho(x)\mathrm dx$. The forward and backward transition probability densities, respectively $p(x,t,\mathrm dy)$ and $p_*(\mathrm dx,t,y)$, satisfy the following 'detailed balance condition':

$$\begin{equation} \rho(x) p(x,t,\mathrm dy) \mathrm dx = p_*(\mathrm dx,t,y) \rho(y) \mathrm dy, \end{equation} \tag{ 52 }$$

where p solves the usual FP equation

$$\begin{equation} \frac{\partial p}{\partial t} = -\frac{\partial}{\partial{y^i}} (b^i(y) p) + \frac{1}{2} \frac{\partial^2}{\partial y^i \partial y^j} (a^{ij}(y) p), \end{equation} \tag{ 53 }$$

with drift vector field $(b^i)$ and positive definite diffusion tensor $(a^{ij})$, and $p_*$ denotes the probability transition to come back around x from a future location y, during the time interval $[0,t]$. Kolmogorov formulated the reversibility as the following equation

$$\begin{equation*} p_*(\mathrm dx,t,y) = p(y,t,\mathrm dx), \end{equation*}$$

and proved that this equation holds if and only if the drift is a gradient:

$$\begin{equation} b^i = \frac{a^{ij}}{2} \frac{\partial (\log\rho)}{\partial x^j}. \end{equation} \tag{ 54 }$$

This stationary gradient form of the drift (54), involved in his equilibrium reversibility condition [37] was therefore generalized via forward and backward transition probabilities (6), (7) to associated drifts describing Schrödinger's time-inhomogeneous processes. Given an initial probability density ρ, the future one follows by integration w.r.t. the elementary solution of equation (53). This means that it is sufficient to solve an initial value FP problem for ρ. The detailed balance condition (52) expresses two equivalent ways to define the two-times joint probability of a reversible homogeneous Markov process.

Comparing relation (52) with (8), it is clear that Schrödinger's nonequilibrium notion of reversibility, needed for his variational problem, is the general form (but local in time) of Kolmogorov's one, formulated five years earlier. As said in the introduction, Schrödinger's context was, apparently, the one of classical statistical physics but his problem was essentially inspired by his doubts about the foundations of quantum mechanics. It is only when the unexpected relation between Nelson's stochastic mechanics [39] (a radical attempt to interpret the quantum wave equation itself as a local in space, regularized, version of Newtonian dynamics) and Schrödinger's 1931 idea was discovered [2] that Kolmogorov's conditions (52) and (54) were suitably generalized (see [2, 16]). The special feature of Schrödinger(-Bernstein) diffusions is that they admit simultaneously two transition probabilities, allowing to preserve his nonequilibrium notions of detailed balance and reversibility.

We shall give a theorem summarizing the key conditions for the construction of the diffusions in $\mathbf{R}^n$ solving Schrödinger's variational problem for the special class of Hamiltonian $\mathcal H$, as in (2)–(3), of the form

$$\begin{equation*} \mathcal H= \mathcal H_0 + U = -\frac{\hbar^2}{2} \Delta + U(x), \end{equation*}$$

with U an additional scalar potential. The positivity of the above free integral kernel h₀ is obvious so the one of $h(x, t-s, y) = (e^{-(t-s) \mathcal H} \delta_x)(y)$ is a condition on the potential U. The theorem involves results of a number of authors, along the years, following Schrödinger's seminal idea (see [18, 40, 41] for more).

Theorem 5.1. Let $\mathcal H = \mathcal H_0 +U$ densely defined in $L^2(\mathbf{R}^n)$, lower bounded self-adjoint Hamiltonian, with $U:\mathbf{R}^n\to\mathbf{R}$ a real measurable function in the Kato class $\mathcal K^n$ [40]. For instance, when n = 3, the harmonic oscillator and Coulomb potential are in $\mathcal K^3$. If $\max(-U,0)\in \mathcal K^n$ and $U \boldsymbol{1}_{B(0,L)} \in \mathcal K^n$, the kernel $h(x,t-s,y)$, $t\geqslant s$ is jointly continuous in the three variables and nonnegative. Let $\rho_{t_0}(x)$ and $\rho_{t_1}(y)$ be two strictly positive probability densities, for $\mathcal H$ as before. Then positive, not necessarily integrable, solution $(\psi,\phi)$ of equation (10) for h as before, exists and is unique. The optimal diffusion X(t) is, therefore, completely determined.

Schrödinger himself checked, in [17], that the simplest diffusion on $\mathbf{R}$ associated with $\mathcal H_0$ (with zero potential $U\equiv 0$), the one-dimensional Brownian motion $X(t) = x+ W(t-s)$, starting from x at time s was compatible with his approach, though as a degenerate, manifestly nonequilibrium example.

In this case we are given that $\rho_{t_0}(x) = \delta_x$ and $\rho_{t_1}(y) = [2\pi (t_1-t_0)]^{-\frac{1}{2}} e^{-\frac{|y-x|^2}{2(t_1-t_0)}}$ and the solution of equation (10) is the trivial one $\phi(y) = 1$, $\psi(x) = \delta_x$. The corresponding solutions of equations (2) and (3), $t\in[t_0,t_1]$ are $w(x,t) \equiv 1$ and $w^*(y,t) = [2\pi (t-t_0)]^{-\frac{1}{2}} e^{-\frac{|y-x|^2}{2(t-t_0)}}$. Using (6) and (7) it is easy to verify that the usual (forward) drift B is indeed the conditional expectation

$$\begin{equation*} B(x,t) = \lim_{\Delta t\to0^+} \mathbf E\left[ \frac{X(t+\Delta t)-X(t)}{\Delta t} \bigg| X(t) = x \right] = \frac{\nabla w}{w} (x,t) = 0, \end{equation*}$$

and the backward one

$$\begin{equation*} B_*(y,t) = \lim_{\Delta t\to0^+} \mathbf E\left[ \frac{X(t)-X(t-\Delta t)}{\Delta t} \bigg| X(t) = y \right] = -\frac{\nabla w^*}{w^*} (y,t) = \frac{y-x}{t-t_0}, \end{equation*}$$

whose singularity in $t = t_0$ expresses our singular initial condition $X(t_0) = x$. Also notice that φ, here, is not integrable. Of course, given $\mathcal H = \mathcal H_0 + U$, each new data of $(\rho_{t_0}, \rho_{t_1})$ as in theorem 5.1, will be associated to a new diffusion. This is the origin of the rich dynamical content of the stochastic theory described here.

In this elementary example, the Hamiltonian function of system (40) reduces, by (39) in this flat one-dimensional case, to

$$\begin{equation*} H(x,p,o) = H_0(x,p) + \frac{1}{2} o, \end{equation*}$$

where $o = \frac{\partial p}{\partial x}$ and $H_0 = \frac{1}{2} |p|^2$ is the free classical (but Euclidean) Hamiltonian. In other words, the logarithm substitution $S = \log w$, for w solving equation (3), provides HJB with final boundary condition (since $p = \nabla S$),

$$\begin{equation*} \frac{\partial S}{\partial t} + H_0(x,\nabla S) + \frac{1}{2} \Delta S = 0. \end{equation*}$$

However, this solution is trivial because of our choice $\phi = w = 1$ for equation (3). The backward HJB counterpart of (41), namely,

$$\begin{equation*} \frac{\partial S^*}{\partial t} + H_0(x,-\nabla S^*) - \frac{1}{2} \Delta S^* = 0 \end{equation*}$$

related with equation (2) by the transformation $S^* = -\log w^*$ and $p = \nabla S^*$, is more interesting since $w^*$ is the elementary solution of the Cauchy problem equation (2). The resulting drift is $B_*$ founded before.

In its traditional (forward) description, with zero drift, the Brownian motion looks hardly 'reversible'. However, as a reciprocal process, its backward drift displays its singular initial boundary condition, an observation known for a long time (see for instance, [42, example 8.3]), but rarely used in applied sciences.

The Lagrangian L₀ associated with H₀ is the free one $L_0(X(t), DX(t)) = \frac{1}{2}|DX(t)|^2$ which, for the Brownian motion, corresponds to a trivial solution of the stochastic Euler–Lagrange equation (43), namely, $DDX(t) = 0$. This example provides us the motivation to elaborate on the classical notion of time reversal and Schrödinger's one.

Traditionally, for any classical observable in phase space $O(x(t),p(t))$, a time reversal operator T_c is introduced so that $\mathcal T_c(x(t),p(t)) = (x(t),-p(t))$. Generally, the classical Hamiltonian observable $H_0(x,p)$ of the system is invariant under $\mathcal T_c$. Then the microscopic reversibility of the system means that if $(x(t),p(t))$ are solutions of Hamilton's equations with initial condition $(x(0),p(0))$, so are $\mathcal T_c(x(t),p(t))$, with time reversed boundary condition $\mathcal T_c(x(0),p(0))$.

For Bernstein's reciprocal process X(t) solving Schrödinger's problem, the notion is a bit more sophisticated, because of the non-differentiability of their paths. Let us summarize it for the special class of Hamiltonian operators $\mathcal H$ of theorem 5.1.

As already said, it follows from the special form (5) of the probability density looked for by Schrödinger that those process X(t) are themselves invariant under time reversal. Regarding the drifts in the simplest example of Brownian motion following theorem 5.1, we introduce the regularized forward time derivative as limit of conditional expectation of difference, as in (35),

$$\begin{equation*} DX(t) = \lim_{\Delta t\to0^+} \mathbf{E}\left[ \frac{X(t+\Delta t)-X(t)}{\Delta t} \bigg| X(t) \right]. \end{equation*}$$

Let us compute this drift for the time reversal $\overleftarrow X(t) = X(-t)$, we find $D\overleftarrow X(t) = -D_*X(-t)$, where $D_*$, the backward derivative is now defined by,

$$\begin{equation*} D_*X(t) = \lim_{\Delta t\to0^+} \mathbf{E}\left[ \frac{X(t) - X(t+\Delta t)}{\Delta t} \bigg| X(t) \right], \end{equation*}$$

generally distinct from DX(t), as illustrated by the above example. This means that our counterpart of the classical time reversal operator $\mathcal T_c$ becomes

$$\begin{equation*} \mathcal T(X(t),DX(t)) = (X(t), - D_*X(t)). \end{equation*}$$

So we can indeed make all computations using traditional tools of (forward) Itô's calculus, but be aware that their results have a backward version and that fully dynamical laws have to be invariant under $\mathcal T$. For instance, our stochastic Hamiltonian equations (37) have a backward counterpart, needed to re-establish the reversibility of our stochastic dynamical system. So, Schrödinger's reversibility, encapsulated in his nonequilibrium detailed balance condition (8), is a consequence of the stochastic deformation of classical microscopic reversibility, founded on the $\mathcal T_c$ operator.

For instance, coming back to our elementary Brownian case and using the stochastic time reversal $\mathcal T$, the backward version of the free Lagrangian L₀ is $L_0(t, X(t), D_*X(t)) = \frac{1}{2} |D_*X(t)|^2$. The backward version of the free Euler–Lagrange equation becomes

$$\begin{equation*} \begin{aligned} D_* D_* X(t) &= D_* B_* = D_* \left( \frac{X(t) - x}{t-s} \right) \\ &= \left( \frac{\partial}{\partial t} + \frac{X(t) - x}{t-s} \nabla_y - \frac{1}{2}\Delta_y \right) \left( \frac{y - x}{t-s} \right)\bigg|_{y=X(t)} = 0. \end{aligned} \end{equation*}$$

Notice that the minus sign in front of the Laplacian term comes from backward Itô's calculus (see [16] and reference therein).

To conclude this part, we observe that, motivated by the product form of the density in equation (4), Schrödinger was interested in Markovian processes. A more general construction, due to B. Jamison, started from any joint probability $\mu_{t_0,t_1}(\mathrm dx,\mathrm dy)$, not necessarily of the form (11) (equation (11) is the only form producing Markovian processes, cf section 2.2). In general, we call Bernstein's reciprocal processes the processes following from an arbitrary choice of joint probability. Those processes have been proved to be relevant as counterparts of quantum statistical systems [43].

More information and examples can be found in [16, 44].

As another illustration of the key relevance of this aspect in the stochastic dynamical theory inspired by Schrödinger, let us mention that the rigorous (Euclidean) versions of Feynman's path integral interpretation of Heisenberg commutation relation requires this time reversal of the drift (see [16]).

Schrödinger's reinterpretation (by anticipation!) of Feynman's path integral commutation relations highlights the deep analogies between his stochastic dynamical approach and quantum mechanics of pure states. They are rooted in the fundamental role of the two adjoint heat equations (2) and (3), and the Euclidean Born's interpretation of (4). Let us summarize some key relations of the resulting Schrödinger's Euclidean Quantum Mechanics, namely between the Markovian processes solving his stochastic dynamics and operators calculus founded on the class of Hamiltonian operators of theorem 5.1.

Let $f(X(t),t)$ be a smooth function of the Markovian process of theorem 5.1. Then its regularized forward derivative (or infinitesimal generator), as in (36), can be expressed, as a differential operator, by

$$\begin{equation*} Df(X(t),t) = \frac{1}{w(X(t),t)} \left( \frac{\partial}{\partial t} - \mathcal H \right) (\,f\,w)(X(t),t). \end{equation*}$$

The proof is an exercise of Itô's calculus using the form of $\mathcal H$ in theorem 5.1 and the fact that w solves equation (2). Explicitly, the l.h.s. describes the dynamics of X(t) and the r.h.s. has close relation with the Hamiltonian operator $\mathcal H$, generator of Euclidean dynamics, for a given Euclidean 'state' w.

The backward version of $Df(X(t),t)$ results from the time reversed expression, for $w^*$ solving equation (2),

$$\begin{equation*} D_*f(X(t),t) = \frac{1}{w^*(X(t),t)} \left( \frac{\partial}{\partial t} + \mathcal H^* \right) (\,fw^*)(X(t),t), \end{equation*}$$

where $\mathcal H^* = \mathcal H$ as $\mathcal H$ is self-adjoint.

It is possible to construct a Hilbert space approach of Schrödinger's Euclidean dynamics, describing how to translate properties of the process X(t) in terms of an operators calculus, in analogy with the quantum model [40]. In it, a number of purely informal results of Feynman's path integral method make perfect mathematical sense (Heisenberg's commutation relations being only one of them).

The expectation of (Euclidean) observables in Hilbert space take a familiar form. For instance, the one of the momentum is

$$\begin{equation*} \langle P \rangle_w = \int w^* P^+ w \mathrm dx = \int ww^* \left( \frac{P^+ w}{w} \right) \mathrm dx, \end{equation*}$$

for $P^+ = \nabla$ so that the forward drift is, indeed $B = \frac{\nabla w}{w}$. Since $P^+$ is not symmetric, the backward one is $\int ww^* \left( \frac{P w^*}{w^*} \right) \mathrm dx$ namely $B^* = -\frac{\nabla w^*}{w^*}$.

Regarding the energy observable, the Hamiltonian $\mathcal H$ of theorem 5.1 is symmetric and

$$\begin{equation*} \langle \mathcal H \rangle_w = \int w^* \mathcal H w \mathrm dx = \int ww^* \left( \frac{\mathcal H w}{w} \right) \mathrm dx, \end{equation*}$$

namely, after using $\frac{\Delta w}{w} = |\frac{\nabla w}{w}|^2 + \nabla \cdot \left( \frac{\nabla w}{w}\right)$, we find the random variable $\frac{1}{2}|B|^2 + \frac{1}{2}\nabla\cdot B - U$, where the sign change of the potential is usual in the Euclidean world and the extra divergence term comes from Itô's deformation or, equivalently, the SO term of HJB equation (41). Indeed, in this Euclidean case, $B = \nabla S$, $H_0 = \frac{1}{2}|p|^2 - U(x)$ and the SO correction $\frac{1}{2} g^{\,jk} (o_{jk} - p_i \Gamma_{jk}^i)$ reduces to $\mathrm{tr}\ o$.