Gradient Flow Structure of a Multidimensional Nonlinear Sixth Order Quantum-Diffusion Equation

The following nonlinear parabolic evolution equation of sixth order is considered:

where $\lambda\geq 0$ is a given parameter. There are at least two different contexts in which (1) plays a role.

The first is the semi-classical approximation of the nonlocal quantum drift-diffusion model by Degond et al [5]. In the formal asymptotic expansion of that equation in terms of the Planck constant $\hbar$, the right-hand side of (1) with $\lambda=0$ appears as the term of order $\hbar^{4}$, after the linear diffusion operator $\Delta u$ at order $\hbar^{0}$ and the operator $-\nabla^{2}:(u\nabla^{2}\log u)$, which is related to the Bohm potential, at order $\hbar^{2}$. More details on the derivation of the model and its formal expansion are given below in Section 2.1.

The other context, more relevant to the paper at hand, is that of gradient flows in the $L^{2}$-Wasserstein distance. Consider the following three functionals, defined on — at the moment for simplicity: strictly positive — probability densities $u:{\mathbb{R}}^{d}\to{\mathbb{R}}$ by

which we shall refer to as — $\lambda$-perturbed if $\lambda>0$ — entropy, Fisher information, and energy, respectively. The celebrated result of [9] is that the gradient flow of $\mathcal{H}_{\lambda}$ in the $L^{2}$-Wasserstein metric is the linear Fokker-Planck equation,

In [8] — see also[1, Example 11.1.10] — it has been shown that the gradient flow of $\mathcal{F}_{\lambda}$ is the so-called quantum drift-diffusion or DLSS equation,

The starting point for our analysis is that (1) is the gradient flow of $\mathcal{E}_{\lambda}$, at least formally, that is, the potential $\Phi$ in (1) is the variational derivative of $\mathcal{E}_{\lambda}$. The reason for considering $\mathcal{E}_{\lambda}$ as potential for a gradient flow is more profound than its formal similarity with $\mathcal{H}_{\lambda}$ and $\mathcal{F}_{\lambda}$. There is an intimate relation between $\mathcal{H}_{\lambda}$, $\mathcal{F}_{\lambda}$ and $\mathcal{E}_{\lambda}$, that has already been the basis for deriving sharp self-similar asymptotics in [13], and that we shall elaborate on in Section 1.3 below. In a nutshell, the dissipation of $\mathcal{H}_{\lambda}$ along the heat flow is $\mathcal{F}_{\lambda}$, the dissipation of $\mathcal{F}_{\lambda}$ along the heat flow is $\mathcal{E}_{\lambda}$, up to a multiple of $\mathcal{F}_{\lambda}$ itself, and this equals the dissipation of $\mathcal{H}_{\lambda}$ along the flow of (4). In this spirit, one may consider (4) and (1), respectively, as fourth and sixth order analogues of the second order linear Fokker-Planck equation (3).

Our analytical results are two-fold. First, we give a proof of existence of weak solutions to the initial value problem for (1) on the whole space ${\mathbb{R}}^{d}$ for initial data with finite entropy and finite second moment. This result is proven with full rigor. Second, we study the long-time asymptotics of solutions using a linearization around the steady state. This part is formal in the sense that we calculate the linearization for sufficiently smooth perturbations of the steady state and discuss the spectral properties of an appropriate closure of the linear operator.

Global existence of non-negative weak solutions to (1) with $\lambda=0$ on the $d$-dimensional torus, i.e., with periodic boundary conditions, has been proven in [11] for $d=1$, and in [3] for $d=2$ and $d=3$. The main technical ingredient of these proofs is a particular regularization of the evolution equation, namely by $\varepsilon\Delta^{3}\log u$. This regularization produces approximations of the true solution that are smooth and have an $\varepsilon$-dependent positive lower bound. Smoothness of both $u$ and $\log u$ then allows to perform all the necessary a priori estimates on the approximation, and these pass to the limit $\varepsilon\downarrow 0$. An extension of this method from the torus to the whole space — if possible at all — is at least not straight-forward, since the intermediate estimates involve simultaneous Sobolev estimates on $u$ and $\log u$, that would contradict each other on unbounded domains.

