Yue Su
email: yue.su@centralesupelec.fr

Nicolas Dupin
email: nicolas.dupin@univ-angers.fr

Sophie N Parragh
email: sophie.parragh@jku.at

Jakob Puchinger
email: jakob.puchinger@centralesupelec.fr

A Column Generation Approach for the Electric Autonomous Dial-A-Ride Problem

Keywords: Dial-A-Ride Problem, Electric Autonomous Vehicles, Labeling Algorithm

). In other words, we can always find a service start time T 1 ⩽ T 2 that consumes not more energy.

Introduction

The Electric Autonomous Dial-A-Ride Problem (E-ADARP) consists in scheduling a fleet of Electric Autonomous Vehicles (EAVs) to provide ride-sharing services for customers that specify their origins and destinations. The E-ADARP was first introduced by [START_REF] Bongiovanni | The electric autonomous dial-a-ride problem[END_REF] as a problem variant of the classical Dial-A-Ride Problem (DARP) defined in [START_REF] Cordeau | The dial-a-ride problem : models and algorithms[END_REF]. Although the E-ADARP shares some constraints of the classical DARP, the E-ADARP differentiates in two aspects : (i) a weighted-sum objective that minimizes both total travel time and total excess user ride time, and (ii) the employment of EAVs and a partial recharging policy ; The first aspects (ii) requires determining minimal-excess-time schedules for a feasible E-ADARP route, while the second aspect complicates the feasibility checking process as en-route recharging needs to be considered. Other E-ADARP features further increase the complexity of solving the problem. These features include (a) the elimination of maximum route duration constraints due to the autonomy of vehicles, (b) a minimum battery level that must be maintained at the end of the route, and (c) limited visits to each recharging station.

This paper proposes a highly-efficient labeling algorithm that can be embedded in the Column Generation (CG) approach to solve the E-ADARP. The original mixed integer programming (MIP) model is decomposed by the Dantzig-Wolfe decomposition into a master problem and pricing subproblems, which are formulated as Elementary Shortest Path Problems with Minimizing Excess Ride Time (hereafter ESPP-MERT). To solve the problem, we first introduce a fragment-based representation. Based on this representation, A novel approach is invoked to abstract fragments to arcs while ensuring excess-time optimality. We then construct a new graph that possesses all feasible routes of the original one, by enumerating all feasible fragments and connecting them with depots and recharging stations in a feasible way. Routes with negative reduced costs are generated by a customized labeling algorithm with strong dominance rules and complex REFs. To demonstrate the performance of the proposed labeling algorithm, we perform our labeling algorithm in CG on all existing E-ADARP instances and compare CG results to best-reported results in [START_REF] Bongiovanni | The electric autonomous dial-a-ride problem[END_REF][START_REF] Su | A Deterministic Annealing Local Search for the Electric Autonomous Dial-a-Ride Problem[END_REF].

The ESPP-MERT Description

The problem is defined on a complete directed graph G = (V, A), where V represents the set of vertices, and A denotes the set of arcs. V can be further partitioned into several subsets, i.e., V = N ∪ S ∪ O ∪ F , where N represents the set of all customers, S is the set of recharging stations, O and F denote the set of origin depots and destination depots, respectively. The set of all pickup vertices is denoted as P = {1, • • • , i, • • • , n} and the set of all drop-off vertices is denoted as

D = {n + 1, • • • , n + i, • • • , 2n}.
The union of P and D is N , i.e., N = P ∪ D. Each customer request is a pair (i, n + i) for i ∈ P and the maximum ride time for users associated with request i is assumed to be m i . A time window is defined on each node i ∈ V , denoted as [e i , l i ], in which e i and l i represent the earliest and latest time at which the vehicle starts its service, respectively. A load q i and a service duration s i are also associated with each node i ∈ V . For pickup node i ∈ P , q i is positive. For the corresponding drop-off node n + i, we have q n+i = -q i . For other nodes j ∈ O ∪ F ∪ S, q j and s j are equal to zero. Vehicles are assumed to be heterogeneous in terms of their maximum vehicle capacities (denoted as C k ) and homogeneous in terms of battery capacities (denoted as Q). The travel time on each arc (i, j) ∈ A is denoted as t i,j , and the battery consumption is denoted as b i,j . The recharging rate at charging facilities is α. To avoid the numerical problem when calculating time and energy, we define h i,j = b i,j /α to convert the battery consumption b i,j on arc (i, j) to the time needed for recharging this amount of energy. Similarly, we can also convert the current energy level to the time needed to recharge to this energy level. Let H denote the time required to recharge from zero to full battery capacity Q. Partial recharging is allowed while a vehicle visits recharging stations, and a minimum battery level γQ must be respected at destination depots, where we analyze three different γ values, namely, γ ∈ {0.1, 0.4, 0.7}. Higher values of γ result in more tightly constrained instances, which are harder to solve. The triangle inequality is assumed to hold for travel times and battery consumption.

The ESPP-MERT consists in finding the minimum-reduced-cost route while satisfying the following E-ADARP constraints :

1) Each route starts at an origin depot and ends at a destination depot ;

2) Pairing and precedence constraints for pickup and drop-off nodes ;

3) Time window constraints ; 4) Maximum user ride time constraints ; 5) Maximum vehicle capacity constraints ; 6) Maximum battery capacity constraints ; 7) Minimum battery level constraints at the destination depot ; 8) At-most-one visit to each recharging station ; 9) No passenger is onboard when the vehicle visits a recharging station ;

The reduced cost cω for a route ω is formulated as :

c ω - i∈V θ i,w λ i (1)
where c ω is the weighted-sum objective function value for a route ω, θ i,w is the binary coefficient that denotes whether node i is on ω, and λ i is the associated dual variable values.

Forward Labeling Algorithm for the ESPP-MERT

The most challenging point of solving the ESPP-MERT is maintaining excess-time optimality in the extension of labels. To handle this issue, we design a labeling algorithm by introducing a fragment-based representation in Section 3.1, which ensures the minimum excess user ride time being non-decreasing during the label extension. Our labeling algorithm includes two steps : From Section 3.2 to 3.3, we construct a new sparser network G sp by abstracting fragments to arcs. In G sp , we guarantee excess-time optimality for all arcs, and we prove that G sp possesses all feasible routes over the original graph G. Section 3.4 generates routes with negative reduced costs by an efficient labeling algorithm with strong dominance rules over G sp .

Representation of partial paths

To minimize the total excess user ride time for a partial path P, we first introduce the battery-restricted fragment as in Definition 1 :

Definition 1. Assuming that F = (i 1 , i 2 , • • • , i k )
is a sequence of pickup and drop-off nodes, where the vehicle arrives empty at i 1 and leaves empty at i k and has passenger(s) on board at other nodes. Then, we call F a battery-restricted fragment if there exists a feasible route R of the form :

