System of continuity equations with Newtonian interactions

Measure solutions to a system of continuity equations driven by Newtonian nonlocal interactions


We prove global-in-time existence and uniqueness of measure solutions of a nonlocal interaction system of two species in one spatial dimension. For initial data including atomic parts we provide a notion of gradient-flow solutions in terms of the pseudo-inverses of the corresponding cumulative distribution functions, for which the system can be stated as a gradient flow on the Hilbert space according to the classical theory by Brézis. For absolutely continuous initial data we construct solutions using a minimising movement scheme in the set of probability measures. In addition we show that the scheme preserves finiteness of the -norms for all and of the second moments. We then provide a characterisation of equilibria and prove that they are achieved (up to time subsequences) in the large time asymptotics. We conclude the paper constructing two examples of non-uniqueness of measure solutions emanating from the same (atomic) initial datum, showing that the notion of gradient flow solution is necessary to single out a unique measure solution.

Key words and phrases:
systems of aggregation equations; Newtonian potentials; uniqueness of solutions; gradient flows; long time asymptotics.
1991 Mathematics Subject Classification:
Primary: 45K05, 35A15; Secondary: 92D25, 35L45, 35L80;
Corresponding author: José A. Carrillo

José Antonio Carrillo

Department of Mathematics, Imperial College London

London SW7 2AZ, United Kingdom

Marco Di Francesco

DISIM - Department of Information Engineering, Computer Science and Mathematics, University of L’Aquila

Via Vetoio 1 (Coppito) 67100 L’Aquila (AQ) - Italy

Antonio Esposito

DISIM - Department of Information Engineering, Computer Science and Mathematics, University of L’Aquila

Via Vetoio 1 (Coppito) 67100 L’Aquila (AQ) - Italy

Simone Fagioli

DISIM - Department of Information Engineering, Computer Science and Mathematics, University of L’Aquila

Via Vetoio 1 (Coppito) 67100 L’Aquila (AQ) - Italy

Markus Schmidtchen

Department of Mathematics, Imperial College London

London SW7 2AZ, United Kingdom

(Communicated by the associate editor name)

1. Introduction

In this work we consider a particular instance of the following nonlocal interaction system for the evolution of two probability measures and on the whole real line


Here model the way any two agents of the same species interact with one another (so-called self-interaction potentials, or intraspecific interaction potentials). Respectively, are called the cross-interaction potentials, or interspecific interaction potentials, as they describe the interaction between any two agents of opposing species.

This model can be easily understood as a natural extension of the well-known aggregation equation (cf. [6, 9, 25, 51, 55]) to two species. This link was first established in the paper [36] as the continuous counterpart of a system of ordinary differential equations. More precisely, for suppose and denote the locations of agents of two different species, each of them with masses and respectively. Then, assuming the velocity of any agent is given as an average of the forces exhibited by all other agents upon that agent, one gets

The choice of the interaction potentials depends on the application or the phenomena of interest. In particular in mathematical biology contexts the potentials are often assumed to be radial, i.e. , for , i.e. they only depend on the relative distance between any two agents. An interaction potential is said to be attractive if and repulsive if . The existence theory developed in [36] covers the case of -potentials , , (with suitable growth conditions at infinity) and provides a semigroup defined in the space of probability measures with finite second moment equipped with the Wasserstein distance, in the spirit of [1]. More specifically, the JKO scheme, [43], can be adopted in the special case . Then Eq. (1) can actually be seen as the gradient flow of the functional


In this case a slightly lower regularity needs to be required on the potentials as long as all of them are convex up to a quadratic perturbation. Thus, uniqueness can be proven via the notion of -convexity along geodesics of , see [50, 1].

Common interaction potentials for the one species case include power laws , as for instance in the case of granular media models, cf. [5, 56]. Another possible choice is a combination of power laws of the form , for where is the space dimension. These potentials, featuring short-range repulsion and long-range attraction, are typically chosen in the context of swarming models, cf. [4, 7, 24, 21, 28, 41, 42, 44]. Other typical choices include characteristic functions of sets or Morse potentials

or their regularised versions where and denote the interaction strength and radius of the attractive (resp. repulsive) part and , cf. [28, 29, 38]. These potentials display a decaying interaction strength, e.g. accounting for biological limitations of visual, acoustic or olfactory sense. The asymptotic behaviour of solutions to one single equation where the repulsion is modelled by non-linear diffusion and the attraction by non-local forces has also received lots of attention in terms of qualitative properties, stationary states and metastability, see [17, 19, 20, 22, 40, 23] and the references therein, as well as its two-species counterparts, cf. e.g. [35, 21, 27, 18], and references therein.