Here, we shall prove existence of solutions in $d\leq 3$ on the whole space ${\mathbb{R}}^{d}$. Our technical device is very different from the one recalled above: we invoke the machinery of metric gradient flows. Our approximations are of lower regularity, and we do not have any information about positivity. To justify the a priori estimates, we need the method of flow interchange with the heat equation [8, 13]. In contrast to the constructions in [11, 3], we do not modify the equation, but perform a variational discretization in time using the minimizing movement scheme. More specifically: given initial probability density $u_{0}$ of finite second moment and entropy, and a time step $\tau>0$, define inductively a sequence $(u_{\tau}^{n})_{n=0}^{\infty}$ by $u_{\tau}^{0}=u_{0}$, and $u_{\tau}^{n}$ being a minimizer of

We recall the $L^{2}$-Wasserstein-distance ${\mathbf{W}}_{2}$ below in Section 2.2, and we prove that the inductive procedure is well-defined. Our result about existence is that the sequences $(u_{\tau}^{n})$ approximate a weak solution to (1).

Theorem 1.1 is proven by time-discrete approximation of the solution $u$ via the celebrated minimizing movement scheme. The main compactness estimate for passing to the time-continuous limit is provided by the dissipation of the (unperturbed) entropy $\mathcal{H}_{0}$, which formally amounts to

with some positive $\kappa>0$. The formal calculations leading to (10) via integration by parts are identical to the ones used in [3]. The justification of these estimates in the whole-space case is technically more involved.

The interpretation of (1) as a Wasserstein gradient flow is essential for our second main result, which concerns the long time asymptotics of $u$. First assume $\lambda>0$, in which case there is an unique equilibrium $U_{\lambda}$ for (1) (with mass equal to 1), given by

Our approach to understanding the dynamics near $U_{\lambda}$ is to formally calculate a linearization of (1) at $U_{\lambda}$ and to determine its spectrum. We wrote a linearization since there are several competing concepts for linearization in nonlinear diffusion equations, providing different pieces of information about the long-time asymptotics. Following the ideas of [7], we study here the “linearization in Wasserstein”, which is given by the so-called displacement Hessian of $\mathcal{E}_{\lambda}$ at $U_{\lambda}$. Very informally, the displacement Hessian of a functional ${\mathcal{U}}$ at a critical point $U_{*}$ is the representation of ${\mathcal{U}}$’s second variational derivative with respect to the scalar product in $H^{-1}({\mathbb{R}}^{d};U_{*}\,\mathrm{d}\mathcal{L}^{d})$. A definition and a more intuitive interpretation in terms of Lagrangian maps is given in Section 4.3.

Displacement Hessians are rarely used for the analysis of long-time asymptotics since the extraction of rigorous analytical information requires a lot of a priori knowledge about regularity from the solution. Particularly when it comes to proving higher order asymptotics, alternative linearizations are often easier to handle; we refer to the discussions in [21, 6]. The most famous application of displacement Hessians concerns the result on self-similar asymptotics for the porous medium equation [17]; further applications, also to higher order diffusion equations, can be found e.g. in [14, 15, 20]. In the situation at hand, we use the Wasserstein linearization because of its compatibility with the special structure of (1) that we outline below.

Since the derivation of the result is probably more interesting than the result itself, we briefly indicate the main ingredient, which is the aforementioned intimate relation between the three functionals in (2): recall that the linear Fokker-Planck equation (3) is the Wasserstein gradient flow of $\mathcal{H}_{\lambda}$. Consider a solution $(w_{r})_{r\geq 0}$ to that flow, i.e.,

It is a well-known fact that $\mathcal{F}_{\lambda}$ is the dissipation of $\mathcal{H}_{\lambda}$ along (12), i.e.,