(o, s i1 , • • • , s iv , F i 1 , i 2 , • • • , i k , s iv+1 , • • • , s im , f ), where s i1 , • • • , s iv , s iv+1 , • • • , s im (v, m ⩾ 0) are recharging stations and o ∈ O, f ∈ F .
In the case of route R is feasible without a recharging station (i.e., v = m = 0), the battery-restricted fragment is equivalent to the one defined in [START_REF] Rist | A new formulation for the dial-a-ride problem[END_REF].

Based on Definition 1, each E-ADARP route can be regarded as the concatenation of batteryrestricted fragments, recharging stations (if required), origin depot, and destination depot. Clearly, this representation guarantees the total excess user ride time for P being non-decreasing. Meanwhile, it eliminates the need to determine the reduced cost at every node along P as we can generalize a sequence of nodes to an arc when extending P by a fragment.

Abstracting fragments as arcs

In the labeling algorithm, a partial path is not fixed until reaching the sink node, which introduces difficulty in ensuring excess-time optimality while maintaining feasibility. This section solves this issue from another perspective : assuming a fragment F = {1, 2, • • • , m} and any excess-time optimal schedule A is presented as a set of service start times over F, i.e.,

A = (A 1 , A 2 , • • • , A m ).
Based on our representation of paths, we only need to determine all possible values of A 1 , A m for any excess-time optimal schedule A, such that we can maintain excess-time optimality by restricting time windows at node 1 and node m. To calculate all possible values of A 1 , A m , we introduce vehicle-waiting-time optimal schedules : Definition 2. A vehicle-waiting-time optimal schedule B for a fragment F is defined as a set of service start times B i , i ∈ F that minimize the total waiting time of the vehicle waits at each node along

F (i.e., m i=2 [B i -(B i-1 + t i-1,i + s i-1 )]
). For a given fragment F, we determine the following two vehicle-waiting-time optimal schedules. We show later in Theorem 1 that these schedules determine all possible values of A 1 , A m . 1. the "latest" vehicle-waiting-time optimal schedule B l ; 2. the "earliest" vehicle-waiting-time optimal schedule B e .

Assuming that B l i is the service start time for schedule B l at node i. Then, the arrival time at each node i is

Arr i = B l i-1 + t i-1,i + s i-1 , 2 ⩽ i ⩽ m.
The waiting time ∆ i at node i is calculated as ∆ i = max{0, B l i -Arr i }. B l i is defined inductively as follows : 1. B l 1 = l 1 ; 2. assuming B l i has been defined for i < v, we define B l v by : (a) if the extension from node (v -1) will not violate the time window constraint at node v (i.e., Arr v ≤ l v ), then we define

B l v = max{e v , Arr v } ; (b) Otherwise, -if min i<v {0, B l i -e i } ≥ Arr v -l v ,
we can update the schedule at nodes 1, • • • , v -1 by moving forward Arr v -l v . Then the extension from node (v -1) to node v will not violate the time window constraint. I.e., we update B l i as

B l i -(Arr v -l v ) for i = 1, • • • , v -1 and define B l v = l v ; -Otherwise,
there is no feasible schedule for {1, 2, • • • , v}, as the time window constraint at node v is always violated.

After determining B l on F, B e is determined by moving forward B l on F by the maximum amount of time that will not change the minimum vehicle waiting time. We denote the maximum amount of time that B l can be moved forward as ← ∆ , which is calculated by taking the minimum value among all the B l i -e i , that is :

← ∆ = min i∈F {B l i -e i }. Then, B e i = B l i - ← ∆ , i ∈ F.
Clearly, for any excess-time optimal schedule A of a fragment F, we can obtain A by moving forward B l by a certain amount of time δ if B e ̸ = B l . In case of waiting time generated, we proved that the service start time of A is the same as B l at the first and last node of the fragment (Theorem 1). The proof is omitted due to the limit of space.

Theorem 1. Assuming that fragment F = {1, 2, • • • , m}, A is an excess-time optimal schedule and B l is the constructed latest vehicle-waiting time optimal schedule over F. Then there exists δ e , δ l ≥ 0 such that :

A 1 = B e 1 + δ e , A m = B e m + δ e ; A 1 = B l 1 -δ l , A m = B l m -δ l .
Based on Theorem restricting time windows at node 1 and node m to [B e 1 , B l 1 ] and [B e m , B l m ] will include all excess-time optimal schedules on fragment F = {1, • • • , m}. Then we can abstract F to an arc (1, m) such that :

1. the total travel time from 1 to m (denoted as

t ′ 1,m ) is B l m -B l 1 ; 2. the time window of node 1 is [B e 1 , B l 1 ] and [B e m , B l m ] for node m ; 3. the energy consumption from 1 to m is m-1 i=1 h i,i+1 ;

Construction of a new sparser graph

In this section, we construct a new sparser graph G sp by enumerating all feasible fragments and abstracting them to arcs. To generate all feasible fragments, we start from each pickup node and extend it in a feasible way, assuming that the initial battery inventory is Q. Then, these "arcs" are connected to depots and recharging stations in a feasible way. Figure 1 shows an example of constructing new arcs in G sp . It should be noted that for two different fragments, even though they have the same start node i+ and end node j-, we must treat them as two different arcs in G sp as they represent different fragments containing different sequences of nodes and have different restricted time windows.

Forward labeling algorithm

We develop a labeling algorithm on G sp by extending the one proposed by [START_REF] Desaulniers | Exact algorithms for electric vehicle-routing problems with time windows[END_REF] considering features of the E-ADARP. We denote L i as the label associated with a partial path ends with node i. The forward labeling algorithm extends labels from a source node o k ∈ O to a sink node f ∈ F . Let a label associated with a partial path P from o k to current vertex i be

L i = {T cost i , (T rchs i ) s∈S , T tM in i , T tM ax i , T rtM ax i , T req i }
, the definition of each resource is described as follows :

1. T cost i : The reduced cost of the partial route until i ; 2. T rchs i : The number of times recharging station s ∈ S is visited along partial path P ; 3. T tM in i : The earliest service start time at vertex i assuming that, if a recharging station is visited prior to i along the partial path, a minimum recharge (ensuring the battery feasibility up to i) is performed ; 4. T tM ax i : The earliest service start time at vertex i assuming that, if a recharging station is visited prior to i along the partial path, a maximum recharge (ensuring the time-window feasibility up to i) is performed ; 5. T rtM ax i : In order to propagate the information along the path, we make the artificial assumption that vehicles can be recharged at all vertices. But in reality, the vehicle will never go to the recharging station when passengers are on board. With this assumption, T rtM ax i denotes the maximum recharging time required to fully recharge at vertex i if a recharging station is visited prior to i along P, a minimum recharge (ensuring the battery feasibility up to i) is performed ; 6. T req i : The set of visited and unreachable requests until i along partial path P. A request is said to be unreachable if time window constraints are violated if visited. We extend a label