As set out earlier we shall study a particular instance of the above system where all interactions are modelled by Newtonian potentials. More precisely, by setting , we consider repulsive Newtonian intraspecific interactions and attractive Newtonian interspecific interactions, i.e. we will deal with the system


Following (2), there is a natural functional that can be associated to system (3), namely


We mention at this stage that this choice of the functional does not fit the set of assumptions in [36], in that the (repulsive) intraspecific parts of are not defined through convex potential (up to a quadratic perturbation).

The corresponding equation for one species has been attracting a lot of interest. In [6] and [9], the authors provide an and an -theory for the aggregation equation , , with initial data in , where and denotes the set of probability measures with bounded second order moments. They consider radially symmetric kernels whose singularity is of order , , at the origin. In particular, the authors prove local well-posedness in for , where . Moreover, when , the exponent is sharp since for any the solution instantaneously concentrates mass at the origin for initial data in . Global well-posedness of solutions with initial data in was proven in [25] for a general class of potentials including in particular . The gradient flow structure was crucial to show a unique continuation after blow-up of solutions to the aggregation equation. Let us also mention [8] that provides a well-posedness theory of compactly supported -solutions for the Newtonian potentials in . The gradient flow structure introduced in [25] in the particular case of in one dimension was further developed in [12], where the authors prove the equivalence of the Wasserstein gradient flow for

with or , and the notion of entropy solution of a scalar nonlinear conservation law of Burgers-type


in the repulsive or attractive case respectively, being . Such a result is relevant in particular in the repulsive case, as it shows that all point particles initial data evolve into an -density on as a simple consequence of the uniqueness of entropy solutions to the corresponding scalar conservation law, [45]. More precisely, a point particle in the repulsive case corresponds to the initial datum for the equation , and the discontinuity can be resolved (in a weak solution sense) either via a stationary Heaviside function or through a rarefaction wave with time-decaying slope connecting the two states and . As the flux is convex, the latter is the only admissible solution in the entropy sense (see e.g. [32]). Therefore, the equivalence result in [12] implies that the distributional derivative is the only gradient flow solution to the repulsive aggregation equation . Notice that such a solution satisfies for all , whereas the initial condition is an atomic measure.

The occurrence of such a smoothing effect in the one-species repulsive case suggests that similar phenomena may be observed in the two-species case, at least in one space dimension. Understanding such an issue is one of the purposes of this work. However, the equivalence to the system of conservation laws


does not provide any useful insights in this case, as we shall discuss in detail in Section 5. We would like to stress at this stage that the persistence of an atomic part for one of the two species in (3) would make the definition of measure solutions rather difficult, as the velocity fields are given by convolutions of the solution with a discontinuous function. In the version (5) this corresponds to the impossibility of e.g. multiplying a discontinuous function by an atomic measure . On the other hand, we shall see that the equivalence to -gradient flows in the pseudo-inverse formalism (see e.g. [13, 30, 12]) provides a natural way to state a suitable notion of solution for (3) with measure initial data giving rise to a unique solution for all times. As for the -regularity of solutions, we mention here that a system similar to (3) with the addition of linear diffusion in both components was studied in the context of semiconductor device modelling, see e.g. [3] in which solutions are shown to maintain the finiteness of the -norms for . As we will show in our paper, the same holds in the one-dimensional diffusion free case (3). However, when initial data feature an atomic part, the attractive part in the functional may inhibit solutions to instantaneously become -densities. This may happen for instance when the two species share an atomic part at the same position initially, see Section 5 below.