The connection between $\mathcal{H}_{\lambda}$ and $\mathcal{E}_{\lambda}$ — and, implicitly, also $\mathcal{F}_{\lambda}$ — is given by

The relation (14) is derived in Section 4.2. It will play the same role in our studies on (1) as (13) has played for the analysis of the DLSS equation (4) in [14]. That is, we use (14) to express the displacement Hessian of $\mathcal{E}_{\lambda}$ in terms of the displacement Hessian of $\mathcal{H}_{\lambda}$.

For the linear Fokker-Planck equation (12), which is the gradient flow of the relative entropy $\mathcal{H}_{\lambda}$ and has $U_{\lambda}$ as a critical point as well, the displacement Hessian has been calculated in [7, Proposition 2]: it is the extension of the formal operator

to $W^{1,2}({\mathbb{R}}^{d};U_{\lambda}\,\mathrm{d}\mathcal{L}^{d})$. The corresponding spectrum of ${\mathbf{L}}_{\lambda}$ is well-known: it is purely discrete with eigenvalues that are precisely the positive integer multiples of $\lambda$. The corresponding eigenfunctions are Hermite polynomials. The derivation of ${\mathbf{L}}_{\lambda}$ might appear a bit weird: first, one re-writes the linear Fokker-Planck equation (12) on the $L^{2}$-Wasserstein space, obtaining a highly non-linear gradient flow, and then calculates its linearization, which is — up to dualization — again the original equation (12). The key point is that the Wasserstein linearization is compatible with the relation (14), i.e., the linearization of the highly non-linear sixth order equation (1) is easily expressible in terms of ${\mathbf{L}}_{\lambda}$, see Theorem 1.3 below. A similar idea has been used in [14], where the authors showed by exploiting the relation (13) below that the displacement Hessian for $\mathcal{F}_{\lambda}$ is formally given by ${\mathbf{L}}_{\lambda}^{2}$, and discussed implications on the long-time asymptotics of the nonlinear fourth order DLSS (4) equations.

Our second main result is:

Theorem 1.3 could be proven by direct calculations, evaluating the second variational derivative of $\mathcal{E}_{\lambda}$ along the transport equation, and then integrating by parts until the desired form is attained. This would be tedious but in principle straight-forward once the desired terminal form ${\mathbf{L}}_{\lambda}^{3}+\lambda{\mathbf{L}}_{\lambda}^{2}$ is known. Here we present a more conceptual approach, based on the relation (14), which leads us to the form of the Hessian in the first place.

As said before, the implications of Theorem 1.3 on the long-time asymptotics of (1) are far from obvious. Naturally, the hope is that the dynamics of all sufficiently smooth solutions $u$ close to equilibrium is approximately the same as that of the linearized equation, and in particular, that the eigenvalues of ${\mathbf{L}}_{\lambda}^{3}+\lambda{\mathbf{L}}_{\lambda}^{2}$ determine the exponential decay of $u$’s “nonlinear modes” in the long-time limit. It is not hard to see that the formal differential operator ${\mathbf{L}}_{\lambda}^{3}+\lambda{\mathbf{L}}_{\lambda}^{2}$ possesses a unique self-adjoint closure in $W^{1,2}({\mathbb{R}}^{d};U_{\lambda}\,\mathrm{d}\mathcal{L}^{d})$, and that the spectrum of the latter is purely discrete with eigenvalues $\lambda^{3}(k^{3}+k^{2})$ for $k=0,1,2,\ldots$

The strong results from [13], where the fully nonlinear exponential stability of the Gaussian steady state for the DLSS equation (4) has been proven on grounds of (13), give rise to a conjecture, namely that the spectral gap $2\lambda^{3}$ in the linearization actually determines the global rate of convergence to equilibrium.

A direct consequence of Conjecture 1 would be $u$’s asymptotic self-similarity in case $\lambda=0$. More precisely, in Section 4.6 we show that if Conjecture 1 were true, then any solution $u$ to (1) with $\lambda=0$ approaches in $L^{1}({\mathbb{R}}^{d})$ with algebraic rate $t^{-1/6}$ the self-similar solution

