Contractions in the 2-wasserstein lenght space and thermalization of granular media

CONTRACTIONS IN THE 2-WASSERSTEIN LENGTH SPACE AND THERMALIZATION OF GRANULAR MEDIA JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI Abstract. An algebraic decay rate is derived which bounds the time required for velocities to equilibrate in a spatially homogeneous flowthrough model representing the continuum limit of a gas of particles interacting through slightly inelastic collisions. This rate is obtained by reformulating the dynamical problem as the gradient flow of a convex energy on an infinite-dimensional manifold. An abstract theory is developed for gradient flows in length spaces, which shows how degenerate convexity (or even non-convexity) — if uniformly controlled — will quantify contractivity (limit expansivity) of the flow. 1. Introduction It has been known since the work of Otto [43] that various familiar diffusion equations can be considered, at least heuristically, to be gradient flows on the space of probability measures, endowed with a manifold structure and local metric whose arc length distance coincides with the quadratic Wasserstein distance, dist2(ρ0, ρ1) = inf ½Z|v−w|2dγ(v, w); γ∈Γ(ρ0, ρ1)¾1/2 ; (1.1) here Γ(ρ0, ρ1) is the set of probability measures on Rd×Rdhaving marginals ρ0and ρ1. Otto showed how to use these heuristics to study the long-time behavior of nonlinear porous-medium type equations. His work has inspired numerous developments, some of which are reviewed in [52]. The present paper deals with applications of this point of view to diffusion equations whose nonlinearities may also present a nonlocal structure, as found in the kinetic models of Benedetto et al for equilibration of velocities in granular media [11, 12, 13]. It is the sequel to our previous work [23], in which we studied these equations by means of entropy methods, or more precisely, convergence of the energy functional towards its infimum as time becomes large. In the present paper we shall pursue two goals. The first of these is to complement our previous study by estimating 1 2 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI rates of convergence in quadratic Wasserstein distance. Although convergence in Wasserstein distance may be weaker than convergence in entropy sense, as explained and illustrated in [23], this approach offers several advantages: 1) Wasserstein distance is the natural distance associated with the gradient flow structure under examination; 2) the assumption of finite Wasserstein distance is much more general than the assumption of finite entropy; 3) most importantly, this approach enables one to directly compare two different solutions, instead of just comparing each solution to the stationary one. Thus much information is gained about the short-time behaviour of the flow as well as its long-time asymptotics. For instance, when we can show that the distance between any two solutions does not grow too quickly as a function of time, uniqueness of solutions and extension of the flow to singular initial data follow immediately. If these distances actually decrease — which is often the case — then existence and uniqueness of a fixed point may also be inferred from contractivity. This investigation will lead us to examine in fine detail the structure of the space of probability measures equipped with the Wasserstein distance. Thus the second goal of our paper is to develop a formal mathematical framework for Otto’s ideas. To do this, we shall study the space of probability measures P2(Rd) with finite second order moments, viewed as a length space. This provides a conceptual setting in which many known results, and some new ones, fall into place naturally. We introduce an additional structure, which we call a Riemannian length space, to axiomatize key ingredients of Otto’s approach and serve as a convenient framework for converting his heuristical arguments more directly into rigorous theorems. This framework is successful in allowing us deal with smooth and positive densities evolving either on bounded domains Ω ⊂⊂ Rdwith no-flux boundary conditions or on Rdbut rapidly decaying as v→ ∞; it still requires the approximations employed by Otto to extend the results to general initial data or evolutions on the whole of Rdwhose phase space decay is not exponential. This paper has been in gestation for quite a long time: its results were announced already at the Azores TMR Summer School on Mass Transportation Methods in Kinetic Theory and Hydrodynamics (4-9 September 2000). The intervening years have seen a number of related and important independent developments illustrating the relevance and power of these ideas. A study of this length space structure was performed by Ambrosio, Gigli & Savar`e [7, 8]. Their construction has a lot to do with ours, even if some of their goals and the tools that they employ are quite different: the authors in [8] studied general results for existence of gradient flows for convex energy functionals on this length space, establishing the 2-uniform contraction rates discussed below as a byproduct of their investigation. For this CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 3 reason, they develop a more general theory, which enables them to handle singular measures, whereas our analysis is largely restricted to probability measures which are absolutely continuous with respect to Lebesgue measure. (Our Riemannian length space structure can be extended to handle singular measures, but the resultant ambiguities in particle labelling then lead to geometrical pathologies). Their focus is on absolutely continuous curves in the ambient length space, whereas ours is on global parameterizations of this space using exponential maps in lieu of Riemann normal coordinates at each point. We refer to [8] for further explanations. Our main theorem gives an explicit estimate on the growth or decay of the distance between any two solutions of a subgradient flow in the Riemannian length setting that we introduce. To apply it in the particular case of the 2-Wasserstein for probability measures and to the partial differential equation models we deal with, we either work with smooth and positive solutions evolving on a bounded domain, or else need to perform a series of approximations to the equations which does not close all the possible cases (see last section for precise open problems). Another approach would be to try to extend the existence theory and contractivity estimates of Ambrosio, Gigli & Savare [8] to energies with degenerate convexity, and then verify the resulting flows coincide with the partial differential evolutions we are interested in. Other authors exploring similar themes include Carlen & Gangbo, who in their investigation of the kinetic Fokker-Planck equation show the length space (P2(Rd),dist2) posseses a conical structure [19] analogous to a warped product in Riemannian geometry, and Cordero-Erausquin, Gangbo & Houdre [27], who establish various expressions of uniform displacement convexity for entropies E:P2(Rd)−→ Rwith respect to more general costs on P2(Rd). When the cost is given by dist2, these relate to our rate of convergence results. The displacement convexity of such entropies — which amounts to convexity along geodesics in (P2(Rd),dist2) — originated in work of McCann, where it was established using a particular geodesic structure without reference to an underlying metric [40]. The application of displacement convexity to rates of convergence in nonlinear evolution equations begun by Otto [43], was recently explored for more general costs associated with different nonlinearities by Agueh [1, 4] and Agueh, Ghoussoub & Kang [2] [3]. Finally Sturm & von-Renesse [49] have shown equivalence of 2-uniform semiconvexity of the Boltzmann entropy (or exponential contractivity of the heat semigroup with respect to 2-Wasserstein distance) to lower bounds for the Ricci curvature on a Riemannian manifold. Sturm [46] [47] [48] and Lott & Villani [38] extended this idea to nonlinear evolution equations, and to a means of defining Ricci curvature bounds in abstract metric-measure spaces. 4 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI Typical equations to which the present considerations apply take the form ∂ρ ∂t =∇·[ρ∇(A0(ρ) + B+C∗ρ)] ,(1.2) where ρ: [0, T]×Ω−→ [0,+∞] is an integrable density, Ω ⊂Rd, and A: [0,∞)−→ Rand B, C :Rd−→ Rare convex potentials. In the one-dimensional models for granular media which motivated our original interest, ρ(t, v) represents a distribution of velocities v∈Ω at each time, and the three potentials model the effects of: A) random interactions of the granules with their environment (a fluid or heat bath), B) friction, and C) inelastic collisions between granules with different velocities — the nonlocal source of nonlinearity. Notice equation (1.2) is appropriate to spatially homogeneous initial conditions — meaning ρ(t, x, v) = ρ(t, v) depends only on the velocity coordinate vin phase space and not the position x— so it would be natural to study the evolution on the entire space tangent space of velocities v∈Ω = Rd. However, for technical reasons, as in [43, 22], it is often convenient to begin by formulating the problem on a bounded convex domain of velocities Ω ⊂⊂ Rdwith no-flux boundary conditions, ρνΩ·∇(A0(ρ) + B+C∗ρ) = 0,on (t, v)∈[0, T]×∂Ω,(1.3) and study the large domain limit Ω →Rdsubsequently. Here νΩdenotes the outer unit normal at vto Ω. Later on in the text, we may use the variables xor yin place of v; regardless of its name, our independent variable always represents a velocity in kinetic models. The notion of solution for which we will obtain our rates of decay results will be detailed in Section 6, but let us announce that rates of decay will be obtained for smooth solutions and thus, generalized for weak solutions constructed by suitable approximation procedures. The spirit of our results is captured by the following examples. We assume power law potentials here for simplicity; more general potentials are addressed in later sections. For a6= 1, we will discuss primarily smooth, positive densities evolving on a bounded domain; we go through the details of extending the conclusions to the entire space Ω ⊂Rdonly for the most relevant situation of linear diffusion (a= 1). Example 1.1. Take B(v) = β|v|b+2/(b+ 2) and C(v) = γ|v|c+2/(c+ 2), and A(%) = ½α%a/(a−1) 1 6=a≥max{d−1 d,d d+2 } α% log % a = 1 (1.4) CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 5 with α, β, γ ≥0and b≥0, c ≥ −1. The Wasserstein L2distance dt:= dist2(ρ1(t), ρ2(t)) between any two solutions of (1.2–1.3) on Ω⊂⊂ Rddecays like dt≤½e−βtd0b= 0 d0(1 + βtb(d0/2)b)−1/b ∼2(βtb)−1/b b > 0(1.5) in the presence of friction β > 0. When friction is negligible, meaning β= 0, the inelastic collisions γ > 0alone yield a decay rate dt≤½e−γtd0c= 0 d0(1 + γtc(d0/√2)c)−1/c ∼√2(γtc)−1/c c > 0,(1.6) provided the center of masses of the two solutions coincide at each point in time; this will be true if, for example, we assume reflection symmetry Ω = −Ωand ρ(0, v) = ρ(0,−v)initially (and hence for all time). In the most interesting cases of interaction potentials C(v) = γ|v|c+2/(c+ 2), we are able to overcome the restriction of reflection symmetry by approximating the solution using very smooth fixed center of mass solutions of the same equation which decay quickly at infinity on all of Rd [23]. Precise statements are given in the last section. Choosing d= 1, a= 1, b= 0, c= 1 produces the one-dimensional granular models of Benedetto, Caglioti, Carrillo, & Pulvirenti [11, 12]. There the presence or absence of friction can mean the difference between exponentially fast and algebraically slow thermalization: indeed, Benedetto, Caglioti & Pulvirenti’s original calculation shows that neither the constants nor the exponent of the algebraic bound (1.6) can be improved when α=β= 0. In this special case, all velocities converge to a single equilibrium value, and the slow convergence results from the rate of collisions dwindling to zero along with the dissipated energy per collision. The mathematical reason for this algebraic rate is collapse of the relative sand grain velocities v−v onto the unique point where the second derivative of the collision potential C(v−v) = |v−v|3/3 vanishes. We eventually showed in our companion paper how exceptional this example is: the algebraic bound (1.6) can be improved to an exponential bound provided α > 0; the presence of a heat bath speeds up thermalization by ensuring that neither the rate of collisions nor the dissipated energy becomes too small. The resulting bound differs from (1.5) however, in that the exponential rate of contraction we derive in this case is not global, but depends on the initial entropy of ρ1(0) and ρ2(0) [23]. For a6= 1, it may also be possible to extract the same results in the large domain limit Ω →Rd, but a complete discussion of most general conditions which permit this would form the subject of separate treatise; the presence 6 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI of friction β > 0 above is sufficient if b≥0. When β= 0 the center of mass condition is necessary for convergence in Wasserstein distance on Ω = Rd: translation invariance implies the average velocities hviρ1(t)and hviρ2(t)from (5.5) do not change; if they differ initially then dt≥ |hviρ1(0) − hviρ2(0)| cannot converge to zero. Compare how barycenter enters explicitly in the inequalities formulated by Agueh, Ghoussoub & Kang [2]. Let us also mention that Wasserstein contraction estimates have been obtained recently by Li & Toscani [37] for the family of one dimensional granular media models introduced in [51]. Their main idea is to use the particular explicit formula of the Wasserstein distance in one dimension. In fact, the optimal transport map in one dimension is always the same for all convex costs and is defined in terms of the inverse distribution functions of the measures involved. A short review of these ideas applied to one dimensional nonlinear diffusion-dominated equations can be found in [25]. Wasserstein contraction estimates play a role in controlling the expansion of the support of solutions for one dimensional nonlinear diffusions as recently pointed out in [21] for the porous medium equation and in [20] for diffusiondominated equations. Similarly, Bolley, Brenier & Loeper have found that scalar conservation laws in one-dimension contract the p-Wasserstein distance between the derivatives of the solutions for all p∈[1,∞] [17]. Finally, let us also point out that a related equation that it is also included in this theory (at least formally) is the one dimensional nonlinear Fokker-Planck equation arising in free probability [15], also called the free Fokker-Planck equation. The linear diffusion term is replaced by the Hilbert transform in this equation. The free Fokker-Planck equation has also a formal gradient flow structure with respect to a logarithmic interaction energy functional. In one dimension, this energy happens to be displacement convex in the sense of McCann [40], as observed and exploited by Blower in the context of random matrix theory [16]. 2. A schematic cartoon of the rate arguments Before attempting to construct an abstract argument in a context fraught with perils of nonsmoothness, infinite dimensions, and degenerate convexity, it is instructive to recall the ideas behind the convergence arguments in their simplest form. The setting will be so simple that not only are the results well-known, they could all be deduced by a good sophomore calculus student. Nevertheless, they serve to contrast the contraction strategy developed hereafter with the Bakry-Emery [10] type entropy production analysis employed in Otto [43] and in our previous work [23]. CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 7 Fix E∈C2(Rd) and consider solutions of the ordinary differential equation dxt dt =−∇E(xt) (2.1) corresponding to steepest descent or gradient flow on the energy (entropy) landscape determined by E. Here Iwill denote the d×didentity matrix. Proposition 2.1 (Contraction / expansion bounds in a semi-convex valley). Fix k∈R. If E∈C2(Rd)satisfies D2E(x)≥kI throughout Rd, and the curves xtand t∈[0,∞)−→ yt∈Rdboth solve the differential equation (2.1), then |xt+t0−yt+t0| ≤ e−kt|xt0−yt0|. Proof. Set f(t) = |xt−yt|2/2. Then f0(t) = −hxt−yt,∇E(xt)−∇E(yt)i =−hxt−yt,Z1 0 D2E[(1 −s)xt+syt] (yt−xt)dsi ≤ −2kf(t)Z1 0 ds. Gronwall’s inequality (integration) implies the desired result: f(t+t0)≤ e−2ktf(t0). ¤ If k > 0, more can be achieved. The convexity of Eis said to be 2- uniform, and we have shown that the solution map x0∈Rd−→ Xt(x0) = xt of the initial value problem (2.1) defines a uniform contraction on Rdfor each t > 0. The C2smoothness of Eensures that the solution map is welldefined locally in space and time; the map is globally defined for all future times since xtis constrained to lie in the level set {x|E(x)≤E(x0)}, whose compactness follows from the coercivity of E(x)≥E(x0) + h∇E(x0), x − x0i+k|x−x0|2/2. Since Rdis complete, the contraction mapping principle dictates that this map has a unique fixed point Xt(x∞) = x∞∈Rd, and each solution curve xt=Xt(x0) must converge to x∞in the long time limit t→ ∞. If we are only interested in the rate of convergence to x∞, an alternative to Proposition 2.1 can be based on the Bakry and Emery entropy production approach. We give that argument here for comparison’s sake. The quantity estimated is the decay rate of the slope |∇E(xt)| → 0; by the analogy discussed at the end of this section, the square of this slope is called the information. Proposition 2.2 (Entropy production and information decay rate).Fix k∈Rsuch that E∈C2(Rd)satisfies D2E(x)≥kI throughout Rd. Then any solution t∈[0,∞)−→ xt∈Rdof (2.1) satisfies |∇E(xt+t0)| ≤ e−kt|∇E(xt0)|. 8 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI Proof. Let f(t) := |∇E(xt)|2/2. Then −f0(t) = −h∇E(xt), D2E(xt) ˙xti =h∇E(xt), D2E(xt)∇E(xt)i ≥2kf(t), and Gronwall’s inequality proves the desired estimate: f(t+t0)≤e−2kt f(t0). ¤ While the conclusions of these two propositions are not immediately comparable, the following consequence (2.2) of 2-uniform convexity relates them. It shows that information dominates the altitude or relative entropy E(x)−E(x∞), which in turn dominates horizontal distance squared. Thus in the limited range of validity k > 0 and yt:= x∞, and apart from constants, Proposition 2.2 trumps Proposition 2.1. On the other hand, (2.3) also shows that if information remains bounded, then convergence in the weakest sense, namely of distance (unsquared), also implies convergence in the stronger sense of relative entropy. Lemma 2.3 (Manifestations of 2-uniform convexity).Let 0≤f∈C2(R) satisfy f(0) = 0 and f00(s)> k > 0for all s∈R. Then ks2≤2f(s)≤ k−1|f0(s)|2and f(s)≤sf0(s)−ks2/2. Proof. Let g(s) := f(s)−ks2/2. Taking two derivatives shows g(s) is convex, so its critical point at the origin must be a minimum: g(s)≥g(0) = 0. This proves the first inequality. Since f(s)≥0 is strictly convex, its minimum f(0) = 0 is its only critical point. Defining h(s) := |f0(s)|2/2−kf(s), we see that f0(s) and hence h0(s) = f0(s)(f00(s)−k) have the same sign as s∈R. Thus h(s) has no critical points apart from a global minimum at zero, and the second inequality is established: h(s)≥h(0) = 0. Finally, let e(s) = sf0(s)−ks2/2−f(s). Then e0(s) = s(f00(s)−k) also has the same sign as s∈R, so its unique critical point is a global minimum at zero: e(s)≥e(0) = 0 to complete the proof of the lemma. ¤ Corollary 2.4 (Log Sobolev, transportation, and HWI inequalities).Suppose E(x∞)≤E(x)∈C2(Rd)and D2E(x)≥kI > 0for all x∈Rd. Then k 2|x−x∞|2≤E(x)−E(x∞)≤1 2k|∇E(x)|2(2.2) and E(x)−E(x∞)≤ |x−x∞||∇E(x)|−k|x−x∞|2/2.(2.3) CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 9 Proof. The conclusions of the lemma continue to hold under the relaxed hypothesis f(s)≥k, as is easily seen by replacing kwith k−1/n and taking a limit n→ ∞. Given x∈Rd, the function f(s) := E(x∞+sx−x∞ |x−x∞|)−E(x∞) satisfies the hypothesis f00(s)≥k. Setting s=|x−x∞|in the conclusion of the lemma, Cauchy-Schwarz yields the desired inequalities (2.2–2.3). ¤ For the reader familiar with Riemannian geometry, it is not hard to extend the results of this section to a C2function E:M−→ Ron a complete Riemannian manifold (M, h,i) satisfying the Hessian bound D2E≥kh,i. For example, (2.3) takes the form E(x)−E(x∞)≤dist(x, x∞)|∇E(x)|−kdist(x, x∞)2/2.(2.4) where dist(x, x∞) denotes arclength (geodesic) distance between xand x∞∈Minduced by the Riemannian metric h,i. Our primary task will be to extend the argument of Proposition 2.1 to the length space M=P2(Rd) of probability measures metrized by quadratic Wasserstein distance, to obtain optimal contraction rates under a range of degenerate convexity assumptions. Following Otto’s work [43], analogs of Propositions 2.2–2.4 have been explored in this setting by Agueh [1], Agueh, Ghoussoub, & Kang [2], Carrillo, J¨ungel, Markowich, Toscani & Unterreiter [22], Cordero-Erausquin, Gangbo & Houdre [27], Otto & Villani [44], and our parallel work [23]. In the classical case of linear diffusion with quadratic confinement (Example 1.1 with α=β=a= 1 and γ=b= 0), the relative energy E(ρ)−E(ρ∞) reduces to the Boltzmann entropy (5.1) of ρ=f2ρ∞, and |∇E(ρ)|2to its Fisher information. As explained in these references, the first inequality in (2.2) becomes Talagrand’s transportation inequality [50], the second the log-Sobolev inequality of Gross and others [35], while (2.4) becomes the HWI inequality of Otto and Villani [44]. In particular, (2.4) indicates how 2-uniform displacement convexity on a neighbourhood of x∞converts convergence in the weak metric dist2to convergence in relative entropy. From there it can often be converted to strong convergence in L1(Rd) via a Csiszar-Kullback inequality, as in [22, 43]. This helps to explain part of the interest in Wasserstein contraction rates. Although we were not aware of it at the time of first writing, analogs of Proposition 2.1 were explored simultaneously and independently by Ambrosio, Gigli and Savare [8], and for the heat equation in a Riemannian setting by Sturm and von Renesse [49], who showed the p-Wasserstein contraction / expansion rate for each p≥1 is given by the sharp lower bound k∈Rfor the Ricci curvature of the manifold. Refinements of this observation have been pursued by Sturm [46] [47] [48] and Lott & Villani [38]. 16 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI since sp ∈ Kx, so the convex function E(σs) has right derivative dE(σs) ds ¯¯¯¯s=0+≥ −h˙u0, pi.(3.16) Similarly, −˙v0∈∂Ev0and σs= expv0(1 −s)qimply dE(σs) ds ¯¯¯¯s=1−≤ h˙v0, qi.(3.17) Using (3.16–3.17) to estimate (3.15) yields d+ dt ¯¯¯¯0 dist2(ut, vt)/2 = dist(u0, v0)d+ dt ¯¯¯¯0 dist(ut, vt) ≤ −dE(σs) ds ¯¯¯¯ s=1− s=0+ ≤ −φ( dist(u0, v0)) dist(u0, v0), by φ-uniform convexity (3.9) of Ealong the geodesicσsof length dist(u0,v0). Time-translation invariance shows the same estimate must hold at any other time t=t0that we derived at t= 0. Thus when φ≥0 (resp. φ≤0), d+ dt ¯¯¯¯t0 Φ( dist(ut, vt)) = 1 φ( dist(ut0, vt0)) d+ dt ¯¯¯¯t0 dist(ut, vt)≤ −1 (resp. ≥) (3.18) holds at each instant t0∈I:= {t∈[0, T)|ut06=vt0}. In case φ≥0, the primitive equation (3.14) defines a continuously increasing function Φ : (0,∞)−→ Rin view of hypothesis (φ0) (3.4), but its limit Φ(0) = −∞ is unbounded due to the Lipschitz continuity of φnear φ(0) = 0 implied by (φ2). Thus the inverse Φ−1: (−∞,Φ(∞)) −→ Ris also a continuously increasing function. If φ≤0, (3.18) is reversed but Φ decreases monotonically from Φ(0) = +∞, and we may need to extend Φ−1(s) to s≤Φ(∞) by setting Φ−1(s) = +∞. In this case the growing bound (3.13) may only remain finite for a short time. Using hypothesis (φ1) (3.5), we obtain that this growth is no larger than exponential and thus, it remains finite for all times. Either way, Gronwall’s inequality completes the proof as long as I= [0, b)⊂[0, T). The only remaining possibility is that the relatively open subset I⊂[0, T) contains a non-empty connected component (a, b)⊂[0, T). We claim this cannot happen. To see why, observe that if φ≥0 then Gronwall’s inequality yields t+ Φ( dist(ut, vt)) non-increasing so s+ Φ( dist(us, vs)) ≥t+ Φ( dist(ut, vt)) (3.19) CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 17 for a<s<t<b. Letting s→ashows Φ(0) ≥t−a+ Φ( dist(ut, vt)), contradicting Φ(0) = −∞. On the other hand, if φ≤0, then (3.18–3.19) are reversed. Taking the limit s→acontradicts Φ(0) = +∞, to conclude the proof of the theorem. ¤ Example 3.13 (Exponential versus algebraic convergence). a) φ(x) = kx with k∈Rimplies Φ(x) = 1 klog xand Φ−1(y) = eky so (3.13) becomes dist(ut, vt)≤e−kt dist(u0, v0).(3.20) b) φ(x)=(k/r)xr+1 with k, r > 0implies Φ(x) = −1 kx−rand Φ−1(y) = (−ky)−1/r so (3.13) becomes dist(ut, vt)≤dist(u0, v0) (1 + tk distr(u0, v0))1/r .(3.21) Remark 3.14 (Rates of expansion).Theorem 3.12 covers semiconvex functionals as well as convex ones. Thus (3.20) with k < 0provides exponential control on the growth of separation between two initial conditions under the subgradient flow. In particular, taking u0=v0shows the time evolution defined by the flow is unique, when it exists. 4. Probability measures form a Riemannian length space As discussed in the introduction, we are interested in the evolution of probability measures verifying certain partial differential equations. Our objective is to formulate this evolution as a subgradient flow on a Riemannian length space. In this section we introduce the relevant Riemannian length space structure on subsets of the space of all Borel probability measures on Rd, i.e., P(Rd). To begin we recall the Kantorovich-Rubinstein-Wasserstein L2distance dist2(ρ, ρ0) [36, 53] between two measures ρ, ρ0∈ P(Rd): its square is defined as an infimum dist2 2(ρ, ρ0) := inf γ∈Γ(ρ,ρ0)ZRd×Rd|x−y|2dγ(x, y) (4.1) over the set Γ(ρ, ρ0) of joint measures γ≥0 on Rd×Rdwith left and right marginals ρand ρ0, respectively. It is not hard to see that dist2 satisfies the triangle inequality and makes P(Rd) a complete metric space [29, 32]. However dist2(ρ, ρ0) = +∞whenever one measure has finite second moment and the other does not, so henceforth we restrict our attention to the connected component P2(Rd) := ½ρ∈ P(Rd)|ZRd|x|2dρ(x)<+∞¾,(4.2) 18 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI itself a complete metric space on which dist2is finite. Let Pac(Rd) denote the set of Borel probability measures on Rdwhich are absolutely continuous with respect to Lebesgue. The intersection P2(Rd)∩ Pac(Rd) is denoted Pac 2(Rd). It is also easy to see that N:= P2(Rd) is a length space: the infimum (4.1) is attained, and the image ρs:= (πs)#γof the optimal joint measure γunder the map πs: (x, y)∈Rd×Rd−→ (1 −s)x+sy ∈Rdtraces out a minimal geodesic in P2(Rd) as s∈[0,1] ranges from zero to one. The notation (πs)#γis defined in (4.5). Although ρsis a measure on Rdand not on Rd×Rd, it can — apart from a dilation depending on s— be visualized as the projection of γonto a d-dimensional subspace {((1 −s)x, sx)|x∈Rd} ⊂ Rd×Rd. What is more subtle are the following facts established in McCann’s thesis [40], where these paths were first introduced and described from a different point of view under the name displacement interpolation: (i) Pac 2(Rd) is geodesically convex; (ii) a minimal geodesic is uniquely determined by its endpoints if either (or both) of them lie in Pac 2(Rd); (iii) in this case, the entire geodesic lies in Pac 2(Rd) except perhaps for its second endpoint. For the present point of view, the most relevant articulation and proof of (ii) is the one given by Carlen & Gangbo [19]. Taking N=P2(Rd) as our complete length space, the subset M= Pac 2(Rd) of absolutely continuous probability measures will carry our Riemannian length space structure. Here we recall the formal Riemannian structure introduced on Pac 2(Rd) by Otto [43], who first realized the connection between this structure and nonlinear diffusions as gradient flows. Although Otto used this connection in a purely formal manner to motivate detailed rate calculations in [43], for the theory developed hereafter it is necessary to state somewhat more precisely the nature of the tangent space, exponential mapping, and structure of M. According to Definition 3.8 we only need to define the exponential mapping over the subset M. Fix ρ∈M. Let spt (ρ) denote the smallest closed subset of Rdcontaining the full mass of ρ, and let Ωρ⊂Rddenote the interior of the convex hull of spt (ρ). We take Hρ:= H1,2(Rd, dρ)⊂C0,1 loc (Ωρ) to consist of those locally Lipschitz functions on Ωρwhose first derivative lies in the weighted space L2(Rd, dρ;Rd), modulo equivalence with respect to the seminorm hψ, ψiρ=ZΩρ|∇ψ|2dρ(x).(4.3) The local Lipschitz condition implies ∇ψ∈L∞ loc(Ωρ) can be interpreted equally well in either the a.e. or the distributional sense, as long as the CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 19 domain of integration in (4.3) is restricted to Ωρinstead of Rd. There is some arbitrariness in this definition; for the purposes which follow we may equally well choose to further restrict the space H1,2(Rd, dρ), also denoted H1,2 ρ, to consist only of those functions which can be expressed locally as a difference of convex functions on Ωρ. At any ρ∈ Pac 2(Rd), the tangent space TρMto Mis identified with the vector space H1,2 ρequipped with the inner product (4.3). The exponential map generates a curve s∈R−→ ρs∈Npassing through ρ0=ρin direction ψ∈ H1,2 ρdefined by imagining a collection of infinitesimally small particles comprising ρ, which evolve freely in time (both future and past) and having velocity profile ∇ψat time s= 0. More precisely, the Borel map F(x) := x+s∇ψ(x) is used to push forward the measure ρon Rdto yield expρsψ := [Id +s∇ψ]#ρ, (4.4) where by definition, the pushed-forward measure F#ρ∈ P(Rd) assigns mass F#ρ[K] := ρ[F−1(K)] (4.5) to each Borel set K⊂Rd. Observe ρs= expρsψ belongs to P2(Rd) by finiteness of the kinetic energy (4.3). Thus expρ:H1,2 ρ−→ Nis well-defined, and surjective as a consequence of the Brenier/McCann theorem [18, 39], which associates to each ρ∈ P2(Rd) and ρ0∈ P2(Rd) a convex function ψ(x) + |x|2/2 on Rd (taking values in R∪ {+∞}) whose gradient pushes ρforward to ρ0[39]. This motivates the identification of the star-shaped set Kρ=½ψ∈ H1,2 ρ|Ψ(x) = 1 2|x|2+ψ(x) convex on Rd,∇Ψ#ρ∈ Pac 2(Rd)¾ (4.6) which allows us to verify the conditions over the exponential map necessary for Pac 2(Rd) to be a Riemannian length space : Proposition 4.1 (Wassersteindistance metrizes a Riemannian lengthspace). The absolutely continuous measures M=Pac 2(Rd)form a Riemannian length space metrized by dist2(ρ, ρ0). In particular, the squared Wasserstein distance is superdifferentiable on the product manifold M×M: letting ρsdenote the minimal geodesic joining ρ0=ρto ρ1=ρ0yields dist2 2(expρtψ, expρ0tψ0)≤dist2 2(ρ, ρ0)−2thψ0,dρs ds ¯¯¯¯1− iρ0−2thψ, dρs ds ¯¯¯¯0+ iρ+ 4t2, (4.7) or equivalently dist2 2(expρtψ, expρ0tψ0)≤dist2 2(ρ, ρ0)+2thψ0, ϕ0iρ0−2thψ, ϕiρ+ 4t2, 20 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI for each pair of unit tangent vectors ψ∈TρMand ψ0∈Tρ0M, where ϕ,ϕ0 are such that ρs= expρsϕ = expρ0(1 −s)ϕ0. Proof. Given ρ, ρ0∈ Pac 2(Rd), let γ0∈Γ(ρ, ρ0) denote the joint measure which achieves the infimum (4.1) defining the Wasserstein distance. This measure can also be expressed in the form γ0= (id ×(∇ϕ+id))#ρ= ((∇ϕ0+id)×id)#ρ0, where the functions ϕ(x) + x2/2 and ϕ0(y) + y2/2 are convex and Legendre transforms, according to Brenier’s theorem [18]; see also McCann [39] and Rachev & R¨uschendorf [45]. (The same theorem has a converse that we also require: every ˜ϕ∈ Kρgives rise to a ˜γ0achieving the Wasserstein distance dist2 2(ρ, (∇˜ϕ+id)#ρ) = h˜ϕ, ˜ϕiρ.) Our prescription for constructing minimal geodesics yields ρs:= [id +s∇ϕ]#ρ= [id + (1 −s)∇ϕ0]#ρ0; (4.8) indeed dist2 2(ρs, ρs+t) = |s−t|2hϕ, ϕiρ=|s−t|2dist2 2(ρ, ρ0)<+∞as in Example 3.2, and ρs∈ Pac 2(Rd) is absolutely continuous according to Proposition 1.3 of [40]. This shows Kρis star-shaped, ϕ∈ Kρ, and the exponential (4.4) maps Kρonto M, taking rays onto minimal geodesics as desired. Also ϕ0∈ Kρ0, and (4.8) shows the geodesic ρs= expρ0(1 −s)ϕ0 can be parameterized from the other end equally well, as required in the Riemannian length space definition 3.8. It remains only to establish (4.7), which will imply (3.11) to complete the proof. Given ψ∈ H1,2 ρand ψ0∈ H1,2 ρ0of unit norm, the map F(x, y)=(x+ t∇ψ(x), y +t∇ψ0(y)) on Rd×Rdcan be used to define a pushed-forward measure γt:= F#γ0via (4.5). Then γt∈Γ(expρtψ, expρ0tψ0), so (4.1) implies dist2 2(expρtψ, expρ0tψ0)≤ZRd×Rd|x0−y0|2dγt(x0, y0) =ZRd×Rd|x−y+t(∇ψ(x)−∇ψ0(y))|2dγ0(x, y) = dist2 2(ρ, ρ0)+ ZR2n 2thx−y, ∇ψ(x)−∇ψ0(y)i+t2|∇ψ(x)−∇ψ0(y)|2dγ0(x, y) ≤dist2 2(ρ, ρ0)+2t ZRd·h−∇ϕ(x),∇ψ(x)idρ(x)+h∇ϕ0(y),∇ψ0(y)idρ0(y)¸+4t2 = dist2 2(ρ, ρ0) + 2thϕ0, ψ0iρ0−2thϕ, ψiρ+ 4t2, yielding the proof of (4.7). ¤ Remark 4.2. We use expρin place of a coordinate chart covering M= Pac 2(Rd). Its star-shaped domain Kρof bijectivity is actually convex, a fact which is not central to our discussion, but is proven below since the convexity CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 21 of certain functions E(expρψ)on Kρis central to the analysis of Ambrosio, Gigli and Savare [7, 8], and can be addressed similarly. Unfortunately, even if the inner product space H1,2 ρhappens to be complete, we cannot define a Hilbert manifold using these atlases because the convex set Kρ⊂ H1,2 ρ of bijectivity for the exponential map is too short in various directions to contain an open neighbourhood of the origin. For example, there are many measures (including the Gaussian) for which ψ(x) = t|x|belongs to Kρif and only if t≥0. Exponentiating from such a measure in the −|x|direction produces a curve which is not length minimizing, even locally. Proof of convexity of Kρ.Given ψ, ψ0∈ Kρ⊂ H1,2 ρ, observe ψt:= (1−t)ψ+ tψ0is in H1,2 ρ. Clearly Ψt(x) = ψt(x)+|x|2/2 is convex for t∈[0,1] since Ψ0 and Ψ1are. We need only deduce the absolute continuity of (∇Ψt)#ρfrom its absolute continuity at the endpoints t= 0 and 1. Let Xt⊂Rddenote the set of Lebesgue points of ρwhere Ψtadmits a second order Taylor expansion in which the Hessian D2Ψt(x) of Aleksandrov is invertible. From Theorem 4.4 of McCann [40], this Borel set carries the full measure of ρwhen t= 0 or 1; the same follows for t∈[0,1] since X0∩X1⊂Xt. Proposition 4.2 of [40] states that the (symmetric) Lebesgue density of the measure (∇Ψt)#ρ at ∇Ψt(x) is finite for x∈Xt, being given by ρ(x)/det[D2Ψt(x)]. On the contrary, almost everywhere with respect to the singular part of (∇Ψt)#ρits symmetric Lebesgue density would be infinite. Since the Borel set ∇Ψt(Xt) has full mass for (∇Ψt)#ρ, we conclude the measure has no singular part, and must be absolutely continuous as desired. ¤ For technical reasons, it is convenient in certain applications to be able to restrict our attention to compactly supported measures. The following corollary to Remark 3.9 shows that probability measures on any convex set Ω⊂Rdalso form a Riemannian length space Pac(Ω) := {ρ∈ Pac 2(Rd)|ρ[Rd\Ω] = 0}. Corollary 4.3 (Geodesic convexity of measures on a convex domain).Let Ω⊂Rdbe convex, and (M, dist) = (Pac 2(Rd),dist2). Then Pac(Ω) forms a geodesically convex subset of M, and hence a Riemannian length space. Proof. Let ρ, ρ0∈ Pac(Ω). Carlen and Gangbo [19, Theorem 2.2] assert the existence of a unique minimal geodesic s∈[0,1] −→ N=P2(Rd) joining ρ0=ρto ρ1=ρ0. We claim that ρs∈ Pac(Ω) ⊂M. The previous proof (4.8) asserts that ρs= [id+s∇ϕ]#ρ∈Mis absolutely continuous and given in terms of a function ϕ∈ Kρ. To see 0 = ρs[Rd\Ω] = ρ[(id +s∇ϕ)−1(Rd\Ω)],(4.9) 22 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI observe that it holds for s= 0 and s= 1 by hypothesis. This means for ρ-a.e. x∈Rdthat x∈Ω and x+∇ϕ(x)∈Ω. Convexity of Ω implies (1 −s)x+s(x+∇ϕ(x)) ∈Ω on the same set where ρhas full measure, thus extending (4.9) to all s∈[0,1]. This proves Pac(Ω) is geodesically convex. That M0:= Pac(Ω) inherits the Riemannian length space structure from Pac 2(Rd) now follows from Proposition 4.1 by Remark 3.9. ¤ Remark 4.4 (Closability and smooth densities).If |log ρ(x)|is bounded throughout Ωρ⊂Rd, then H1,2 ρmay form a Hilbert space. For more general densities, closability of the Dirichlet form (4.3) can be a delicate question [5] which we prefer to avoid; this is why we do not require that the inner product space H1,2 ρbe a Hilbert space generally. 4.1. Differentiable curves on M.Finally, we make contact with Otto’s formalism [43] by pointing out that the charts described above correspond to normal coordinates around the point ρ∈ Pac 2(Rd), in the sense that the metric assumes the canonical form (4.3). If one chooses to parameterize Pac 2(Rd) by some other set of coordinates near ρ, a corresponding linear transformation is induced on the components ψ(x) of each tangent vector. In particular, the linear transformation ψ∈ H1,2(Rd, dρ)−→ div [ρ∇ψ]∈ D0(Ωρ) (4.10) of the tangent space has a distinguished role, since formally at least, the geodesic path ρtdefined by (4.4) satisfies the instantaneous transport equation ∂ρt ∂t ¯¯¯¯t=0 + div [ρ∇ψ] = 0; (4.11) for ψsemiconvex, (4.11) is actually a consequence at a.e. point where ρis differentiable of the Monge-Amp`ere equation ρ0(x) = ρt(x+t∇ψ(x)) det[I+ tD2ψ(x)] [40]. Given any non-geodesic path ρt∈Mthrough ρ0=ρsmooth enough that the elliptic problem (4.11) has a solution ψ(x)∈ H1,2(Rd, dρ), Otto asserts that the solution ψrepresents a tangent vector to the curve in normal coordinates. Notice that ∂ρ/∂t has zero total mass, so any boundary conditions for the elliptic problem must ensure no net flux. Since ρ∈ Pac 2(Rd) is assumed fixed, the solution ψthe depends linearly on (∂ρ/∂t)t=0 as desired. On a smooth bounded domain Ω ⊂⊂ Rd, the following lemma gives sufficient conditions for differentiability of such a curve, and identifies its tangent vector. The outward unit normal to the domain boundary is denoted by νΩ(x) at x∈∂Ω. Lemma 4.5 (Tangent to a smooth curve; c.f. [43]).Fix Ω⊂⊂ Rdbounded smooth domain. Suppose a C2smooth function ψt(x) := ψ(t, x)∈Rand a CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 23 smooth C2curve of probability densities ρt(x) := ρ(t, x)≥0on [0, T]×Ω are related by the transport equation and no-flux (Neumann) condition ∂ρ ∂t + div [ρ∇ψ] = 0 on [0, T]×Ω, ∇ψ(t, x)·νΩ(x) = 0 on [0, T]×∂Ω.(4.12) Then t→ρtis a differentiable curve in Pac 2(Rd). A tangent vector ˙ρtto the curve at t∈[0, T]is given by ψt∈ H1,2 ρt. Proof. Without loss of generality, we’ll establish right differentiability of the curve at t= 0, and show ψ0∈ H1,2 ρ0is a tangent vector. To compare ρtwith the geodesic ˜ρt:= expρ0t∇ψ0, integrate dXt(x) dt =∇ψt(Xt(x)) X0(x) = x(4.13) to find the one-parameter family of diffeomorphisms Xt: Ω −→ Ω generated by ∇ψt. The Wasserstein distance between ρtand ˜ρt= (Yt)#ρ0is estimated using the joint measure γt:= (Xt×Yt)#ρ0constructed from Xtand Yt(x) := x+t∇ψ0(x). Note that Taylor’s theorem and (4.13) imply Xt(x) = Yt(x)+ O(t2); the C2smoothness of ψ(t, x) allows the error term to be estimated uniformly in x∈Ω as t→0. By definition (4.1), dist2 2(ρt,˜ρt)≤ZRd×Rd|x−y|2dγt(x, y) =ZΩ|Xt(z)−Yt(z)|2dρ0(z) =O(t4), with is more than the definition (3.11) of differentiability requires. [Mere continuity of dX/dt on [0, T]×Ω is enough to yield dist2(ρt,˜ρt) = o(t)]. ¤ For certain applications, we will also be interested in proving differentiability of paths of measures defined on the whole space Rd. Lemma 4.6 (Tangent to a smooth curve in Rd).Suppose a C2smooth function ψt(x) := ψ(t, x)∈Rand a smooth C2curve of probability densities ρt(x) := ρ(t, x)≥0on [0, T]×Rdare related by the transport equation ∂ρ ∂t + div [ρ∇ψ] = 0 on [0, T]×Rd. Assume that |∇ψt(x)| ≤ C0(1 + |x|) (4.14) 24 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI and ZRd|x|2ρ(t, x)dx < C0(4.15) for each t∈[0, T ]. Then t→ρtis a differentiable curve in Pac 2(Rd). A tangent vector ˙ρtat t∈[0, T ]is given by ψt∈ H1,2 ρt. Proof. Let us follow the notation of the previous lemma. We shall denote by Cvarious positive constants only depending on C0and T. Since ∇ψtis of class C1and linearly growing at ∞, standard classical results of ODE’s ensure the global existence in [0, T], uniqueness and regularity of the solutions of the initial value problem (4.13). Therefore, the family Xtof C1diffeomorphisms is well defined for any t∈[0, T] and there is no difficulty in deducing from the transport equation that Xt#ρ0=ρt for any t∈[0, T]. Again, the Wasserstein distance between ρtand ˜ρt=Yt#ρ0is estimated using the joint measure γt:= (Xt×Yt)#ρ0constructed from Xtand Yt(x) := x+t∇ψ0(x). Note that Xt−Yt=Zt 0 [∇ψs(Xs)−∇ψ0(X0)] ds, and together with our bounds on ∇ψthis implies in particular (1) |Xt−Yt| ≤ Ct(1 + |x|), (2) (Xt−Yt)/t converges towards 0 as t→0, for all x. By Lebesgue’s dominated convergence theorem, it follows that lim t→0 1 t2Z|Xt(z)−Yt(z)|2dρ0(z) = 0, which is what the definition (3.11) of differentiability requires. ¤ Remark 4.7 (Differentiable curves defined by gradient flows).The previous lemma remains valid under less stringent conditions on the growth of ∇ψin x, provided extra integrability assumptions on ρtare satisfied. For instance, (4.14) can be replaced by any hypothesis implying a well-defined flow map for the ODE system (4.13) in the whole interval [0, T]for any T > 0. We refer to this as global existence for (4.13). Linear growth of the function defining an ODE system is the simplest assumption implying global existence of (4.13). The use of a Liapunov functional L(x)is one of the standard tools for proving global existence for (4.13). In particular, any autonomous gradient-flow, i.e., dXt(x) dt =∇ψ(Xt(x)) CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 25 has a Liapunov functional given by L(x) = −ψ(x). Coercivity of L(x), i.e. boundedness of its sublevel sets, is enough to ensure a well-defined family of diffeomorphisms Xtfor any t∈[0, T ]. Therefore, in this case the growth of L(x) = −ψ(x)when |x| → ∞ need not be restrictive. Nonetheless, hypothesis (4.15) needs to be strengthened by suitably bounded moments related to the growth of ∇ψ(x)at infinity. In particular, a set of hypotheses for autonomous gradient-flows ensuring the conclusion of the previous lemma is ψ(x)≃ −A|x|k, when |x|→∞,|∇ψ(x)| ≤ C(1 + |x|k−1)with k≥2and uniform time estimates on the 2(k−1)th-moment of the densities ρt. In the case of non-autonomous gradient-flow ODE’s systems, i.e., dXt(x) dt =∇ψt(Xt(x)), the conditions on L(t, x) = −ψt(x)which imply global existence are −dψt(x) dt −|∇ψt(x)|2≤0 (4.16) and −ψt(x)≥ −˜ ψ(x)for any t∈[0, T]with −˜ ψ(x)coercive. Therefore, a set of hypotheses for non-autonomous gradient-flows ensuring the thesis of previous lemma can also be written in the same spirit as for the autonomous case by adding to (4.16) uniform bounds in time for the gradient of ψt(x) and suitable uniform time estimates of moments of ρt. However, these assumptions are difficult to meet in applications. 5. Energy functionals on M In this section we turn to the model for granular media which motivates the foregoing theory. The energy functional E(ρ) that we consider is a sum of three terms: E(ρ) = A(ρ) + B(ρ) + C(ρ) =ZRd A(ρac(x)) dx + ZRd B(x)dρ(x) +1 2ZRd×Rd C(x−y)dρ(x)dρ(y), (5.1) which can be defined on Pac(Rd), though we only need it on M=Pac 2(Rd). Here ρac denotes the Radon-Nikodym derivative of ρwith respect to Lebesgue measure. Let us first clarify the assumptions over each of these three terms. (A1) The internal energy A(%) is lower semicontinuous,A(0) = 0 and λ7−→ λdA(λ−d) is convex nonincreasing on λ∈(0,∞). (5.2) It follows that A(%) is proper, continuous and convex throughout [0,∞). Also, in terms of the pressure function P(%) := A0(%)%−A(%), (5.2) becomes 32 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI Let X⊂Ω denote the set where Ψ(x) can be approximated to second order by a quadratic polynomial; this set has full measure according to Aleksandrov’s theorem. Define vs(x) = det[(1 −s)I+sD2Ψ(x)] at x∈X. For each s < 1 the monotone change of variables theorem [40, Theorem 4.4] yields E(expρsψ) = ZX⊂Ω Aµρ(x) vs(x)¶vs(x)dx. (5.13) We shall shortly justify interchange of the integral Eρ(ψ) = lim s→0ZX A(ρ(x)/vs(x)) vs(x)−A(ρ(x)) sdx (5.14) with the limit lim s→0 A(ρ/vs)vs−A(ρ) s= [A(ρ)−ρA0(ρ)]∂vs ∂s ¯¯¯¯s=0 =−P(ρ)∂vs ∂s ¯¯¯¯s=0 . Let us first assume the legitimacy of this interchange, to complete the proof. Note that ∂vs(x)/∂s|s=0 = tr D2ψ(x) for each x∈X. Now the convexity of A(%) with A(0) = 0 yield A(%)≤%A0(%), and convexity of Ψ implies the distributional Laplacian ∆Ψ is a non-negative Radon measure on Ω with tr D2Ψ|Xas its absolutely continuous part. Thus Eρ(ψ) = −ZX P(ρ)tr D2ψ dx ≥ −ZΩ P(ρ)∆ψ dx =ZΩh∇P(ρ),∇ψidx −Z∂Ω P(ρ)∂ψ ∂ν dHd−1(x) ≥ZΩh∇A0(ρ(x)),∇ψ(x)iρ(x)dx =: hϕ, ψiρ. Here the last inequality follows from P(ρ)≥0, (5.12), and the identity P0(ρ) = ρA00(ρ). Finally, let us justify the exchange of the integral with the limit in (5.14). As in the proof of [40, Theorem 2.2], hypothesis (5.2) implies the integrand of (5.13) is convex as a function of s∈[0,1]. It follows that the integrand in (5.14) is non-decreasing on s∈[0,1]. Therefore, the corresponding integrands of (5.14) form a non-increasing sequence of functions as sdecreases to 0+verifying for 0 < s < 1/2 −P(ρ)tr D2ψ= [A(ρ)−ρA0(ρ)][(tr D2Ψ) −d] CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 33 ≤A(ρ/vs)vs−A(ρ) s ≤2A¡ρ/v1/2¢v1/2−2A(ρ) ≤2A¡2dρ¢/2d−2A(ρ) for all x∈X. Since A(ρ(x)) and P(ρ(x)) are C1(Ω) on a compact domain, tr D2ψ≥ −d, and ZX tr D2ψ dx ≤ZΩ ∆ψ dx =Z∂Ω ∂ψ ∂ν dHd−1(x) <Hd−1(∂Ω) sup x,y∈Ω|x−y|. We now have L1(Ω) bounds above and below throughout X⊂Ω and thus, we deal with non-increasing sequences of integrable functions with bounded integrals. Lebesgue’s dominated convergence theorem completes the justification. ¤ Let us verify the same result for the other two terms of the energy functional. Lemma 5.8 (Friction subgradient).Given ρ∈ Pac(Ω) on Ω⊂⊂ Rdand B:Rd−→ Rsemiconvex, ϕ:= B∈ H1,2 ρ(Ω) is a subgradient of the potential energy: ϕ∈∂Bρ⊂TρM. Proof. Semiconvexity of the integrand Bimplies lower semicontinuity and displacement semiconvexity of the functional E(ρ) := B(ρ) on Pac 2(Rd), by Remark 5.4 and Lemma 5.5. Thus it suffices to show ϕ∈(∂Eρ)0, according to Lemma 5.6. Semiconvexity also implies Bis Lipschitz since is locally finite, so both the function and its derivative are uniformly bounded on the bounded domain Ω. A straightforward application of Lebesgue’s dominated convergence theorem proves that Eρ(ψ) = lim s→0ZΩ B(x+s∇ψ(x)) −B(x) sdρ(x) =ZΩ lim s→0½B(x+s∇ψ(x)) −B(x) s¾dρ(x) =ZΩh∇B(x),∇ψidρ(x) =: hϕ, ψiρ. Thus ϕ:= B∈(∂Eρ)0.¤ 34 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI Lemma 5.9 (Collision subgradient).Given ρ∈ Pac(Ω) on Ω⊂⊂ Rdand C:Rd−→ Rsemiconvex, ϕ:= ρ∗C∈ H1,2 ρ(Ω) is a subgradient of the interaction energy: ϕ∈∂Cρ⊂TρM. Proof. Semiconvexity of the integrand Cimplies lower semicontinuity and displacement semiconvexity of the functional E(ρ) := C(ρ) on Pac 2(Rd), by Remark 5.4 and Lemma 5.5. Thus it suffices to show ϕ∈(∂Eρ)0, according to Lemma 5.6. Again Cis locally Lipschitz, so its derivative is uniformly bounded in the domain Ω ⊂⊂ Rd, and it is straightforward to prove using Lebesgue’s dominated convergence theorem that ϕ:= ρ∗C∈W1,∞(Ω) with ∇ϕ=ρ∗∇Cand Eρ(ψ) = = lim s→0 1 2ZΩ×Ω C[x−y+s(∇ψ(x)−∇ψ(y))] −C(x−y) sdρ(x)dρ(y) =1 2ZΩ×Ω lim s→0½C[x−y+s(∇ψ(x)−∇ψ(y))] −C(x−y) s¾dρ(x)dρ(y) =ZΩh∇C(x−y),∇ψ(x)idρ(x)dρ(y) =: hϕ, ψiρ. ¤ 5.3.2. Subdifferentiability of energies in Rd.The treatment of subdifferentiability for the whole space problem poses new challenges based on the need to control the behaviour of the densities and tangent vectors at +∞. Lemmas 5.8 and 5.9 are easily generalized to Rdprovided we restrict to a suitable geodesically convex subset of Pac 2(Rd). In fact, let us add hypotheses which further restrict the confinement and interaction potentials: (B2)B:Rd−→ Ris semiconvex on Rdand |∇B(x)| ≤ RB(1 + |x|αB); (C2)C:Rd−→ Ris semiconvex on Rdand |∇C(x)| ≤ RC(1 + |x|αC). (5.15) with αB, αC≥1 and RB, RC>0. Note that the convexity of |x|2αBand |x|2αCfor αB, αC≥1 implies that the functionals ˜ B(ρ) := ZRd|x|2αBdρ(x) and ˜ C(ρ) := ZRd×Rd|x−y|2αCdρ(x)dρ(y), CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 35 are displacement convex defined on Pac 2(Rd) by Lemma 5.2. Therefore, the set M0:= {ρ∈ Pac 2(Rd)|˜ B(ρ)<+∞and ˜ C(ρ)<+∞} is a geodesically convex subset of Pac 2(Rd) and hence, it inherits the structure of length space and the star-shaped set K0 ρof the tangent space is restricted to those tangent vectors joining measures in M0and therefore, lying inside the subset M0. Let us point out that hypotheses (5.15) are verified for the most relevant cases in applications: namely, power-like radial potentials. Lemma 5.10 (Friction subgradient in Rd).Given ρ∈M0,ϕ:= Bis a subgradient of the potential energy: ϕ∈∂Bρ⊂TρM0. Proof. Following the same proof as in Lemma 5.8, we need just to justify the application of Lebesgue’s dominated convergence theorem to interchange the limit s→0 and the integral over Rd. Since Bis locally Lipschitz, we can estimate the integrand as follows ¯¯¯¯ B(x+s∇ψ(x)) −B(x) s¯¯¯¯≤ ≤RB|∇ψ(x)|max 0≤s≤1(1 + |x+s∇ψ(x)|αB)≤ ≤RBmax{1 + |x|αB,1 + |x+∇ψ(x)|αB}|∇ψ(x)| by convexity of |x|αB. Since both ends of the geodesic lie on M0then the right-hand side of the inequality is integrable with respect to ρand thus, we have L1(Rd) control uniformly in s.¤ Subdifferentiability of the collision functional on all of Rdfollows a similar argument to previous lemma, so we omit the proof. Lemma 5.11 (Collision subgradient in Rd).Given ρ∈M0,ϕ:= ρ∗Cis a subgradient of the interaction energy: ϕ∈∂Cρ⊂TρM0. For the subdifferentiability of the entropy functional we will require some additional smoothness hypotheses on the density ρ∈M0and the energy density A. Assume (A2) A∈C2(0,∞)∩C[0,∞) satisfies A(0) = 0, A00(ρ)>0 and P(%)≥0 and P0(%) %1/2is integrable in (0,1).(5.16) Thus, the primitives P(%) and Q(%) of the differential equations P0(%) = %A00(%) and Q0(%) = %1/2A00(%) define diffeomorphisms on (0,∞). Moreover, 36 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI assumption (5.16) allows us to normalize, so that P(0) = Q(0) = 0. Indeed P(%) := %A0(%)−A(%) and Q(%) := Z% 0 s1/2A00(s)ds =Z% 0 P0(s) s1/2ds, (5.17) where the last integral converges due to (5.16). For example if m > 1/2 then A(%) = (%m−%)/(m−1), P(%) = %m, Q(%) := 2m 2m−1ρ(2m−1)/2. Lemma 5.12 (Integration by parts in the whole space for positive densities).Fix ρ∈ Pac 2(Rd)∩C∞(Rd)a positive density ρ > 0in Rd. Assume ϕ, ψ ∈ H1,2 ρ(Rd)where Ψ(x) = ψ(x) + |x|2/2is a convex function on Ωρ, ϕ(x) := A0(ρ(x)) ∈L2(Rd)and P(ρ)2ρ−1∈L1(Rd). If Asatisfies (A2) and D2ψ(x)denotes the Hessian of ψin the a.e. sense of Aleksandrov, then ZRd P(ρ(x))tr D2ψ(x)dx ≤ −ZRdh∇P(ρ(x)),∇ψ(x)idx. (5.18) Proof. Take the sequence of bump functions χn(x) = χ(x/n), n∈N, where χ(x)∈C∞ c(Rd) is supported in the unit ball being unity on B1/2(0). Thus, we have ∇χn(x) bounded uniformly in nin L∞(Rd) and being unity on Bn/2(0) with support inside Bn(0). Since ρ > 0 and smooth, the same is true for P(ρ), and we can use that 0 ≤tr D2Ψ≤∆Ψ where ∆ denotes the distributional Laplacian to get ZRd P(ρ)χntr D2Ψdx ≤ ≤ZRd P(ρ)χn∆Ψ dx =−ZRd χnh∇P(ρ),∇Ψidx −ZRd P(ρ)h∇χn,∇Ψidx. Since by hypotheses we have that ϕ, Ψ∈ H1,2 ρ(Rd) and we can rewrite the first term of the right-hand side as h∇P(ρ),∇Ψi=h∇Q(ρ), ρ1/2∇Ψi= hρ1/2∇A0(ρ), ρ1/2∇Ψi, then the dominated convergence theorem proves that ZRd χnh∇P(ρ),∇Ψidx →ZRdh∇P(ρ),∇Ψidx when n→ ∞. To pass to the limit in the second term, we first notice that ∇χn(x) is bounded uniformly in nin L∞(Rd) and converges pointwise to zero when n→ ∞. Moreover, by the assumptions and H¨older’s inequality |∇Ψ|P(ρ)∈L1(Rd), and thus, the second term vanishes when n→ ∞ by Lebesgue dominated convergence theorem. Application of Fatou’s lemma finally results in the desired inequality for Ψ. CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 37 To produce the same inequality for ψ= Ψ −|x|2/2, it is enough to check that (5.18) becomes an equality when Ψ(x) = |x|2/2. The only source of possible inequality is the relation between the Aleksandrov and distributional Laplacians, but in this particular case, they are equal tr D2Ψ(x) = ∆Ψ(x) = d. This concludes the proof of the lemma. ¤ Finally, the proof of Lemma 5.7 can be applied without any change to obtain the same conclusion based on the integration by parts inequality proved in previous Lemma. Lemma 5.13 (Entropy subgradient in Rd).Assume Asatisfies (A2) and ρ∈M0satisfies the additional hypotheses of Lemma 5.12. Moreover assume that A(2dρ)∈L1(Rd). Then ϕ(x) := A0(ρ(x)) ∈ H1,2 ρ(Ω) is a subgradient of the entropy (5.1):ϕ∈∂Aρ⊂TρM0. Let us remark that an additional assumption like A(2dρ)∈L1(Rd) — or at least A(ηρ)−ηA ◦ρ∈L1(Rd) for some η > 1 — is needed to have the L1bound from above in the proof of Lemma 5.7. When Aobeys an Orlicz condition, as in the homogeneous case, this hypothesis is implied by finiteness of the entropy. In more suggestive notation, we have shown that the variational derivative δE/δρ ∈ H1,2 ρgiven by δE δρ (ρ(x)) = A0(ρ(x)) + B(x)+(ρ∗C)(x) is a subgradient δE/δρ ∈∂Eρat any ρ∈M0with the additional smoothness assumptions of Lemmas 5.10–5.13. 6. Application to granular media The goal of this final section is to apply the strategy developed above to obtain contractivity properties of weak solutions of the family of PDEs: ∂ρ ∂t =∇··ρ∇µA0µdρ dx¶+B+C∗ρ¶¸, under assumptions on A,Band C, to be specified below. In order to apply the methodology of Section 3 to the framework of probability densities endowed with the Euclidean Wasserstein distance, we are forced to work in a smooth setting due to the differentiability structure we imposed in Section 4 and the analysis of the differentiability of curves from Lemmas 4.5 and 4.6. Therefore, our approach in this section is to show these contractivity properties in smooth situations (Subsection 6.1) by a direct application of our main theorem of decay rates for subgradient flows in length spaces, Theorem 3.12. These results will be generalized for weak solutions of this family of PDEs by approximation arguments. Therefore, a final step to obtain our general rates theorems needs to be done with our 38 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI strategy, which depends on particular cases of the general family of PDEs (1.2). In fact, we elect to write our theorems in an abstract framework by using a concept of approximable solutions (Subsection 6.2) in order to clarify the assumptions required from the approximations. As a consequence, we will use several approximation procedures from bounded velocity domains already existing in the literature to generalize the contractivity properties of the distances to weak solutions (Subsection 6.3). We have already seen that several technical problems appear for nonnegative solutions in whole of Rdas compared to strictly positive solutions on bounded domains, particularly, in the concept of differentiability of curves in Section 4 and additional hypotheses on the potentials (5.15) to have a well defined subgradient flow. In the most relevant case of linear diffusion, we will perform a different approximation procedure in order to overcome the challenge of ensuring the approximate solutions have fixed center of mass, as required for uniform displacement convexity of C(ρ) (Subsection 6.4). Before proceeding to this program in the next subsections, let us clarify the notion of weak solution [23] of the family of PDEs (1.2) we will deal with. The basic assumptions we make on the potentials and the diffusion are conditions (5.2) and (5.3). Let R+:= (0,∞) and R+ 0:= [0,∞). Definition 6.1 (Weak solution: linear diffusion).In the linear diffusion case, we will say that ρ∈C(R+ 0;P2(Rd)) is a weak solution of (1.2) if ∇C∗ρ∈L∞ loc(R+×Rd), such that for all T > 0and smooth, compactly supported test-functions η∈ D(Rd), ZRd η dρT−ZRd η dρ0=ZT 0ZRd ∆η dρtdt−ZT 0ZRd∇η·∇(B+C∗ρt)dρtdt. (6.1) Definition 6.2 (Weak solution: nonlinear diffusion).In the nonlinear diffusion case, given an initial data ρ0∈L1(Rd)we require ρ∈L∞ loc(R+;L1(Rd)) such that ∇C∗ρ∈L∞ loc(R+×Rd),∇P(ρ)∈L1 loc(R+×Rd,Rd), and ZRd η(0, x)ρ0(x)dx +Z∞ 0ZRd ∂η ∂t ρ(t, x)dx dt =Z∞ 0ZRd∇η·[∇P(ρt) + ∇(B+C∗ρt)ρt(x)] dx dt, (6.2) holds for all compactly supported test-functions η∈ D([0,∞)×Rd). Corresponding notions of weak solutions to equation (1.2) can be written in bounded domains Ω with no-flux boundary conditions assuming that that CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 39 the test functions are in η∈ D(R+×Ω) with supp(η)∩({0} × ∂Ω) = ∅. Weak solutions of (1.2) were constructed for particular cases under suitable additional hypotheses on A,Band Cin [43, 22, 23]. We will say that the diffusion is degenerate if A0(0) >−∞ and non-degenerate otherwise. 6.1. Smooth settings. We will give rate of convergence results for smooth solutions in two distinct situations: bounded velocity domains and the whole Rd. These results are direct consequences of the theory developed in the previous sections. Theorem 6.3 (Smooth setting: bounded velocity domain).Let A,Band Csatisfy conditions (5.2) and (5.3),B, C ∈C2(Rd),A∈C2(0,∞). Let ρ1, ρ2strictly positive C2probability density solutions of (1.2) in a smooth bounded convex domain Ω⊂Rdwith no flux boundary conditions ρ∇(A0(ρ) + B+C∗ρ)·νΩ= 0 on ∂Ω. Then, (1) If B:Rd−→ Rand C:Rd−→ Rare (semi)convex, say D2B(x)≥ βI and D2C(x)≥γI for a.e. x∈Rd, some β∈Rand γ≤0, then dist2(ρ1(t), ρ2(t)) ≤e−(β+γ)tdist2(ρ1(0), ρ2(0)) (6.3) holds for all t≥0. If the hypotheses are strengthened by insisting γ > 0, the stronger conclusion 6.3 will be true provided the centers of mass hxiρ1(t)=hxiρ2(t)= 0 remain equal for all t≥0. (2) Let φ(s) = (k/r)sr+1,k, r > 0, and assume the potentials satisfy one of the following two conditions: (i) B(x)is φ-uniformly convex on Rd, or (ii) C(x)is √2φ(·/√2)-uniformly convex on Rd, and the center of mass remains fixed hxiρn 1(t)=hxiρn 2(t)= 0 for all t≥0. Then for all t≥0the solutions ρ1(t)and ρ2(t)verify (3.21): dist2 2(ρ1(t), ρ2(t)) ≤dist2 2(ρ1(0), ρ2(0)) (1 + tk distr 2(ρ1(t), ρ2(0)))2/r . Proof. Given the smooth velocity fields ui=−∇ψiwith ψi=A0(ρi) + B+C∗ρi verifying ui·νΩ= 0 on the boundary of Ω for i= 1,2, we deduce that ∂ρi ∂t +∇·[ρiui] = 0, with ui·νΩ= 0 on the boundary for i= 1,2. Therefore, Lemma 4.5 ensures that both solutions are differentiable curves in M=Pac(Ω) with tangent 40 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI vectors given by ˙ρ1=−ψ1and ˙ρ2=−ψ2respectively. From Lemmas 5.7– 5.9 of subsection §5.3, we deduce that both curves are subgradient flows for the energy functional E(ρ) since ψi=A0(ρi) + B+C∗ρi=δE δρ (ρi(x)) for i= 1,2 and thus, −˙ρi∈∂Eρi. We have checked that our two curves on Mare differentiable and subgradient flows with respect to the energy functional E(ρ). A direct application of Theorem 3.12 imply the results stated. In fact, careful comparison of Examples 3.7 and 3.13 of Theorem 3.12 yield the conclusions of Example 1.1. For instance, if c > 0 then C(x) = γ|x|c+2/(c+ 2) is √2φ(s/√2)-uniformly convex with φ(s) = γ2−c/2sc+1. In the absence of friction β= 0 we need net momentum to vanish hxiρ1(t)=hxiρ2(t)= 0 since Lemma 5.2 shows C(ρ) to be φ-uniformly convex on Pac 0(Ω) ⊂ Pac 2,0(Rd) but not generally on Pac(Ω) ⊂ Pac 2(Rd). ¤ Remark 6.4 (Particular cases). •In case C= 0, smooth solutions are generic for strictly positive smooth initial data. Moreover, solutions are smooth for non-degene- rate diffusions for all positive times. Therefore, the hypotheses of previous theorem are not restrictive at all in the non-degenerate cases. •In the case of degenerate diffusions and C= 0, global weak solutions for initial probability densities in ρ(0) ∈L1(Ω) ∩L∞(Ω) were constructed in [14]. These weak solutions were obtained by approximating the degenerate diffusion by non-degenerate ones. •The fixed center of mass hypotheses mentioned in the previous theorem is satisfied for symmetric bounded velocity Ω = −Ωdomains and symmetric initial data ρ(0, x) = ρ(0,−x). In the whole space Rd, it will sometimes be convenient to work with smooth solutions of (1.2) which decay quickly as |x| → ∞. This is caused by the need to center the mass of the approximating solutions in order to extend the B= 0 cases of previous theorem to the whole space Rd. Theorem 6.5 (Smooth setting: Rd).Let A,Band Csatisfy conditions (5.2),(5.15) and (5.16),B, C ∈C2(Rd). Let ρ1, ρ2be curves of strictly positive C2probability densities decaying rapidly in Rdwhich satisfy (1.2), such that the integrability properties on ρ1, ρ2from Lemmas 5.12–5.13 are satisfied, and such that the velocity fields grow at most linearly as |x| → ∞. Then, the conclusions of Theorem 6.3 hold. CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 41 The preceding result follows analogously to Theorem 6.3 by using Lemma 4.6 to ensure the differentiability of curves in M=Pac 2(Rd) and Lemmas 5.10–5.13 to deduce the subgradient flow structure. 6.2. Approximable solutions: general theorems. In order to extend these rate of convergence results to weak solutions for non-degenerate and degenerate diffusions, we recall that weak solutions are typically constructed by taking the limit of solutions to a sequence of better behaved approximating problems involving non-degenerate diffusions on smooth bounded convex domains [14, 43, 22, 23]. Since the asymptotic rates of the preceding section apply directly to the approximating problems, we need only decide when the approximations are good enough for the asymptotic rate to survive the limit procedure. Sufficient conditions for this are summarized in the the notion of approximable solution introduced below. These conditions are not typically as strong as those which need to be imposed to construct weak solutions by the approximation method in the first place, so our conclusions apply a fortiori to the solutions constructed by Bertsch & Hilhorst, Otto, Carrillo et al, and in our earlier paper. Let us denote by ˆρthe extension of the function ρto Rdby setting ˆρ(x) = 0 for x6∈ Ω. Definition 6.6 (Approximable Solution).We say that ρis a approximable solution of (1.2) with initial data ρ0∈ Pac 2(Rd)if: (1) ρis a weak solution to (6.1)–(6.2). (2) There exists a sequence of smooth C2positive functions ρn, which are solutions of regularized equations of the form ∂ρn ∂t =∇·[ρn∇(A0 n(ρn) + Bn+Cn∗ρn)] ,(6.4) either in bounded smooth convex domains Ωn⊂Rdwith no flux boundary conditions ρn∇(A0 n(ρn) + Bn+Cn∗ρn)·νΩn= 0 on ∂Ωn, where An,Bnand Cnsatisfy the hypotheses of Theorem 6.3, or else in Rdbut decaying rapidly as |x| → ∞ and with velocity fields growing at most linearly as |x| → ∞ and with An,Bnand Cn satisfying the hypotheses of Theorem 6.5. Moreover, we assume that Bnand Cnare respectively φB n- and φC n-convex potentials with φB nand φC nconverging uniformly on compact subsets of of R+ 0to the respective moduli φBand φCof convexity of the potentials B and C. (3) ˆρnconverges towards the weak solution ρat least verifying ˆρn(t)→ρ(t)weakly in L1(Rd)a.e. t > 0.(6.5) 48 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI From Jensen’s inequality it results that Mα:= Rρ∞(1 + |x|α) satisfies Mα≤CM1−2/α α, in particular Mα≤Cα/2. Let h=ρ/ρ∞. Since ∇(log ρ) = ∇(log ρ∞) + ∇(log h), and since ∇(log ρ∞) = −∇(C∗ρ∞) satisfies the desired bound, it is sufficient for us to prove that |∇(log h)| ≤ R(6.11) and |∇C∗(ρ−ρ∞)| ≤ R. (6.12) Let ∂C =∂iCfor some index i. Since ρand ρ∞have the same mass and the same center of mass, we can write ∂C ∗(ρ−ρ∞) = ZRd ∂C(x−y) (ρ−ρ∞)(y)dy =ZRd [∂C(x−y)−∂C(x)−∇∂C(x)·y] (ρ−ρ∞)(y)dy. By Taylor’s formula and the uniform bound on D2∂C, we can bound this expression by RZRd|y|2|ρ−ρ∞|(y)dy, which is bounded by a uniform constant and (6.11) is proved. Let us proceed to estimate h. We will use the notation C=C∗ρand C∞=C∗ρ∞. Some tedious but easy computations lead to the equations ∂th= =∆h+(2∇log ρ∞+∇C)·∇h+ +(∆ log ρ∞+|∇log ρ∞|2+∇log ρ∞·∇C+ ∆C)h = ∆h+ (∇C−2∇C∞)·∇h+ (|∇C∞|2−∆C∞−∇C∞·∇C+ ∆C)h, then, with u= log h, ∂tu= ∆u+|∇u|2+(∇C−2∇C∞)·∇u+(|∇C∞|2−∆C∞−∇C∞·∇C+∆C). Let b:= ∇C−2∇C∞and c:= |∇C∞|2−∆C∞−∇C∞·∇C+∆C. Another calculation yields ∂t|∇u|2 2= ∆|∇u|2 2−kD2uk2+ +∇u·∇|∇u|2+ (2∇u+b)·∇|∇u|2 2+h∇b·∇u, ∇ui+∇c·∇u, CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 49 where k·kstands for the Hilbert-Schmidt norm. Let a:= 2∇u+b, we find (∂t−∆−a·∇)·|∇u|2 2≤(k∇bk+ 1)|∇u|2+|∇c|2. Our goal is to prove that |∇u|remains bounded on each interval [0, T ], knowing that it is bounded at time t= 0. If we manage to prove that both k∇bkand |∇c|are bounded, then the conclusion will follow by maximum principle. From our assumptions, ∇bis bounded. Let us estimate ∇c: the terms ∇∆Cand ∇∆C∞are bounded, so we only have to estimate ∇h|∇C∞|2−∇C∞·∇Ci=∇h∇C∞·∇C∗(ρ−ρ∞)i. And in view of the bounds |∇C∞| ≤ C(1 + |x|), |D2C∞| ≤ Cand (6.11) we only have to prove |D2C∗(ρ−ρ∞)| ≤ R 1 + |x|. Similarly to (6.11), if ∂2C=∂2 ijCfor some indices iand j, we can write |∂2C(x)| ≤ RZRd |y|3 1 + min(|x|,|x−y|)|ρ(y)−ρ∞(y)|dy. Assume |x| ≥ 1. The contribution of those y’s such that |y|≤|x|/2 to the integral above is bounded by R 1 + |x|ZRd|y|3|ρ(y)−ρ∞(y)|dy. On the other hand, by Chebyshev’s inequality, the contribution of those y’s such that |y| ≥ |x|/2 is bounded by R |x|ZRd |y|4 1 + min(|x|,|x−y|)|ρ(y)−ρ∞(y)|dy ≤R |x|ZRd|y|4|ρ(y)−ρ∞(y)|dy. We conclude that indeed |D2C(x)−D2C∞(x)| ≤ R/(1 + |x|), as was announced. ¤ The previous lemma allows us to apply Lemma 4.6 showing that the solution ρt(x) := ρ(t, x) is a differentiable curve on Pac 2(Rd) and therefore, that we have an approximation procedure in case Csatisfies the additional hypotheses: Cuniformly convex, D2Cbounded from above and |D3C| ≤ R/(1 + |x|) for a given constant R. Subdifferentiability of the energy functional results directly from Lemmas 5.10–5.13 due to the smoothness of the solution ρ(t) and the hypotheses on C. 50 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI Let us finally remark that this approximation procedure in Rdovercomes the difficulty of fixing the center of mass in the sequence of approximations on bounded domains. However, as a trade-off we need to face a new challenging problem, that is, to show the existence of a well-defined global-in-time flow map for the velocity field ∇ψt. In order to do so, we needed to impose hypotheses on the interaction potential Cfor which we are able to prove linear growth in xof the velocity field. We finally generalize the class of interaction potentials by a further approximation. Given a general interaction potential of the form C(x) = |x|c+2 with c > 0, we approximate it by a sequence of smooth interaction potentials Cnwith quadratic behavior at ∞, in such a way that the modulus of convexity φnof Cnconverges uniformly in compact subsets of R+ 0to the modulus of convexity φof C. This can be accomplished in this radial case by radial approximating functions obtained by smoothly truncating the second radial derivative near zero and outside a large interval [0, n]. In this way, we obtain potentials Cnsatisfying the quadratic growth at ∞for which: Cnuniformly convex, D2Cnbounded from above and |D3Cn| ≤ Rn/(1 + |x|). Therefore by Lemma 6.21, we ensure that our evolution defines smooth enough curves for the differentiability structure we need. We skip all the details since most of the work has already been done either in subsection 6.2 or in our companion paper [23] for the properties of the solutions and approximation. ¤ Let us finally remark that even in the presence of linear diffusion we have not been able to show exponential convergence towards equilibrium with degenerately convex interaction potential. This was done by the entropy method in our companion paper [23] and it remains an open problem to derive this result by means of measuring the convexity of the involved functionals in the approach just presented. Feasibility of the latter approach was explored in collaboration with NSERC summer undergraduate research assistant Tim Capes at the University of Toronto, who showed that an apriori bound on kρ(t)kL∞allows 2-uniform convexity of the entropy to be quantified, since the bound keeps us far away from the Dirac measures δx0where the convexity degenerates. Acknowledgements: The authors are pleased to thank an anonymous referee for providing a careful and thorough critique of the manuscript. JAC and CV acknowledge support from the European IHP network “Hyperbolic and Kinetic Equations: Asymptotics, Numerics, Applications” HPRNCT-2002-00282. JAC acknowledges the support from the Spanish DGIMCYT/FEDER project BFM2002-01710. RJM gratefully acknowledges the support of Natural Sciences and Engineering Research Council of Canada CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 51 grant 217006-03 RGPIN and NSF grants DMS 0074037 and 0354729. The authors also thank the Georgia Institute of Technology [RJM & CV], the Fields [JAC] and Pacific Institutes for the Mathematical Sciences [JAC & RJM], the Centre Recerca Matem`atica de Catalunya [RJM], and the Universities of Texas [JAC, RJM & CV] and Marne-la-Vall´ee [RJM] for the kind hospitality and excellent research environments they provided during various stages of this work. References [1] Agueh, M. Existence of solutions to degenerate parabolic equations via the Monge- Kantorovich theory. PhD Thesis, Georgia Institute of Technology, 2002. [2] Agueh, M., Ghoussoub, N., and Kang, X. Geometric inequalities via a general comparison principle for interacting gases. Geom. Funct. Anal. 14, (2004), 215-244. [3] Agueh, M., Ghoussoub, N., and Kang, X. The optimal evolution of the free energy of interacting gases and its applications. C. R. Math. Acad. Sci. Paris, Ser. I 337, (2003), 173-178. [4] Agueh, M. Asymptotic behavior for doubly degenerate parabolic equations. C. R. Math. Acad. Sci. Paris, Ser. I 337, (2003), 331-336. [5] Albeverio, S., R¨ ockner, M., Classical Dirichlet froms on topological vector spaces–closability and a Cameron-Martin formula. J. Funct. Anal. 88, (1990), 395- 436. [6] Alt, H. W., Luckhaus, S., Quasilinear elliptic-parabolic differential equations. Math. Z. 183, (1983), 311-341. [7] Ambrosio, L.A., Gigli, N., and Savar´ e, G., Gradient flows with metric and differentiable structures, and applications to the Wasserstein space. To appear in the proceedings of the meeting ”Nonlinear Evolution Equations” held in the Academy of Lincei in Rome. [8] Ambrosio, L.A., Gigli, N., and Savar´ e, G., Gradient flows in metric spaces and in the space of probability measures. Lectures in Mathematics, Birkh¨auser, 2005. [9] Ball, K., Carlen, E.A., and Lieb, E.H. Sharp uniform convexity and smoothness inequalities for trace norms. Invent. Math. 115, (1994), 463–482. [10] Bakry, E. and Emery, M. Diffusions hypercontractives in Sem. Probab. XIX LNM 1123 Springer, New York, 1985, pp 177–206. [11] Benedetto, D., Caglioti, E., and Pulvirenti, M. A kinetic equation for granular media. RAIRO Mod´el. Math. Anal. Num´er. 31, 5 (1997), 615–641. [12] Benedetto, D., Caglioti, E., Carrillo, J. A., and Pulvirenti, M. A nonmaxwellian steady distribution for one-dimensional granular media. J. Stat. Phys. 91, (1998), 979-990. [13] Benedetto, D., Caglioti, E., Golse, F., and Pulvirenti, M. A hydrodynamic model arising in the context of granular media Comput. Math. Appl. 38, (1999), 121–131. [14] Bertsch, M., and Hilhorst, D. A density dependent diffusion equation in population dynamics: stabilization to equilibrium. SIAM J. Math. Anal. 17 (1986), 863-883. [15] Biane, P., Speicher, R. Free diffusions, free entropy and free Fisher information. Ann. Inst. H. Poincar Probab. Statist. 37, (2001), 581-606. 52 JOS´ E A. CARRILLO, ROBERT J. MCCANN , C´ EDRIC VILLANI [16] Blower, G. Displacement convexity for the generalized orthogonal ensemble. J. Statist. Phys. 116, (2004), 1359-1387. [17] Bolley F., Brenier, Y. and Loeper, G. Contractive metrics for scalar conservation laws. To appear in Journal of Hyperbolic Differential Equations. [18] Brenier, Y. Polar factorization and monotone rearrangement of vector-valued functions. Comm. Pure Appl. Math. 44, (1991), 375–417. [19] Carlen, E. and Gangbo, W. Constrained steepest descent in the 2-Wasserstein metric. Annals Math. 157, (2003), 807-846. [20] Carrillo, J. A., Fellner, K. Long time asymptotics via entropy methods for diffusion dominated equations. Asymptotic Analysis 42, (2005), 29-54. [21] Carrillo, J. A., Gualdani, M. P., and Toscani, G. Finite speed of propagation for the porous medium equation by mass transportation methods. C. R. Math. Acad. Sci. Paris, Ser. I 338, (2004), 815-818. [22] Carrillo, J. A., J¨ ungel, A., Markowich, P. A., Toscani, G., and Unterreiter, A. Entropy dissipation methods for degenerate parabolic systems and generalized Sobolev inequalities. Monatsh. Math. 133 (2001) 1–82. [23] Carrillo, J. A., McCann, R. J., and Villani, C. Kinetic equilibration rates for granular media and related equations: entropy dissipation and mass transportation estimates. Rev. Matem´atica Iberoamericana 19 (2003), 1-48. [24] Carrillo, J. A., and Toscani, G. Asymptotic L1-decay of solutions of the porous medium equation to self-similarity. Indiana Univ. Math. J. 49, (2000), 113-141. [25] Carrillo, J. A., and Toscani, G. Wasserstein metric and large–time asymptotics of nonlinear diffusion equations. New Trends in Mathematical Physics, (In Honour of the Salvatore Rionero 70th Birthday), World Scientific, 2005. [26] Cordero-Erausquin, D. Some applications of mass transport to Gaussian-type inequalities. Arch. Rational Mech. Anal. 161, (2002),257–269. [27] Cordero-Erausquin, D., Gangbo, W., and Houdre, C. Inequalities for generalized entropy and optimal transportation. Recent advances in the theory and applications of mass transport, 73–94, Contemp. Math., 353, Amer. Math. Soc., Providence, RI, 2004. [28] Cordero-Erausquin, D., Nazaret, B., and Villani, C. A mass transportation approach to sharp Sobolev and Gagliardo-Nirenberg inequalities. Adv. Math. 182, (2004), 307–332. [29] Dudley, R.M. Probabilities and metrics - Convergence of laws on metric spaces, with a view to statistical testing. Universitet Matematisk Institut, Aarhus, Denmark, 1976. [30] Evans, L.C., and Gariepy, R.F. Measure Theory and Fine Properties of Functions. CRC Press, Boca Raton, 1992. [31] Gangbo, W., and McCann, R. J. Shape recognition via Wasserstein distance. Quart. J. Appl. Math. 4, (2000), 705-737. [32] Givens, C.R., and Shortt, R.M. A class of Wasserstein metrics for probability distributions. Michigan Math. J. 31, (1984), 231–240. [33] Gromov, M. Structures m´etriques pour les vari´et´es riemanniennes. J. Lafontaine and P. Pansu, eds. Cedic/Fernand Nathan, Paris, 1981. [34] Gromov, M. Metric Structures for Riemannian and non-Riemannian Spaces. J. Lafontaine and P. Pansu, eds. With appendices by S. Semmes. Birkhauser, Boston, 1999. [35] Gross, L. Logarithmic Sobolev inequalities. Amer. J. of Math. 97, (1975), 1061- 1083. CONTRACTIONS AND THERMALIZATION OF GRANULAR MEDIA 53 [36] Kantorovich, L.V., and Rubinstein, G.S. On a space of completely additive functions. Vestnik Leningrad. Univ. 13, (1958), 52–59. [37] Li, H., and Toscani, G. Long–time asymptotics of kinetic models of granular flows. Arch. Rat. Mech. Anal. 172, (2004), 407–428. [38] Lott J. and Villani C. Ricci curvature for metric-measure spaces via optimal transport. Preprint at http://www.math.lsa.umich.edu/∼lott/ [39] McCann, R. J. Existence and uniqueness of monotone measure-preserving maps. Duke Math. J. 80, (1995), 309–323. [40] McCann, R. J. A convexity principle for interacting gases. Adv. Math. 128, 1 (1997), 153–179. [41] McCann, R. J. Equilibrium shapes for planar crystals in an external field. Comm. Math. Phys. 195 (1998), 699–723. [42] McCann, R. J. Polar factorization of maps on Riemannian manifolds. Geom. Funct. Anal. 11 (2001) 589–608. [43] Otto, F. The geometry of dissipative evolution equations: the porous medium equation. Comm. Partial Differential Equations 26 (2001) 101–174. [44] Otto, F., and Villani, C. Generalization of an inequality by Talagrand and links with the logarithmic Sobolev inequality. J. Funct. Anal 173 (2001) 361–400. [45] Rachev, S. T., and R¨ uschendorf, L. Mass transportation problems. Probability and its Applications. Springer-Verlag, New York, 1998. [46] Sturm, K.-T. Convex functionals of probability measures and nonlinear diffusions on manifolds. To appear in J. Math. Pures Appl. [47] Sturm, K.-T. Generalized Ricci bounds and convergence of metric measure spaces. To appear in C. R. Acad. Sci. Paris S´er. I Math. [48] Sturm, K.-T. On the geometry of metric measure spaces. SFB Preprint #203, Bonn. http://www-wt.iam.uni-bonn.de/∼sturm/en/index.html [49] Sturm, K.-T. and von Renesse, M.-K. Transport inequalities, gradient estimates, entropy and Ricci curvature. To appear in Comm. Pure Appl. Math. [50] Talagrand, M. Transportation cost for Gaussian and other transport measures. Geom. Func. Anal. 6, (1996), 587-600. [51] Toscani, G. One-dimensional kinetic models of granular flows. RAIRO Mod´el. Math. Anal. Num´er. 34, 6 (2000), 1277–1291. [52] Villani, C. Topics in optimal transportation. Graduate Studies in Mathematics Vol. 58, Amer. Math. Soc, Providence, 2003. [53] Wasserstein, L.N. Markov processes over denumerable products of spaces describing large systems of automata. Problems of Information Transmission 5, (1969), 47–52. Jos´e A. Carrillo, ICREA and Dept. de Matem`atiques Universitat Aut`onoma de Barcelona, E-08193 Bellaterra, SPAIN. Robert J. McCann, Department of Mathematics, University of Toronto, Toronto, Ontario, M5S 3G3 CANADA C´edric Villani, U.M.P.A., Ecole Normale Sup´erieure de Lyon, F-69364 Lyon Cedex 07, FRANCE