Another interesting issue related to (3) is its asymptotic behaviour for large times. The one species case features total collapse, i.e. the formation of one point particle in a finite time in the attractive case with all the mass of the system, see [25], and large time decay to zero in the repulsive case (as a consequence of the results in [12] and on classical results on the large time behaviour for scalar conservation laws, see e.g. [48]). The asymptotic behaviour in the case (3) is much more complex. Finite time concentration for smooth initial data cannot happen in the view of the non-expansiveness of the -norms proven in the present work. Similarly the case for the large time decay to zero is impossible, as we will construct explicit solutions featuring a steady atomic part for all times. We shall prove that the -limit in a suitable topology for (3) is a subset of , which also coincides with the minimising set of the corresponding functional . The rest of this paper is organised as follows:

  • Section 2 contains preliminary concepts on gradient flows in Wasserstein spaces and about the one-dimensional case in particular.

  • Section 3 deals with the existence and uniqueness of solutions. We first prove it in Subsection 3.1, for the notion of solution introduced in Definition 7. In Subsection 3.2 we consider the case of densities as initial conditions, more precisely in for some , and we show that a suitable notion of gradient flow solution in the Wasserstein sense (see Definition 9) can be achieved via the Jordan-Kinderlehrer-Otto scheme, which also allows to prove that the -regularity is maintained. In addition a uniform-in-time control of the second moment is obtained. Moreover, we prove that our solutions also satisfy Definition 7 given the additional regularity. All the results on the absolutely continuous case are collected in Theorem 12.

  • Section 4 contains a detailed study of the steady states for (3), as well as of the minimisers of (4). A characterisation of the steady states is provided in Proposition 9. A consequent result concerning the asymptotic behaviour is provided in Theorem 14.

  • Section 5 describes two relevant examples of atomic initial data. In both cases, non uniqueness of weak measure solutions is shown, and the relevant gradient flow solution is singled out as well. These two examples allow to conclude interesting properties related with the occurrence or not of the smoothing effect (or lack thereof) of initial atomic parts, see Remarks 9 and 10. The link with the hyperbolic system (5) is described in detail in Subsection 5.3 leading to several open problems.

2. Preliminaries

This section is devoted to setting up the framework to show existence and uniqueness of solutions to the system

with Newtonian interactions, . Throughout the paper, and will be considered as time dependent curves with values on the set of probability measures on . System (6a) is equipped with initial data

for some . Moreover, we write to denote the set of probability measures with finite second moment, i.e.

In the following we shall use the symbol referring to elements of which are absolutely continuous with respect to the Lebesgue measure. Consider a measure and a Borel map . We denote by the push-forward of through , defined by

for all Borel sets . We refer to as the transport map pushing to . Next let us equip the set with the -Wasserstein distance. For any measures it is defined as


where is the class of transport plans between and , that is,

where , , denotes the projection operator on the component of the product space . Setting as the class of optimal plans, i.e. minimisers of (7), the Wasserstein distance becomes

for any . The set equipped with the -Wasserstein metric is a complete metric space which can be seen as a length space, see for instance [1, 54, 57, 58].

Remark 1.

Given two measures , by using the elementary inequality and the above definition of the 2-Wasserstein distance, one can easily show the inequality

which will be used later on.

Next we introduce the notion of the Fréchet sub-differential in the space of probability measures.

Definition 1 (Fréchet sub-differential in ).

Let be a proper and lower semicontinuous functional, and let . We say that belongs to the Fréchet sub-differential at , denoted by , if

Moreover, if we denote by the element of minimal -norm in .

This definition will play a crucial role when introducing the notion of gradient flow solutions to system (6) later on, cf. Section 3.2.

A curve is a constant speed geodesic if for all . Due to [1, Theorem 7.2.2], a constant speed geodesic connecting and can be written as

where and thus and . In the literature is also known as McCann interpolation, cf. [50]. Next, we introduce a modified notion of convexity, the so-called -geodesic convexity, which is of paramount importance in the study of gradient flows in the metric space .

Definition 2 (-geodesic convexity).

Let . A functional is said to be -geodesically convex in if for every there exists such that

for any .

It is necessary to recall that the -geodesic convexity is strictly linked with the concept of k-flow.

Definition 3 (-flow).

A semigroup is a -flow for a functional with respect to if, for an arbitrary , the curve is absolutely continuous on and satisfies the evolution variational inequality (E.V.I.)


for all , with respect to any reference measure such that .

As already mentioned, the previous concepts of -flow and -convexity are closely intertwined. Indeed, a convex functional possesses a uniquely determined flow for , and if a functional possesses a flow, then it is convex with , cf. Refs. [1, 37, 49], for further details. The notion of k-flow will be of great help in the use of the flow interchange technique, cf. Subsection 3.2.

Now, let be an absolutely continuous curve in . We can define the metric derivative of as

which is well-defined almost everywhere since is an absolutely continuous curve, cf. [1, 54].

As system (6) describes the evolution of two interacting species, it is necessary to work on the product space equipped with the 2-Wasserstein product distance, which is defined in the canonical way