at least if $u$ is already sufficiently close to self-similarity initially.

We sketch the derivation of (1) from a quantum mechanical model. In [5] — building on [4] — the following non-linear and non-local quantum analogue of the classical Fokker-Planck equation (12) has been derived:

Here $u$ is the macroscopic density of quantum particles whose dynamics aims at minimizing the ensemble’s relative von Neumann entropy $\mathtt{H}_{\lambda}$. The precise definition of $A[u]$, sometimes referred to as $u$’s quantum logarithm, is intricate; for the sake of completeness, we mention one possible definition of $A[u]$ as the uniquely determined potential $A:{\mathbb{R}}^{d}\to{\mathbb{R}}$ such that

for all test functions $\varphi\in C^{\infty}_{c}({\mathbb{R}}^{d})$. Here $\operatorname{Tr}$ denotes the trace over $L^{2}({\mathbb{R}}^{d})$, and $\exp$ is the exponential of a self-adjoint operator.

Under the specific hypotheses made in [5], the von Neumann entropy can be expressed in terms of the macroscipic density $u$ alone, and is given by

We remark that $\mathtt{H}_{\lambda}(u)$’s variational derivative is $A[u]+\frac{\lambda}{2}|x|^{2}$, and thus equation (18) has the formal structure of a gradient flow in ${\mathbf{W}}_{2}$. In the semi-classical limit $\hbar\to 0$, the von Neumann entropy $\mathtt{H}_{\lambda}(u)$ formally approaches the Boltzmann entropy $\mathcal{H}_{\lambda}(u)$, and the variational derivative $A[u]$ formally approaches the pointwise logarithm $\log u$. Consequently, (18) turns into the classical Fokker-Planck equation.

The existence analysis for the full equation (18) goes far beyond classical parabolic theory, and has been carried out just recently, and only in one space dimension [18]. Already the rigorous definition of $A[u]$ as solution to an inverse problem has been challenging [16]. A way to approximate (18) by local equations is via the expansion of $A[u]$ in powers of the small parameter $\hbar$. In [2, Appendix] the following asymptotic expansion up to $\mathcal{O}(\hbar^{6})$ has been computed (for $\lambda=0$):

As mentioned above, a reduction of $A[u]$ in (18) to the leading order term $A_{0}[u]$ yields the linear Fokker-Planck (or rather: heat — since $\lambda=0$) equation (12),

Replacing $A[u]$ by $A_{1}[u]$ leads to the Derida-Lebowitz-Speer-Spohn (DLSS) equation

Finally, since $A_{2}[u]$ is identical to the functional $\Phi(u)$, the equation $\partial_{t}u=\mathrm{div}(u\nabla A_{2}[u])$ coincides with the sixth-order equation (1) under consideration here.

In this section we briefly review some basics about the theory of optimal transport and $L^{2}$-Wasserstein gradient flows, but only as far as it is needed later. For a more profound introduction to these topics, we refer to the text books [1, 19, 22]. By $\mathcal{P}_{2}(\mathbb{R}^{d})$ we denote all probability measures with finite second moment,

We shall frequently identify absolutely continuous (with respect to the Lebesgue measure $\mathcal{L}^{d}$) probability measures $\mu$ on ${\mathbb{R}}^{d}$ with their densities $u\in L^{1}({\mathbb{R}}^{d})$. A sequence $(\mu_{n})_{n\in\mathbb{N}}\subset\mathcal{P}(\mathbb{R}^{d})$ of probability measures converges narrowly to $\rho\in\mathcal{P}(\mathbb{R}^{d})$ if

holds for all bounded, continuous functions $f:\mathbb{R}^{d}\rightarrow\mathbb{R}$. The $L^{2}$-Wasserstein distance between two measures $\mu,\nu\in\mathcal{P}_{2}(\mathbb{R}^{d})$ is defined via