L i = {T cost i , (T rchs j ) s∈S , T tM in i , T tM ax i , T rtM ax i
, T req i } along arc (i, j) of G sp using Resource Extension Functions (REFs) in [START_REF] Desaulniers | Exact algorithms for electric vehicle-routing problems with time windows[END_REF]. The label feasibility is checked via the following proposition : Proposition 1. The route R is feasible if and only of ∀j ∈ R, the label L j satisfies :

T tM in j ⩽ B l j , T tM in j ⩽ T tM ax j , T reqp j ⩽ 1, ∀p ∈ P, T rtM ax j ⩽ (1 -γ)H, j ∈ F H, otherwise
If j is a recharging station, then T rchj j ⩽ 1 must be satisfied. In case of a feasibility violation, the corresponding label will be discarded. The following dominance rules are applied for two labels associated with the same node.

Definition 3. Let L k = {T cost i , (T rchs j ) s∈S , T tM in i , T tM ax i , T rtM ax i , T req i }, k ∈ 1, 2
be two labels. Assume that the partial path associated to L 1 and L 2 are P 1 and P 2 , respectively, and P 1 , P 2 end at the same node. L 1 dominates L 2 if and only if :

T r 1 ⩽ T r 2 , ∀r ∈ {cost, rch, tM in} (2)

Results and Discussion

The proposed labeling algorithm has been embedded in the CG approach to solving the E-ADARP. The obtained CG results are compared to the best-reported results in [START_REF] Bongiovanni | The electric autonomous dial-a-ride problem[END_REF] and [START_REF] Su | A Deterministic Annealing Local Search for the Electric Autonomous Dial-a-Ride Problem[END_REF]. We summarize all results in Table 1, where #opt, #bestlb, and #bestub denote the number of times the respective algorithm provides optimal solutions, the best lower bound/integer solution of all considered algorithms, respectively. B/F denotes the ratio between the number of best solutions obtained and the number of feasible solutions obtained. LB% is the gap between the obtained lower bound and the best-obtained integer solution among all considered algorithms.

On small-to-medium-sized instances (type-a and -u instances), Our CG algorithm optimally solves 50 out of 84 instances and generates high-quality lower bounds with a 0.31% average deviation to the best-known solutions. We obtain 40 equal lower bounds and 24 better lower bounds than those reported in [START_REF] Bongiovanni | The electric autonomous dial-a-ride problem[END_REF]. The overall quality of the lower bound is improved by 1.33% on average. In addition, we provide 10 new optimal solutions on previously solved and unsolved instances. As [START_REF] Bongiovanni | The electric autonomous dial-a-ride problem[END_REF] does not test type-r instances (large-scale instances with up to 8 vehicles and 96 requests), we compare CG results to the best-reported heuristic results of [START_REF] Su | A Deterministic Annealing Local Search for the Electric Autonomous Dial-a-Ride Problem[END_REF]. We report 14 better solutions (5 of them are optimal) as well as 17 lower bounds. Hence, we prove the performance of our labeling algorithm.

type-a and -u instances results summary CG with labeling algorithm

Best-reported results [ Our labeling algorithm offers new insights into designing an exact algorithm for solving a practical version of the electric DARP. One important "by-product" of our labeling algorithm is an exact and efficient schedule optimization method that can determine the excess-time optimal schedule for a given E-ADARP route. This is the first time that excess user ride time minimization is handled exactly in the E-ADARP. This schedule optimization method can also be applied to DARPs with multiple objectives, where total excess user ride time is minimized.

FIG. 1 -

 1 FIG. 1 -Example of constructing arcs in the new graph G sp

  1] scenario #opt #bestlb LB% B/F avg.time(s) #opt #bestlb LB% B/F avg.time(s) #bestlb #bestub B/F avg.time(s) #opt #bestlb #bestub B/F avg.time(s)

	γ = 0.1 21	23	0.10% 25/28 1154.12	23	23	1.17% 26/28 1502.83
	γ = 0.4 20	22	0.13% 26/28 1382.65	19	19	1.84% 21/27 1701.36
	γ = 0.7 9	20	0.71% 20/26 2454.04	11	15	1.93% 14/19 3965.85
	Overall 50	64	0.31% 71/82 1663.60	53	57	1.64% 61/74 2390.01
			type-r instances results summary
		CG with labeling algorithm		Best-reported results [5]
	scenario #opt γ = 0.1 3	10	8	8/10 9423.78	0	NA	2	2/10	981.36
	γ = 0.4 2	7	6	6/8	9974.08	0	NA	3	3/9	1368.13
	Overall 5	17	14	14/18 9698.93	0	NA	5	5/19 1174.74

TAB. 1 -Results summary of CG algorithm on all benchmark instances


Simone Di Marino
email: simone.dimarino@unige.it

Augusto Gerolin
email: agerolin@uottawa.ca

Luca Nenna
email: luca.nenna@universite-paris-saclay.fr

UNIVERSAL DIAGONAL ESTIMATES FOR MINIMIZERS OF THE LEVY-LIEB FUNCTIONAL