for all belonging to . Now, let us introduce another crucial tool in our paper. For a given its cumulative distribution function is given by


Since is a non-decreasing, right continuous function such that

we may define the pseudo-inverse function associated to , by


for any . It is easy to see that is right-continuous and non-decreasing as well. Having introduced the pseudo-inverse, let us now recall some of its important properties. First we notice that it is possible to pass from to as follows


For any probability measure and the pseudo-inverse , , associated to it, we have


for every bounded continuous function . Moreover, for , the Hoeffding-Fréchet theorem [53, Section 3.1] allows us to represent the 2-Wasserstein distance, , in terms of the associated pseudo-inverse functions as


since the optimal plan is given by , where is the Lebesgue measure on the interval , cf. also [57, 30]. We have seen that for every we can construct a non-decreasing according to (10), and by the change of variables formula (12) we also know that is square integrable. We now recall that this mapping is indeed a distance-preserving bijection between the space of probability measure with finite second moments and the convex cone of non-decreasing -functions

Proposition 1 ( is an isometry).

The map


mapping probability measures onto the convex cone of non-decreasing -functions is an isometry.

Let us introduce the notion of sub-differential for functions in .

Definition 4 (Fréchet sub-differential in ).

For a given proper and lower semi-continuous functional on we say that belongs to the sub-differential of at if and only if

as , with the notation . The sub-differential of at is denoted by , and if then we denote by the element of minimal -norm of .

Remark 2 (Mass normalisation).

In the most general situation possible, the two species, and , have different masses . The change of variables

allows to rewrite system (6a) (by dropping the tildes) as

and the gradient flow structure in the product Wasserstein metric would be lost because the two interspecific potentials are different. However, this problem can be overcome by using a weighted version of the product distance of the form

as done in [36]. As these multiplying constants do not bring significant technical difficulties (while making the notation much heavier), for the sake of convenience we shall assume throughout the whole paper that .

3. Existence and Uniqueness

In this section we provide the mathematical theory for system (6). In the first subsection we will deal with the case of general probability measures in as initial conditions. In the second subsection we shall restrict ourselves to the case of measures that are absolutely continuous with respect to the Lebesgue measure. In the former we will provide a notion of solutions that is linked to the concept of gradient flows in Hilbert space a-la Brézis [15], working with the pseudo-inverse formulation of the problem. In the latter, a better regularity can be achieved and the theory is developed in the framework of gradient flows in Wasserstein space, [1]. Before entering the details, let us recall the definition of interaction energy functional in (4): for all we set

which is well defined due to the control on the second order moment.

3.1. General Measures Initial Data

In this first subsection we will use the concept of -gradient flow by studying system (1) in terms of the pseudo-inverse functions and defined in Section 2. Throughout the rest of this section we set and to simplify the notation. Hence, system (6) (formally) becomes


for and , cf. [11, 30, 47] for similar computations. In order to give a meaning to the above system in the case of or having atoms, we use the convention . In terms of the pseudo inverses and , the functional becomes


In the remainder of this section we shall see that (16) is the -gradient flow associated to (an extended version of) the energy functional (17).

Remark 3.

Later on in the paper we need to distinguish between the self-interaction part of and its cross-interaction part. Thus, let us rewrite as

where is the energy functional arising from the self-interactions and is associated to the cross-interaction.

Following the procedure of [12], since we are dealing with distribution of particles, we have to ensure that the flow remains in the set , with defined in (14), see also [10, 13]. Hence, for given, we introduce the the indicator function of , defined by,


Thus we consider the extended functional


In [12, Proposition 2.8] the authors proved that the self-interaction part of , i.e. , is actually linear when restricted to . Let us recall this result in the next proposition.

Proposition 2.

Let . Then

As a trivial consequence of Proposition 2 we have the following result.

Proposition 3.

The functional is convex on .


The proof is trivial since , the cross-interaction part of , is convex due to the convexity of the Newtonian potential . Moreover, from [12, Proposition 2.9] we argue that the remaining part is convex. ∎

We now present the definition of -gradient flow solutions to system (16).

Definition 5.

Let . An absolutely continuous curve , , is a gradient flow for the functional if is a Lipschitz function on , i.e. (in the sense of distributions) and if it satisfies the sub-differential inclusion


for every with .

We observe that the assumption is natural as (respectively ) is the pseudo-inverse of the cumulative distribution of the initial measure (respectively ). We also observe that this assumption easily implies .