where $\Pi(\mu,\nu)$ denotes the set of all transport plans between $\mu$ and $\nu$, that is all $\pi\in\mathcal{P}(\mathbb{R}^{d}\times\mathbb{R}^{d})$ with respective marginals $\mu$ and $\nu$. The Wasserstein distance metrizes narrow convergence on $\mathcal{P}_{2}(\mathbb{R}^{d})$ and is itself lower semi-continuous with respect to that convergence.

We are not going to define metric gradient flows on $(\mathcal{P}_{2}(\mathbb{R}^{d}),\textnormal{W}_{2}(\cdot,\cdot))$ in general. Here we just need a particularly nice subclass.

Solutions to (1) will be constructed via discrete-in-time approximation by means of the minimizing movement scheme, i.e., a variational Euler method, see Section 3.2 below. Inductively, the approximation $\mu_{\tau}^{n}$ at time $t=n\tau$ is obtained from $\mu_{\tau}^{n-1}$, the one at $t=(n-1)\tau$, as minimizer in

For passage to the continuous limit, a priori estimates independent of the time step $\tau$ are essential. The key is to give a rigorous meaning to a priori estimates related to dissipations $-\frac{\mathrm{d}}{\,\mathrm{d}t}{\mathcal{V}}(\mu_{t})$ of auxiliary functionals ${\mathcal{V}}$ already on the time-discrete level. For this, we shall use the so-called flow interchange method [13, Theorem 3.2], where the minimizer $\mu_{\tau}^{n}$ is perturbed along the $\alpha$-flow $\mathcal{S}^{\mathcal{V}}_{(\cdot)}$ of the auxiliary functional ${\mathcal{V}}$.

We begin by giving a proper definition of the energy functional.

The energy $\overline{\mathcal{E}}_{\lambda}$ defined above is part of a family of second-order functionals that have been studied in [8, Section 3] in connection with the DLSS equation (4). We recall and adapt two results here.

It is important to remark here that [8, Corollary 3.2] indeed applies to the case $\Omega={\mathbb{R}}^{d}$: the derivation of (29) only involves an integration by parts with the vector field $\mathbf{v}=|\nabla\sqrt[4]{u}|^{2}\nabla\sqrt{u}$, which is integrable on ${\mathbb{R}}^{d}$ since $\sqrt{u}\in W^{2,2}({\mathbb{R}}^{d})$ and $\sqrt[4]{u}\in W^{1,4}({\mathbb{R}}^{d})$. We refer to Lemma C.1 and to the subsequent discussion in the Appendix.

Let a time step size $\tau>0$ be fixed. For a given $\nu\in\mathcal{P}_{2}(\mathbb{R}^{d})$, define the Yosida-penalized energy functional $\overline{\mathcal{E}}_{\lambda,\tau}(\cdot;\nu)$ by

As a consequence of Lemma 3.5, the minimizing movement scheme for $\overline{\mathcal{E}}$ is well-defined for every initial condition $\mu_{0}\in\mathcal{P}_{2}(\mathbb{R}^{d})$. That is, starting from $\mu_{\tau}^{0}:=\mu_{0}$, one can define a sequence $(\mu_{\tau}^{n})_{n=0}^{\infty}$ inductively by choosing for $\mu_{\tau}^{n}$ as a minimizer of $\overline{\mathcal{E}}_{\tau}(\cdot;\mu_{\tau}^{n-1})$ for $n=1,2,\ldots$ For notational convenience, we also introduce the usual time-discrete “interpolation” $\bar{\mu}_{\tau}:{\mathbb{R}}_{\geq 0}\to\mathcal{P}_{2}(\mathbb{R}^{d})$ by

for $n=1,2,\ldots$, and $\bar{\mu}_{\tau}(0)=\mu_{0}$. The sequence $(\mu_{\tau}^{n})_{n=0}^{\infty}$ and its interpolation $\bar{\mu}_{\tau}$ satisfy a variety of energy estimates, that directly follow from the construction via minimization.