estimate of the probability of the particles being in positions (x 1 , . . .

Introduction

Density Functional Theory attempts to describe all the relevant information about a many-body quantum system at ground state in terms of the one electron density ρ. Following the Levy and Lieb's approach [START_REF] Levy | Universal variational functionals of electron densities, first-order density matrices, and natural spin-orbitals and solution of the v-representability problem[END_REF][START_REF] Elliott | Density functionals for Coulomb systems[END_REF] the ground state energy can be rephrased as the following variational principle involving only the electron density

E 0 [v] = inf √ ρ∈H 1 (R 3 ) ρ(R 3 )=N R 3 v(x) dρ<+∞ F LL,ε [ρ] + R 3 v(x) dρ ,
where v is an external potential and the Levy-Lieb functional F LL,ε is defined as

F LL,ε [ρ] := min ψ∈H N ψ →ρ R 3N ε|∇ψ| 2 (x) + v ee (x)|ψ| 2 (x) dx , (1)
where H N = N i=1 H 1 (R 3 ; C) is the fermionic space, v ee (x 1 , . . . , x N ) = i<j

1 |x i -x j |
is the Coulomb interaction potential between the N electrons and ψ → ρ means that the one-body density of ψ is ρ, that is ρ = N R 3(N -1) |ψ| 2 .

The Levy-Lieb functional is indeed the lowest possible (kinetic plus interaction) energy of a quantum system having the prescribed density ρ. This universal functional is the central object of Density Functional Theory, since knowing it would allow one to compute the ground state energy of a system with any external potential v. For a complete review on it we refer the reader to [START_REF] Mathieu Lewin | Universal functionals in density functional theory[END_REF].

Although the electrons are fermions, in this article we will treat only the case of bosonic wave-functions, i.e. consider in the constraint search in [START_REF] Benamou | A numerical method to solve multi-marginal optimal transport problems with coulomb cost[END_REF] wavefunctions ψ that are symmetric rather than anti-symmetric. Notice that the ground-state energy of bosonic systems are generally higher than fermionic ones. In our analysis, however, the bosonic case is not very restrictive since we are looking at the regime ε small.

Our approach interprets the Levy-Lieb functional as a (Fisher-information regularized) multi-marginal optimal transport problem.

Connection with Optimal Transportation Theory: It has been recently shown that the limit functional as ε → 0 corresponds to a multi-marginal optimal transport problem [START_REF] Bindini | Optimal transport with Coulomb cost and the semiclassical limit of density functional theory[END_REF][START_REF] Cotar | Density functional theory and optimal transportation with Coulomb cost[END_REF][START_REF] Cotar | Smoothing of transport plans with fixed marginals and rigorous semiclassical limit of the Hohenberg-Kohn functional[END_REF][START_REF] Mathieu Lewin | Semi-classical limit of the Levy-Lieb functional in Density Functional Theory[END_REF] (see also the seminal works in the physics and chemistry literature [START_REF] Buttazzo | Optimal-transport formulation of electronic density-functional theory[END_REF][START_REF] Seidl | The strictly-correlated electron functional for spherically symmetric systems revisited[END_REF][START_REF] Seidl | Strong-interaction limit of density-functional theory[END_REF][START_REF] Seidl | Strictly correlated electrons in density-functional theory: A general formulation with applications to spherical densities[END_REF][START_REF] Seidl | Strictly correlated electrons in densityfunctional theory[END_REF]): rather than wave-functions, one has now enlarged the constrained search in [START_REF] Benamou | A numerical method to solve multi-marginal optimal transport problems with coulomb cost[END_REF] to minimize among probability measures on R 3N having ρ as marginal, that is

F 0 [ρ] := inf P∈Π N (ρ) R 3N v ee (x 1 , . . . , x N ) dP(x 1 , . . . , x N ) , (2)
where Π N (ρ) denotes the set of probability measures on R 3N having ρ/N as marginals.

The multi-marginal optimal transport with Coulomb cost (2) has garnered attention in the mathematics, physics and chemistry communities and the literature on the subject is growing considerably. Recent developments include results on the existence and non-existence of Monge-type solutions minimizing (2) (e.g., [START_REF] Bindini | On Seidl-type maps for multi-marginal optimal transport with Coulomb cost[END_REF][START_REF] Buttazzo | Optimal-transport formulation of electronic density-functional theory[END_REF][START_REF] Colombo | Multimarginal optimal transport maps for one-dimensional repulsive costs[END_REF][START_REF] Colombo | Equality between Monge and Kantorovich multimarginal problems with Coulomb cost[END_REF][START_REF] Colombo | Counterexamples in multimarginal optimal transport with coulomb cost and spherically symmetric data[END_REF][START_REF] Cotar | Density functional theory and optimal transportation with Coulomb cost[END_REF]18,[START_REF] Friesecke | A simple counterexample to the Monge ansatz in multimarginal optimal transport, convex geometry of the set of Kantorovich plans, and the Frenkel-Kontorova model[END_REF]), structural properties of Kantorovich potentials (e.g., [START_REF] Buttazzo | Continuity and estimates for multimarginal optimal transportation problems with singular costs[END_REF][START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF][START_REF] Di Marino | Optimal Transportation Theory with Repulsive Costs[END_REF][START_REF] Gerolin | Duality theory for multi-marginal optimal transport with repulsive costs in metric spaces[END_REF]), grand-canonical optimal transport [START_REF] Di Marino | Grand-canonical optimal transport[END_REF], efficient computational algorithms (e.g., [START_REF] Benamou | A numerical method to solve multi-marginal optimal transport problems with coulomb cost[END_REF][START_REF] Coyaud | Approximation of optimal transport problems with marginal moments constraints[END_REF][START_REF] Friesecke | Genetic column generation: Fast computation of high-dimensional multi-marginal optimal transport problems[END_REF][START_REF] Khoo | Semidefinite relaxation of multimarginal optimal transport for strictly correlated electrons in second quantization[END_REF][START_REF] Di | An optimal transport approach for the schrödinger bridge problem and convergence of sinkhorn algorithm[END_REF]) and the design of new density functionals (e.g., [START_REF] Chen | Pair densities in density functional theory[END_REF][START_REF] Gerolin | Kinetic correlation functionals from the entropic regularisation of the strictly-correlated electrons problem[END_REF][START_REF] Lani | The adiabatic strictly-correlated-electrons functional: kernel and exact properties[END_REF][START_REF] Mirtschink | Energy density functionals from the strong-coupling limit applied to the anions of the he isoelectronic series[END_REF]). The first order expansion around the limit ε → 0 of the Levy-Lieb functional was obtained in [START_REF] Colombo | First order expansion in the semiclassical limit of the Levy-Lieb functional[END_REF].

We refer to the surveys (and references therein) [START_REF] Di Marino | Optimal Transportation Theory with Repulsive Costs[END_REF][START_REF] Friesecke | The strong-interaction limit of density functional theory[END_REF] for a self-contained presentation on multi-marginal optimal transport approach in Density Functional Theory as well as the review article [START_REF] Vuckovic | Density functionals based on the mathematical structure of the stronginteraction limit of dft[END_REF] for a the recent developments from a chemistry standpoint.

Main result of this paper:

In [START_REF] Buttazzo | Continuity and estimates for multimarginal optimal transportation problems with singular costs[END_REF][START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF][START_REF] De | Optimal transport with Coulomb cost. Approximation and duality[END_REF] it is shown that the support supp(P * ) of a solution P * of the limiting problem (2) is uniformly bounded away of the diagonal, i.e. one has always |x i -x j | ≥ δ > 0 for any x i , x j ∈ supp(P * ). In other words, the electrons are always at a certain distance away from each other, which is the expected behaviour since we are in a classical framework.

In the sequel we will denote with D δ the enlarged diagonal

D δ = {(x 1 , . . . , x N ) ∈ R 3N : ∃ i = j s.t. |x i -x j | ≤ δ}.
In particular the result in [START_REF] Buttazzo | Continuity and estimates for multimarginal optimal transportation problems with singular costs[END_REF][START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF] can be rephrased saying that the solution to the multi-marginal optimal transport problem is concentrated on the complement of D δ for some δ. An important feature of the results is that δ depends only on concentration properties of ρ. In fact defining κ(ρ, r) := sup x∈R 3 ρ(B(x, r))/N, the authors in [START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF] prove that if κ(ρ, β) < 1 2(N -1) then one can choose δ = β 2N . Our main result is to extend this property also for ε > 0 small. In particular we do not expect to have ψ ε = 0 on D δ but we show that the probability of having the electrons in position x ∈ D δ is very small (3).

Theorem 1.1 (Exponential off-diagonal localization for Coulomb). Let ψ ε be a minimizer for (1) where v ee (x 1 , . . . , x N ) = i<j

1 |x i -x j | . Let us consider β such that κ(ρ, β) ≤ 1 4(N -1) then, let α ≤ β
32N , and suppose εN 2 α/16. Then, for

P ε (x) = |ψ ε | 2 (x) we have Dα P ε (x) dx ≤ e -1 24 √ α ε . (3)
In the proof we actually work with a general repulsive pairwise potential v ee , which satisfy the hypothesis [START_REF] Buttazzo | Optimal-transport formulation of electronic density-functional theory[END_REF], stated in the next section. The result in general is the following one: Theorem 1.2 (Exponential off-diagonal localization for general interaction cost). Let ψ ε be a minimizer for [START_REF] Benamou | A numerical method to solve multi-marginal optimal transport problems with coulomb cost[END_REF] where v ee satisfies [START_REF] Buttazzo | Optimal-transport formulation of electronic density-functional theory[END_REF]. Let us consider β such that κ(ρ, β) ≤ 1 4(N -1) then, let α such that m(2α) ≤ 8(N -1)M (β/2), and suppose εN 2 α 2 m(2α). Then, for

P ε (x) = |ψ ε | 2 (x) we have Dα P ε (x) dx ≤ e -1 6 α 2 m(2α) 8ε . (4)
Notice that in [START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF] the diagonal estimate is proven also in the weaker (and sharper) hypotesis κ(ρ, β) < 1 N : while we believe that also in that case a similar generalization in the case ε > 0 holds true, the proof will be more technical and not so trasparent. For the same reason the inequality κ(ρ, β) ≤ 1 4(N -1) is used instead of κ(ρ, β) < 1 2(N -1) in order to have more transparent estimates in the end.

Remark 1.1 (Natural assumption on ε). We want here to justify the hypotesis εN 2 α/16. Notice that we expect |ψ ε | = 0 on the diagonal since v ee = +∞ there; if |ψ ε (x)| ∼ 1 on a point in D α then we then expect |∇ψ ε | ∼ α -1 and so the kinetic energy locally is ε α 2 while the potential energy is 1 α . Since the localization phenomenon outside the enlarged diagonal is a classical one, we expect that it holds in the regime where the kinetic energy is negligible with respect to the potential energy, and that is precisely when ε α. Then there is a factor N which is constant, and we believe it is actually technical.

Organization of the paper: In Section 2 we introduce the notations we are going to use throughout all the paper. In Section 3 we give some estimates concerning kinetic energy term in the Levy-Lieb functional. Section 4 is then devoted to the construction of a competitor for the Levy-Lieb functional; finally in Section 5 we derive the diagonal estimates for the wave-function and, thus, prove Theorem 1.1 and Theorem 1.2 via the iteration of a decay estimate.

Notation

Consider a subset I ⊆ {1, . . . , N }, with cardinality k = |I|, defined as

I = {i 1 , . . . , i k }, with 1 ≤ i 1 < i 2 < • • • < i k . Then, the I-projection π I : R 3N → R 3k , π I ((x 1 , . . . , x N )) = (x i 1 , x i 2 , . . . , x i k ).
Sometimes we will denote x I = π I (x) and if I = J c , then x = (x I , x J ). With a slight abuse of notation, for a function P ∈ L 1 (R 3N ), I ⊆ {1, . . . , N } and J = I c we denote

(π I ) (P)(x I ) = P(x I , x J ) dx J ,
which on density of measures act precisely as the push-forward through the projection function π I .

As we have already mentioned above, we denote by Π N (ρ) the set of probability measures on R 3N having the N one body marginals equal to ρ.

In the following we will consider an electron-electron pair interaction repulsion potential, that v ee with the following form:

v ee (x 1 , . . . , x N ) = i<j c(x i , x j ), where m(|x -y|) ≤ c(x, y) ≤ M (|x -y|) ∀x, y ∈ R 3 for some m, M : (0, ∞) → [0, ∞) decreasing such that lim t→0 + m(t) = +∞ (5)
Moreover, with a slight abuse of notation, we will denote by

P ∈ P(R 3N ) → v ee (P) := R 3N v ee (x 1 , . . . , x N ) dP(x 1 , . . . , x N ) (6)
Notice that we will often identify a measure P with its density. Finally, given an open set Ω ⊆ R 3N , then for every r > 0 we denote

Ω -r := {x ∈ R 3N : B(x, r) ⊆ Ω}. (7)

Estimate for the kinetic energy

In this section we give some preliminary estimates for the kinetic energy term of the Levy-Lieb functional. Denoting L 1 + the cone of positive L 1 functions, we define E kin : L 1 + (R 3N ) → R the Kinetic energy associated to some absolutely continuos N -probability measure h

E kin (h) := R 3N N i=1 |∇ i h(x 1 ,...,x N )| 2 h(x 1 ,...,x N ) dx 1 , . . . , dx N if √ h ∈ H 1 (R 3N ) +∞ otherwise (8)
When it will be clear from the context we will also abbreviate E kin (h) = |∇h| 2 h dx. Notice that the kinetic energy functional is also know as the Fisher information.

Moreover if ψ ∈ H 1 (R 3N ; R), then |∇ψ| 2 dx = E kin (|ψ| 2 ) = E kin (P ψ ),
where P ψ = |ψ| 2 is the joint probability associated to the wave-function ψ. The string of equalities above is thus true when ψ is a minimizer for the bosonic case. The following Lemma summarises some results concerning the homogeneity, subadditivity (which is a consequence of theorem 7.8 in [START_REF] Lieb | Analysis[END_REF]) and the decomposability under projection of the kinetic energy (a similar result also appears in [START_REF] Michael | The Hartree limit of Born's ensemble for the ground state of a bosonic atom or ion[END_REF][START_REF] Rougerie | Théorèmes de de Finetti, limites de champ moyen et condensation de Bose-Einstein[END_REF]). Lemma 3.1. Let E kin defined as in [START_REF] Colombo | Equality between Monge and Kantorovich multimarginal problems with Coulomb cost[END_REF]. Then (i) E kin is 1-homogeneous, that is E kin (λP) = λE kin (P) for every λ > 0;

(ii) given P 1 , . . . , P k ∈ L 1 (R 3N ), we have E kin (P 1 + . . . + P k ) ≤ E kin (P 1 ) + E kin (P 2 ) + . . . + E kin (P k );

(iii) Let P ∈ L 1 + (R 3N
). Given I, J ⊆ {1, . . . , N } two nonempty disjoint sets such that I = J c , we denote by P I = (π I ) P and P J = (π J ) P. Then we have (here N I = I and N J = J)

E N kin (P) ≥ E N I kin (P I ) + E N J kin (P J )
, where the equality holds if and only if P(x) = P I (x I )P J (x J )/λ, where λ = P. In particular if P is the density of a probability measure, we have that the equality happens if and only if x I and x J are independent under the probability P.

Proof. (i) The 1-homogeneity is obvious. (ii) For the subadditivity it is sufficient to prove it for k = 2; then for every x, by Cauchy-Schwarz inequality we have

(P 1 (x) + P 2 (x)) |∇P 1 (x)| 2 P 1 (x) + |∇P 2 (x)| 2 P 2 (x) ≥ (|∇P 1 (x)| + |∇P 2 (x)|) 2 ,
which, after using the triangular inequality and dividing by P 1 + P 2 can be rewritten as

|∇(P 1 + P 2 )| 2 P 1 + P 2 ≤ |∇P 1 | 2 P 1 + |∇P 2 | 2 P 2
, which integrated gives us the conclusion.

(iii) As for the last point we fix x J and we use the Cauchy-Schwarz inequality with respect to the measure dx I :

P(x I , x J )dx I • |∇ J P(x I , x J )| 2 P(x I , x J ) dx I ≥ |∇ J P(x I , x J )|dx I 2 ≥ ∇ J P(x I , x J )dx I 2 ,
where in the last passage we used the triangular inequality and we took the derivative out of the integral. Now we recognize P J (x J ) = P(x I , x J )dx I and so we can write this as

|∇ J P(x I , x J )| 2 P(x I , x J ) dx I ≥ |∇ J P J (x J )| 2 P J (x J )
.

Integrating this with respect to dx J and doing a similar computation for x I , we obtain the conclusion, that is

|∇P(x I , x J )| 2 P(x I , x J ) dx I dx J ≥ |∇ J P J (x J )| 2 P J (x J ) dx J + |∇ I P I (x I )| 2 P I (x I ) dx I .
From the equality cases in C-S and triangular inequality combined we get ∇ J P(x I , x J ) = v(x J )P(x I , x J ) for some vector field v; by a simple integration we actually get v = ∇(P J )/P J ; this can be seen as ∇ J log(P) = ∇ J log P J ; similarly we can get ∇ I log(P) = ∇ I log P I . Summing up this two equalities we get ∇(P(x)/P I (x I )P J (x J )) = 0.

The following lemma is a straightforward adaptation of Theorem 3.2 in [START_REF] Cycon | Schrödinger operators with application to quantum mechanics and global geometry[END_REF] giving the IMS localization formula; we have added a short proof for sake of completeness. Lemma 3.2. Let η 1 , η 2 , η 3 : R 3N → [0, 1] be C 1 functions such that η 1 +η 2 +η 3 ≡ 1. Then, for every function P ∈ L 1 + (R 3N ) we have

E kin (Pη 1 )+E kin (Pη 2 )+E kin (Pη 3 ) = E kin (P)+ |∇η 1 | 2 η 1 + |∇η 2 | 2 η 2 + |∇η 3 | 2 η 3 P dx.
Proof. For every i = 1, 2, 3 pointwisely we have:

|∇(Pη i )| 2 Pη i = |η i ∇P + P∇η i | 2 Pη i = η 2 i |∇P| 2 + 2η i P∇P • ∇η i + P 2 |∇η i | 2 Pη i = η i |∇P| 2 P + 2∇P • ∇η i + P |∇η i | 2
η i Adding them up and using that η i = 1 and

∇η i = 0, we get i |∇(Pη i )| 2 Pη i = |∇P| 2 P + P i |∇η i | 2 η i ,
which integrated, gives us the desired identity.

New trial state: swapping particles and estimate for the potential

The scope of this subsection is to create a competitor for the minimization of the functional

F LL,ε (P) = εE kin (P) + v ee (P) if P ∈ Π N (ρ) +∞ otherwise,, (9)
where E kin is defined in ( 8) and v ee satisfies [START_REF] Buttazzo | Optimal-transport formulation of electronic density-functional theory[END_REF]. The idea is to try to mimic what it is done in [START_REF] Buttazzo | Continuity and estimates for multimarginal optimal transportation problems with singular costs[END_REF][START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF][START_REF] De | Optimal transport with Coulomb cost. Approximation and duality[END_REF], in the semiclassical case ε = 0: in that case we take two points y, z ∈ R 3N and substitute them with ỹ, z where we have interchanged their first compenent, that is ỹ = (z 1 , y 2 , . . . , y n ) and z = (y 1 , z 2 , . . . , z n ).

In order to do so for the n-particle distribution P, we will consider two small bumps centered around y and z

η 1 (x) = λ 1 η x -y r 1 and η 2 (x) = λ 2 η x -z r 2 , (10)
for some λ 1 , λ 2 , r 1 , r 2 to be chosen later and some η ∈ C 1 c (B(0, 1)), η ≥ 0. First of all we assume that supp(η 1 ) ∩ supp(η 2 ) = ∅, which can be granted as long as

r 1 + r 2 < |y -z|, (11)
and then we assume η 1 P = η 2 P = m which can be accomplished again by choosing the appropriate λ i , r i . Let us then define

ρ i 1 (x 1 ) = (π {1} ) (η i P), ρ i 1(x 2 , x 3 , . . . , x N ) = (π {1} c ) (η i P), (12)
P 1 = 1 m ρ 2 1 ρ 1 1, P 2 = 1 m ρ 1 1 ρ 2 1, (13)
where ρ i 1 and ρ i 1 are the marginals of η i P and P 1 and P 2 are densities concentrated around ỹ = (z 1 , y 2 , . . . , y n ) and z = (y 1 , z 2 , . . . , z n ) respectively. We then finally consider

P := P -Pη 1 -Pη 2 + P 1 + P 2 , ( 14
)
which will be the competitor for a minimizer P of the functional F LL,ε .

Remark 4.1. In order to have that P is a competitor, we still have to check that P ≥ 0 and that it has the correct marginals. For the positivity, notice that if λ 1 and λ 2 are small enough we have also η 1 +η 2 ≤ 1 and in particular P-η 1 P 1 -η 2 P 2 , which will guarantee that P ≥ 0.

For the marginal constraint, notice that by [START_REF] Cotar | Density functional theory and optimal transportation with Coulomb cost[END_REF] and (13) we have that η 1 P+η 2 P and P 1 + P 2 have the same marginals, in particular also P and P share the same marginals.

Lemma 4.1. Let P be such F LL,ε (P) < +∞. Given y, z ∈ R 3N , let η 1 , η 2 , P 1 , P 2 , P defined by [START_REF] Colombo | First order expansion in the semiclassical limit of the Levy-Lieb functional[END_REF],( 11), ( 12), ( 13) and [START_REF] Coyaud | Approximation of optimal transport problems with marginal moments constraints[END_REF]. Then E kin ( P) ≤ E kin (P) + P(x)

|∇η 1 | 2 η 1 + |∇η 1 | 2 η 2 + |∇η 1 + ∇η 2 | 2 1 -η 1 -η 2 dx; v ee ( P) = v ee (P) -P(η 1 + η 2 ) i>1 c(x 1 , x i )dx + (P 1 + P 2 ) i>1 c(x 1 , x i )dx.
Proof. Let us consider η 3 = 1 -η 1 -η 2 . Then we have P = η 3 P + P 1 + P 2 . Using Lemma 3.1, in particular the subadditivity and the exact energy split in case of independent variables for E kin , we get (by ( 13))

E kin ( P) ≤ E kin (η 3 P) + E kin (P 1 ) + E kin (P 2 )

= E kin (η 3 P) + E kin (ρ 2 1 ) + E kin (ρ 1 1) + E kin (ρ 1 1 ) + E kin (ρ 2 1); (15)
we then recall [START_REF] Cotar | Density functional theory and optimal transportation with Coulomb cost[END_REF] and the inequality for the split energy (Lemma 3.1 (iii)) to get

E kin (ρ i 1 ) + E kin (ρ i 1) ≤ E kin (η i P) (16)
and so we conclude using ( 15), ( 16) and then Lemma 3.2.

For the estimate with the potential, it is clear that v ee ( P) = v ee (P) -P(η 1 + η 2 )v ee (x)dx + (P 1 + P 2 )v ee (x)dx;

Since v ee (x) = i<j c(x i , x j ) we just need to show that the contribution due to c(x i , x j ) whenever 1 < i < j cancels out in the last two integrals. In fact in both integrals we can integrate out the first variable: denoting I = {1} and J = I c for example we have

Pη 1 c(x i , x j ) dx I dx J = c(x i , x j ) Pη 1 dx I dx J = c(x i , x j )ρ 1 1(x J ) dx J = c(x i , x j )ρ 1 1(x J ) ρ 1 2 (x I ) m dx I dx J = c(x i , x j )P 1 dx.
In a similar way we can show that Pη 2 c(x i , x j ) dx = P 2 c(x i , x j ) dx.

In the sequel we will denote

C 1 (x) = N i=2 c(x 1 , x i ) 5.

Diagonal estimates for the wave-function

We devote this last section to derive the diagonal estimates for the bosonic wavefunction which minimizes the Levy-Lieb functional proving in particular Theorem 1.1 and Theorem 1.2. Lemma 5.1. Let ρ be a marginal distribution and let β > 0 be such that κ(ρ, β) ≤ 1 4(N -1) . Then, for every P ∈ Π N (ρ) and y ∈ R 3N , for every r 1 , r 2 such that r 1 + 2r 2 < β and δ > 0, there exists z ∈ R 3N such that, defining η 1 , η 2 , P 1 , P 2 , m as in [START_REF] Colombo | First order expansion in the semiclassical limit of the Levy-Lieb functional[END_REF], ( 12) and (13)

(i) C 1 (P 1 + P 2 )dx ≤ 2(N -1)M (β -r 1 -2r 2 )m;
(ii) z is a (1 + δ, 1/2)-doubling point at scale r 2 for P, that is

B(z,r 2 ) P dx ≤ 2(1 + δ) 3N B(z,r 2 /(1+δ)) P dx.
Proof. For γ > 0, let us consider the set

Ω = y ∈ R 3N : |y 1 -y i | ≥ γ - δ 1 + δ r 2 and |y 1 -y i | ≥ γ - δ 1 + δ r 2 , ∀i = 2, . . . , N .
We know that if z ∈ Ω we will have of course

C 1 (y 1 , z 2 , . . . , z N )+C 1 (z 1 , y 2 , . . . , y N ) ≤ 2(N -1)M (γ-r 1 -2r 2 ) ∀y ∈ B(y, r 1 ), z ∈ B(z, r 2 ),
which in particular implies C 1 (P 1 + P 2 ) dx ≤ 2(N -1)M (γ -r 1 -2r 2 )m. Now we want to see that there exists a 1/2 doubling point in Ω; in order to do that, it is easy to see that

Ω -δ r 2 1+δ ⊆ {y ∈ R 3N : |y 1 -y i | ≥ γ and |y 1 -y i | ≥ γ, ∀i = 2, .
. . , N } And now a similar computation to what is done in [START_REF] Buttazzo | Continuity and estimates for multimarginal optimal transportation problems with singular costs[END_REF][START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF] will give us

Ω -δ r 2 1+δ P(x) dx ≥ 1 -2(N -1)κ(ρ, γ),
where κ(ρ, r) = sup x∈R 3 ρ(B(x, r)). Now if we consider γ = β we have κ(ρ, β) ≤ 1 4(N -1) , and so we can apply Lemma 5.2 with r = r 2 1+δ get the existence of a (1 + δ, 1/2)-doubling point at scale r 2 in Ω. Lemma 5.2 (Existence of doubling points). Let P ∈ L 1 + (R 3N ) be the density of a probability measure and let δ > 0. Let us consider an open set Ω ⊆ R 3N ; then for every r > 0 we denote C r := Ω -r P(x) dx, where Ω -r is defined as in [START_REF] Colombo | Multimarginal optimal transport maps for one-dimensional repulsive costs[END_REF]. Then, whenever C δr > 0, there exists y ∈ Ω, such that B(y,(1+δ)r)

P(x) dx ≤ (1 + δ) 3N
C δr B(y,r)

P(x) dx,
that is, the measure P(x) dx is doubling at the point y at scale r, with doubling constant (1+δ) 3N C δr . Proof. Suppose on the contrary that for every y ∈ Ω the reversed inequality holds . Then there exists α 0 = α(β, ε) such that if P minimizes (9), we have that for every y ∈ D α such that α ≤ α 0 , and every r 1 ≤ α/2, we have

B(y,r 1 /(1+δ)) P(x) dx ≤ 1 δ 2 r 2 1 m(2α) 2(1+δ) 2 ε + 1 B(y,r 1 ) P(x) dx, (17)
whenever δ > 0 is such that m(2α) > 256εC(δ)/β 2 , where

C(δ) := (1 + δ) 2 • (2(1 + δ) 3N -1) δ 2 . ( 18
)
An implicit choice for α 0 is for example m(2α 0 ) > 8 max{(N -1)M (β/2), 850εN 2 /β 2 }.

Proof. Let y ∈ D α and without loss of generality we can assume that |y 1 -y 2 | < α; let z given by Lemma 5.1. We then consider r 1 , r 2 , η 1 , η 2 , P 1 , P 2 , P defined by ( 10),( 11), ( 12), ( 13) and ( 14); being P ∈ Π(ρ), we get, by the minimality of P, F LL,ε ( P) ≥ F LL,ε (P), εE kin ( P) + v ee ( P) ≥ εE kin (P) + v ee (P); now we can use the estimates in Lemma 4.1 in order to conclude that

ε P |∇η 1 | 2 η 1 + |∇η 2 | 2 η 2 + |∇η 1 + ∇η 2 | 2 1 -η 1 -η 2 dx ≥ P(η 1 +η 2 )C 1 dx-(P 1 +P 2 )C 1 dx. Now we make the choice η(x) = min (1+δ)(1-|x|) + δ , 1 2
. In particular 0 ≤ η ≤ 1 and η ≡ 1 if |x| < 1 1+δ , and moreover |∇ √ η| ≡ 1+δ δ if 1 1+δ ≤ |x| ≤ 1 and 0 otherwise. We thus have

P |∇η 1 | 2 η 1 = (1 + δ) 2 δ 2 r 2 1 B(y,r 1 )\B(y,r 1 /(1+δ)) Pλ 1 dx = (1 + δ) 2 δ 2 r 2 1 B(y,r 1 ) Pλ 1 dx - B(y,r 1 /(1+δ)) Pλ 1 ≤ (1 + δ) 2 δ 2 r 2 1 B(y,r 1 ) Pλ 1 dx -m
In a similar way we have

P |∇η 2 | 2 η 2 ≤ (1 + δ) 2 δ 2 r 2 2 B(z,r 2 ) Pλ 2 dx -m ≤ (1 + δ) 2 • (2(1 + δ) 3N -1) δ 2 r 2 2 •m = C(δ) r 2 2 •m.
where in the last steps we used Lemma 5.1 (ii) and the definition of C(δ) (18). Now we use that C 1 Pη 1 dx ≥ m(α + 2r 1 ) • m and the estimates we have for C 1 (P 1 + P 2 ) dx to get

P(η 1 + η 2 )C 1 dx -(P 1 + P 2 )C 1 dx ≥ [m(α + 2r 1 ) -2(N -1)M (β -r 1 -2r 2 )] • m.
Putting everything together we have

ε(1 + δ) 2 δ 2 r 2 1 B(y,r 1 ) λ 1 P(x) dx ≥ m(α + 2r 1 ) -2(N -1)M (β -r 1 -2r 2 ) - εC(δ) r 2 2 + ε(1 + δ) 2 δ 2 r 2 1 •m.
We can now use use m ≥ B(y,r 1 /(1+δ)) λ 1 P(x) dx and, dividing λ 1 , we can write the inequality as B(y,r 1 /(1+δ))

P(x) dx ≤ 1 δ 2 r 2 1 F (r 1 ,ε,α) (1+δ) 2 ε
+ 1 B(y,r 1 )

P(x) dx, (19)
where

F
In this way we have F (r 1 , ε, α) ≥ m(2α)/2: plugging this estimate in [START_REF] Di Marino | Grand-canonical optimal transport[END_REF] we get precisely B(y,r 1 /(1+δ))

P(x) dx ≤ 1 δ 2 r 2 1 m(2α)
2(1+δ) 2 ε + 1 B(y,r 1 ) P(x) dx.

In order to understand for which α and δ this inequality holds, we have to ensure that the two conditions [START_REF] Friesecke | A simple counterexample to the Monge ansatz in multimarginal optimal transport, convex geometry of the set of Kantorovich plans, and the Frenkel-Kontorova model[END_REF]

We will now iterate the estimate in Proposition 5.1

Theorem 5.1. Let us consider ρ and β such that κ(ρ, β) < 1 4(N -1) . Then let us consider α < α 0 (as in Proposition 5.1) and suppose A := α 2 m(2α) 8ε N 2 . Then if P minimizes (9) we have that Proof. Let us consider δ such that δ 2 A = e 2 . By the hypothesis on A we have δ 1/N ; in particular, by (18) we can estimate C(δ) ≤ 2 δ 2 , and then it is easy to see that m(2α) > 256εC(δ)/β 2 and thus we can apply Proposition 5.

BProposition 5 . 1 ( 1 4

 511 (y,(1+δ)r)P(x) dx >(1 + δ)3N C δr B(y,r)P(x) dx.Then we can integrate this inequality on the whole Ω Ω B(y,(1+δ)r)P(x) dx dy > (1 + δ) 3N C δr Ω B(y,r) P(x) dx dy. Now we can use Fubini and get ω 3N • (1 + δ)r 3N = R 3N B(y,(1+δ)r) P(x) dx dy ≥ Ω B(y,(1+δ)r) P(x) dx dy C δr ω 3N • r 3N = Ω -δr P(z)|B(z, r)| dz = Ω B(y,r)∩Ω -δr P(x) dx dy ≤ Ω B(y,r) P(x) dx dy, and so we get a contradiction. One step decay). Let us consider ρ and β such that κ(ρ, β) < (N -1)

(r 1 2 :

 12 , ε, α) := max m(α + 2r 1 ) -2(N -1)M (β -r 1 -2r 2 ) -εC(δ) r 2 r 2 > 0 .Now we can choose r 1 ≤ α/2 ≤ β/4 and r 2 = β/8, and then choose α < α 0 such that

D α/ 2 P

 2 (x) dx ≤ e -

1 with r 1 = 1 δ 2 ( 1 + δ) 2k+2 e 2 B

 11212 α k = α 2 (1 + δ) -k to obtain for every y ∈ D α B(y,α k+1 ) P(x) dx ≤ α 2 m(2α) 8(1+δ) 2k+2 ε + 1 B(z,α k ) P(x) dx ≤ (y,α k ) P(x) dx.

  notice that α 0 can be characterized as maximal α for which there exists some δ for which (21) is satisfied that is when C(δ) as small as possible, which is approximately achieved for δ = 2 3N . With this choice we have C(2/(3N )) ≤ 26N 2 and thus m(2α 0 ) ≥ 8 max (N -1)M

	are satisfied, that is
	m(2α) ≥ max 8(N -1)M	β 2	, 256	εC(δ) β 2	;	(21)
		β 2	, 850	εN 2 β 2 .

Acknowledgement. Part of this work has been developed during a Research in Pairs program where the authors were hosted at the MFO (Mathematisches Forschungsinstitut Oberwolfach) in January 2017. S.D.M. is a member of GNAMPA (INdAM). A.G. acknowledges partial support of his research by the Canada Research Chairs Program and Natural Sciences and Engineering Research Council of Canada. L.N. is partially on academic leave at Inria (team Matherials) for the year 2022-2023 and acknowledges the hospitality if this institution during this period. His work was supported by a public grant as part of the Investissement d'avenir project, reference ANR-11-LABX-0056-LMH, LabEx LMH and from H-Code, Université Paris-Saclay. The authors want to thank M. Lewin for comments (and references) on a preliminary draft of the paper.

We can now iterate the estimate for k = 0, . . . , k 0 where (1 + δ) 2k 0 +2 ≤ e 2 ≤ (1 + δ) 2k 0 +4 . At that point we have

Integrating this inequality for y ∈ D α we get

In particular

P(y) dy;

notice that since A N 2 we have

Proof. (Theorem 1.1 and Theorem 1.2) First we notice that if ψ ε is a minimizer for (1) then P ε = |ψ ε | 2 is a minimizer for [START_REF] Colombo | Continuity of multimarginal optimal transport with repulsive cost[END_REF]. Then we notice that if m(2α) ≤ 8(N -1)M (β/2) and εN 2 α 2 m(2α), we have also α < α 0 and so we can apply Theorem 5.1. From that we finish using that P ε is a probability density and so D 2α P ε (y) dy ≤ 1. The conclusions for Theorem 1.1 are then implied by using m(t) = M (t) = 1/t.