Remark 4.

The gradient flow notion defined in Definition 5 is taken from the book [15, Theorem 3.1]. Actually, in [15, Theorem 3.1] the following extra condition is required, namely


According to [15, Theorem 3.1] a solution in the sense of Definition 5 in conjunction with (21) directly verifies the following properties:

  1. admits a right derivative for every and

    for every ;

  2. the function is right continuous and the function is non-increasing;

  3. if and are two solutions to system (20), then there holds

    for all .

On the other hand, [39, Section 9.6, Theorem 3] shows that condition (21) can be avoided in order to prove uniqueness. In fact, the estimate (21) can be proven as a consequence of the properties stated in Definition 5.

Since we know that is a proper, lower semi-continuous, and convex functional on the Hilbert space , it easy to show that is a maximal monotone operator. Thus we can apply the theory of Brézis [15, Theorem 3.1] combined with [39, Section 9.6, Theorem 3] in order to prove existence and uniqueness of an absolutely continuous curve satisfying the differential inclusion above.

Theorem 6.

Let . There exists a unique solution in the sense of Definition 5 with initial datum .

Now, let us go back to system (6) and state our definition of solution.

Definition 7.

Let . An absolutely continuous curve is a gradient flow solution to system (6) if the pseudo-inverses of the space cumulative distribution functions associated to are a solution to system (16) in the sense of Definition 5 with initial datum .

According to Definition 7 the following theorem is then a consequence of the isometry (15) and Theorem 6.

Theorem 8.

Let . There exists a unique solution to the system (6) in the sense of Definition 7.

So far we assumed that the link between (6) and (16) is somewhat natural and we just referred to similar situations in the literature. However, the theory developed in this subsection would be somewhat meaningless if we did not show that the concept of solution in Definition 5 extends a more classical notion of solution for (6). The following subsection is dedicated to establishing exactly this link.

3.2. Absolutely Continuous Initial Data

In this subsection we consider the case of densities as initial data. Following the approach of [1] combined with the results from [11, 12, 25, 26, 36], we pose system (6) as the gradient flow of the interaction energy functional (4) (that we recall here for the reader’s convenience)

for all , and , for all .

Definition 9.

Given any , an absolutely continuous curve is a gradient flow for if and solve the following system in the distributional sense


with initial datum and the velocity field such that

for and

for a.e. .

Note that it is easy to check that the element of minimal norm in is given by

Using the results obtained in [12, 25, 26, 36], we easily get the following proposition.

Proposition 4.

The functional is -geodesically convex on for all . Moreover, for all the vector field


is the unique element of the minimal Fréchet sub-differential of , where


The geodesic convexity of on is the consequence of two observations. First, the cross-interaction part is geodesically convex as the interspecific interaction potential is given by , a convex function. Second, the geodesic convexity of the intraspecific self-interaction part can be proven using a nice monotonicity property of the transport map between two measures in in one dimension. Indeed, in [26, Lemma 1.4] the authors prove that, given , the transport map is essentially non-decreasing, i.e. it is non-decreasing except on a -null set. Thus, when interpolating measures, it is sufficient to consider the restriction of on , which is convex. Hence, we can prove that is geodesically convex. Formula (23) can be proven as in [25, Proposition 2.2]. ∎

Remark 5.

We highlight that in the presence of atomic parts for or the sub-differential may be empty, as shown in [12].

Recall that, for , the slope of a functional on is defined as

and it can be written as

under certain conditions, cf. [1, Chapter 10].

Definition 10.

An absolutely continuous curve is a curve of maximal slope for the functional if the map is an absolutely continuous function and the following inequality holds

for all .

In order to construct a solution to system (6) in the sense of Definition 9 we follow the strategy proposed in [1] and used in [25, 36]. First, we prove the existence of a curve of maximal slope by means of the so-called “Minimizing Movement Scheme”, cf. [1, 34], or Jordan-Kinderlehrer-Otto scheme, cf. [43]. Then, we prove the limit curve of the scheme is absolutely continuous w.r.t. the Lebesgue measure provided the initial datum is in for some .

Let be a fixed time step and let be a fixed initial datum such that . We define a sequence recursively. We set and, for a given with , we choose as


Note that (24) is well-posed arguing as in [25, Lemma 2.3 and Proposition 2.5] for each component. Next we define the piecewise constant interpolation of the sequence