Aggregation-diffusion equations with non-linear mobility
Seminar, NTNU, Trondheim
Collaborators
\[ \newcommand{\diff}{d} \newcommand{\dr}{\diff r} \newcommnand{\R}{{\mathbb R}} \newcommand{\Rd}{{\mathbb R^d}} \newcommand{\diver}{\nabla \cdot} \newcommand{\calF}{\mathcal{F}} \newcommand{\calU}{\mathcal{U}} \newcommand{\defeq}{\overset{\mathrm{def}}{=}} \newcommand{\ee}{\varepsilon} \newcommand{\mob}{\mathrm{m}} \newcommand{\Wass}{\mathfrak{W}} \]
Aggregation-Diffusion
\[\begin{equation*} \label{eq:ADE} \tag{ADE} \frac{\partial \rho}{\partial t} = \diver \Big( \mob(\rho) \nabla ( U'(\rho) + V + W * \rho ) \Big). \end{equation*}\] The mobility \(\mob\) is often linear, i.e., \(\mob(\rho) = \rho\).
Modelling
For the introductory part we will follow [Gómez-Castro, 2024].
Continuity equations
Let \(\rho\) be a density and \(\omega \subset \Rd\) any control volume, if \(\mathbf F\) is the out-going flux \(\frac{\diff }{\diff t} \int_\omega \rho \diff x =- \int_{\partial \omega} \mathbf F \cdot \mathbf n \diff S= -\int_{ \omega} \diver \mathbf F \diff x\).
Hence we arrive at the \(\frac{\partial \rho}{\partial t} = -\diver \mathbf F.\)
Heat equation: Fourier’s law to model the flux \(\mathbf F = -D \nabla \rho\) \[ \frac{\partial \rho}{\partial t} = D\Delta \rho . \]
Porous-Medium Equation: \(\mathbf F = -\rho \mathbf v\); Darcy’s law \(\mathbf v = - \frac{k} \mu \nabla p\);
general state equation \(p = \phi(\rho)\). For perfect gases \(\phi (\rho) = p_0 \rho^\gamma\). \[ \frac{\partial \rho}{\partial t} = \Delta \rho^m. \]
Particle systems
Non-interacting particles
Consider \(N\) particles moving in a velocity field \(\frac{\diff X_i}{\diff t} = \mathbf v(t, X_i(t))\)
Define the empirical distribution of equal masses \(1/N: \quad\) \(\mu_t^N = \sum_{j=1}^N \frac 1 N \delta_{X_j(t)}\)
This is a weak solution to the continuity equation
\[ \partial_t \mu^N + \diver (\mu^N \mathbf v ) = 0. \]
Interacting particles: Aggregation and Confinement
Non-local interactions: \(\frac{\diff X_i}{\diff t} = -\sum_{\substack{ j=1 \\ j \ne i}}^N \frac 1 N \nabla W (X_i - X_j)- \nabla V (X_i)\).
We recover Aggregation-Confinement Equation
\[ \partial_t \mu^N = \diver (\mu^N \nabla ( W * \mu^N + V ) ). \]
No-flux conditions
In a bounded domain \(\Omega\), we preserve mass if we set \(\mathbf F \cdot \mathbf n = 0\) on \(\partial \Omega\), i.e., \[\begin{equation*} \tag{BC} \label{eq:BC Omega} \rho \nabla ( U'(\rho) + V + W*\rho ) \cdot \mathbf n = 0 \qquad \text{ for } t > 0 \text{ and } x \in \partial \Omega \end{equation*}\] (and set \(\rho = 0\) in \(\Rd \setminus \Omega\) for the convolution).
Special cases
Linear diffusion
Heat equation \(\quad \partial_t \rho = \Delta \rho_t\)
It corresponds to \(\eqref{eq:ADE}\): \(U = \rho \log \rho\) and \(V = W = 0\).
In \(\Rd\) admits the Gaussian solution \(\rho_t = K_t * \rho_0,\) \[ \rho_t \to 0, \qquad \|\rho_t\|_{L^1} = M , \qquad \| \rho_t - M K_t \|_{L^1} \to 0 \quad\text{ as }\, t \to \infty. \] A good survey on the matter can be found in [Vázquez, 2017].
Linear Fokker-Planck. A self-similar change of variable leads to the linear Fokker-Planck equation \[ \partial_t u = \Delta u_t + \nabla \cdot (x u_t) \]
ADE: \(U_1 (\rho) = \rho \log \rho, V = \frac{|x|^2}{2} \text{ and }W = 0.\)
Asymptotic behaviour: \(u_t \overset{t \to \infty}{\longrightarrow} G\).
Porous-medium equation \(\quad \frac{\partial \rho}{\partial t} = \Delta \rho^m\)
The range \(m \in (0,1)\) is sometimes called Fast Diffusion Equation.
ADE: \(U_m (\rho) = \frac{\rho^m}{m-1}, \quad V = W = 0.\)
Self-similar solution \[ B_t = t^{-\alpha} (C - k|x|^2t^{-2\beta})^{\frac 1 {m-1}}_+, \] where \(\alpha = \tfrac{d}{d(m-1)+2},\, \beta = \tfrac{\alpha}{d},\, k = \tfrac{\alpha(m-1)}{d}\).
(see [Vázquez, 2006b])
We get the self-similar solution from a change to self-similar variables \[ \frac{\partial u}{\partial t} = \Delta u^m + \nabla\cdot(x u). \]
Keller-Segel model
The Keller-Segel [Keller & Segel, 1970] (or Patlak-Keller-Segel model [Patlak, 1953]) model for describes the motion of cells by chemotactical attraction by means of the coupled system \[\begin{equation*} \tag{KSE} \label{eq:Keller-Segel un-normalised} \begin{dcases} \frac{\partial \rho}{\partial t} = \Delta \rho - \diver (\rho \nabla V_t ) \\ -\Delta V_t= \rho_t . \end{dcases} \end{equation*}\] This model was first studied for \(d = 2\).
ADE using Fourier transform: \(V_t = W_{\rm N} * \rho\) where \(W_{\mathrm N} = \mathrm{F}^{-1}[|\xi|^{-2}].\)
Some authors replace \(\Delta \rho\) by \(\Delta \rho^m\):
[Calvez & Carrillo, 2006], [Bedrossian, Rodríguez & Bertozzi, 2011],[Luckhaus & Sugiyama, 2007], [Sugiyama & Kunii, 2006], [Sugiyama, 2007a], [Sugiyama, 2007b].
Chapman-Rubinstein-Schatzman
[Chapman, Rubinstein & Schatzman, 1996] introduced a mean-field model from Gizburg-Landau in the context of superconductivity \[\begin{equation*} \frac{\partial \rho}{\partial t} = \diver( |\rho| \nabla V_t ) \qquad \qquad -\Delta V_t + V_t = \rho_t . \end{equation*}\]
As long as \(\rho \ge 0\) it is ADE: \(U, V = 0, W = \mathrm{F^{-1}}[ (1 + |\xi|^2)^{-1} ]\).
[Huang & Svobodny, 1998], [Schätzle & Styles, 1999] [Ambrosio & Serfaty, 2008], [Ambrosio, Mainini & Serfaty, 2011]
Newtonian vortex problem
In [Lin & Zhang, 1999] taking a different limit from Gizburg-Landau equations \[\begin{equation*} \frac{\partial \rho}{\partial t} = \diver (\rho \nabla V_t ) \qquad \qquad -\Delta V_t = \rho. \end{equation*}\]
[Masmoudi & Zhang, 2005], [Bertozzi, Laurent & Léger, 2012].
Caffarelli-Vázquez
[Caffarelli & Vazquez, 2011] introduces a non-local porous medium-type equation given by \[\begin{equation*}
\begin{dcases}
\frac{\partial \rho}{\partial t} = \diver (\rho \nabla V_t ) \\
(-\Delta)^s V_t = \rho ,
\end{dcases}
\end{equation*}\] where \((-\Delta)^s\) denotes the fractional Laplacian. ADE: Riesz potential \(W = \mathrm F^{-1}[|\xi|^{-2s}]\).
[Caffarelli & Vázquez, 2011], [Biler, Imbert & Karch, 2015], [Caffarelli, Soria & Vázquez, 2013],[Caffarelli & Vázquez, 2015], [Serfaty & Vázquez, 2014], [Lisini, Mainini & Segatti, 2018].
Porous-Medium Equation with non-local pressure
\[ \frac{\partial \rho}{\partial t} = \diver (\mob(\rho) \nabla V_t ) \qquad \text{where } \mob(\rho)=\rho^{m-1} \]
[Stan, Del Teso & Vázquez, 2016], [Stan, Del Teso & Vázquez, 2019],
[del Teso & Jakobsen, 2023]
The power-type family of Aggregation-Diffusion Equations
Many authors have devoted their attention to the power-type family of non-linearities given by \[ \begin{gathered} U = \frac{\rho^m}{m-1}, \qquad \qquad V_\lambda (x) = \begin{dcases} \frac{|x|^\lambda}{\lambda} & \text{if } \lambda \ne 0 \\ \log|x| & \text{if } \lambda = 0 , \end{dcases}, \\ W_k (x) = \begin{dcases} \frac{|x|^k}{k} & \text{if } k \in (-d,d) \setminus \{0\} \\ \log|x| & \text{if } k = 0 , \end{dcases} \end{gathered} \]
In bounded domains we can replace \(W*\rho\) by \(\int_{\Omega} K(x,y) \rho(y) dy\).
More examples
- McKean-Vlasov and Kuramoto models
- Neural networks of machine learning [Fernández-Real & Figalli, 2022]
See [Gómez-Castro, 2024].
Numerics
Write linear-mobility ADE \(\mob(\rho) = \rho\) as \[ \partial_t \rho = - \nabla \cdot F, \qquad F = \rho v, \qquad v = -\nabla \xi, \qquad \xi = U'(\rho) + V + \int_\Omega K(\cdot, y)\rho(y)dy. \]
Consider the implicit up-winding finite-volume scheme [Bailo, Carrillo & Hu, 2020] \[ \begin{aligned} \frac{\rho_i^{ n+1} - \rho_i^{n}}{\Delta t} &= - \frac{F^{}_{i+\frac 1 2} (\rho^{n+1}) - F^{}_{i-\frac 1 2} (\rho^{ n+1})}{\Delta x} , \qquad i \in \{1, \cdots, N\}, n \in \mathbb N , \\ F_{i+\frac 1 2}^{} (\rho) &= \rho_i (v_{i+\frac 1 2}^{} (\rho))^+ + \rho_{i+1} (v_{i+\frac 1 2}^{}(\rho))^- , \\ v_{i+\frac 1 2}^{} (\rho) &= -\frac{\xi_{i+1}^{}(\rho) - \xi_{i}^{}(\rho)}{\Delta x} , \\ \xi_i^{}(\rho) &= U_{}'(\rho_i) + V (x_i) + \sum_{j} K_{ij} \rho_j, \\ F_{\frac 1 2}^{} (\rho) &= F_{N + \frac 1 2}^{} (\rho) = 0 . \end{aligned} \]
Key properties
- Well-posed and mass preserving
- Free energy that decays [Bailo, Carrillo & Hu, 2020]
- Is convergent under some conditions [Bailo, Carrillo, Murakawa & Schmidtchen, 2020]
- When \(W = 0\) the scheme is monotone, and therefore there is comparison principle
[Carrillo, Fernández-Jiménez & Gómez-Castro]
(Some) Non-linear mobilities
When \(F = \mob(\rho) v\) and up-winding is more difficult.
When \(\mob(\rho) = \mob_{\uparrow}(\rho) \mob_{\downarrow}(\rho)\) with \(\mob_{\uparrow}\) non-decreasing and \(\mob_{\downarrow}\) non-increasing.
[Bailo, Carrillo & Hu, 2023] proposed:
\[ F_{i+\frac 1 2}^{} (\rho) = \mob_{\uparrow}(\rho_i) \mob_{\downarrow} (\rho_{i+1}) (v_{i+\frac 1 2}^{} (\rho))^+ + \mob_{\uparrow} (\rho_{i+1}) \mob_{\downarrow} (\rho_{i}) (v_{i+\frac 1 2}^{}(\rho))^- . \\ \]
Asymptotic behaviour
Long-time asymptotic behaviour: \(\rho_t \to \widehat \rho\) as \(t \to \infty\)
Blow-up: \(\rho_t \rightharpoonup \widehat \mu \not \in L^1(\Omega)\) as \(t \nearrow T \in (0,\infty]\).
- We say there is concentration if \(\widehat \mu\) contains a Dirac
Extinction: \(\rho_t \searrow 0\) as \(t \nearrow T \in (0,\infty]\)
- If there is \(B_t\) such that \(\frac{\| \rho_t - B_t \|_{L^p}}{\|\rho_t\|_{L^p}}\) we say there are intermediate asymptotics
Formal steady states
ADE: \(\qquad \partial\rho_t = \nabla \cdot( \rho \nabla (U'(\rho) + V + W*\rho) )\)
Steady state: \(\qquad \rho \nabla (U'(\rho) + V + W*\rho) = 0.\)
This suggests that, for steady states, \[ U'(\rho) + V + W * \rho = C_D \qquad \forall D \subset \mathrm{supp} \rho \text{ open and connected}. \]
The constant can naturally differ between components.
Example. Linear Fokker-Planck. \(\mathrm{supp} \rho=\Rd\)
Example. PME with \(m > 1\) in self-similar variables. Barenblatt profile.
Intermediate asymptotics
Heat equation:
in \(\Rd\): Gaussian self-similar solution [Vázquez, 2017]
In bounded domains: from spectral decomposicion \(a_1 e^{-\lambda_1 t} \varphi_1\).
Porous-Medium Equation: Barenblatt solution \[\begin{equation} \label{eq:intermediate asymptotics} \frac{\| \rho_t - B_t \|_{L^p}}{\| \rho_t \|_{L^p}} \to 0, \qquad \text{as } t \to \infty. \end{equation}\] See [Vázquez, 2006b], [Vázquez, 2006a].
\(\frac{\partial \rho}{\partial t} = \diver(\rho \nabla V )\): blow-up
We can look at the examples \(V = \frac{|x|^\alpha}{\alpha}\) with \(\alpha > 0\), i.e. \(v = -\nabla V = - |x|^{\alpha-2} x\).
Then, \(\rho_t = (X_t)_\sharp \rho_0\) with characteristic field the solutions to \[ \frac{\partial X_t}{\partial t} = -|X_t|^{\alpha-2} X_t, \qquad \qquad X_t(0) = y. \]
This can be explicitly solved as \[ X_t (y) = \begin{dcases} y e^{-t} & \text{if }\alpha=2, \\ ( |y|^{2-\alpha} - (2-\alpha) t )^{\frac 1 {2-\alpha}} \frac{y}{|y|} & \text{if } \alpha \ne 2. \end{dcases} \]
For \(\alpha < 2\) the characteristics arrive at \(0\) at finite time.
We get a Dirac delta at \(0\) at a certain finite time for any \(\rho_0 \not\equiv 0\).
Keller-Segel problem: conditional blow-up
Some authors noticed by direct techniques that the Keller-Segel model has finite-time blow-up for initial data with large enough mass:
[Jäger & Luckhaus, 1992] showed for \(\eqref{eq:Keller-Segel un-normalised}\) in a bounded domain with no-flux
\(\exists M^*\) such that if \(\| \rho_0 \|_{L^1} > M^*\), at \(T^* < \infty\) we have then \(\rho_{T^*} = a \delta_0 + \cdots\) .In \(\Rd\), [Herrero & Velázquez, 1996] showed \(M^* = 8 \pi\).
In more generality in [Dolbeault & Perthame, 2004] \[ \frac{\diff}{\diff t} \int_{\mathbb R^2} |x|^2 \rho(t,x) \diff x = 4 M \Big( 1- \frac{M}{8\pi }\Big), \qquad \qquad M = \int_{\mathbb R^2} \rho_0 . \] [Blanchet, Carrillo & Laurençot, 2009] studied \(d > 2\)Conversely, for \(\| \rho_0 \|_{L^1} < M^*\) global existence for \(\eqref{eq:Keller-Segel un-normalised}\) is known.
In \(d = 2\) it was shown in [Blanchet, Dolbeault & Perthame, 2006].
Free energy and gradient-flow
Free energy: a Lyapunov functional
Consider the free energy \[\begin{equation*} \tag{FE$^*$} \label{eq:free energy star} \mathcal F(\rho) =\int_\Omega U(\rho) + \int_\Omega V \rho + \frac 1 2 \iint_{\Omega \times \Omega } K(x,y) \rho(x) \rho(y) \diff x \diff y . \end{equation*}\] Its first variation in the pressure field we usually encounter \[\begin{equation} \label{eq:first variation} \frac{\delta \mathcal F}{\delta \rho} [\rho] \defeq U'(\rho) + V + \int_\Omega K(\cdot,y)\rho(y) \diff y. \end{equation}\]
We write ADE as \[\begin{equation*} \partial_t \rho = \diver\left( \mob(\rho) \nabla \frac{\delta \mathcal F}{\delta \rho} [\rho]\right) \end{equation*}\]
Then, the free energy is a Lyapunov functional \[\begin{equation*} \frac{d}{dt} \mathcal F[\rho_t] = \int \frac{\delta \mathcal F}{\delta \rho} [\rho] \partial_t \rho_t = - \int \mob(\rho_t) \left|\nabla \frac{\delta \mathcal F}{\delta \rho} [\rho] \right|^2. \end{equation*}\]
Gradient-flow structure
When \(\mob(\rho) = \rho\) there a nice gradient-flow theory:
Otto calculus (see, e.g., [Ambrosio, Brué & Semola, 2021]):
ADE is formally the \(2\)-Wasserstein gradient flow \[ \frac{\partial \rho}{\partial t} = - \nabla_{\Wass_2} \mathcal F[\rho_t]. \]There is a good notion of convexity (McCann’s displacement convexity).
Under certain condition we have \(\Wass_2(\mu_t, \widehat \mu_t) \le e^{-\lambda t} \Wass_2(\mu_0, \widehat \mu_0).\)
For general \(\mob(\rho)\) there is suitable Wasserstein distance
[Dolbeault, Nazaret & Savaré, 2009], [Carrillo, Lisini, Savaré & Slepčev, 2010]
This is only a distance if \(\mob\) is concave.
It is (almost) never displacement convex in our settings.
Under convexity, the asymptotic behaviour are usually the local minimisers.
Local minimisers: Euler-Lagrange condition
If \(\calF(\widehat \rho) \le \calF(\rho)\) for all \(\rho \ge 0\) such that \(\int \rho = \int \widehat \rho\) we have that, for some \(C \in \mathbb R\)
\[\begin{equation*} \begin{dcases} \frac{\delta \calF}{\delta\rho} [\widehat \rho] = C & \text{if } \rho > 0, \\ \frac{\delta \calF}{\delta\rho} [\widehat \rho] \ge C & \text{if } \rho = 0. \\ \end{dcases} \end{equation*}\]
Recall \(\frac{\delta \calF}{\delta\rho} [\rho] = U'(\rho) + V + W*\rho\).
The free-energy minimiser is the steady state with same constant in all connected components.
\[ \widehat \rho(x) = \max\Bigg\{0, (U')^{-1} \big(C-V(x) - W*\widehat \rho\big)\Bigg\}. \]
Notice that when \(W = 0\) we the solution is parametrised by \(C\).
Many steady states, only one free-energy minimiser
Global minimizer may not be the global attractors
Newtonian vortex with
power-type non-linear mobility
\[\begin{equation} \tag{NVE} \begin{dcases} \frac{\partial \rho}{\partial t} = \diver (\rho^\alpha \nabla V_t ) \\ -\Delta V_t = \rho. \end{dcases} \end{equation}\]
[Carrillo, Gómez-Castro & Vázquez, 2022a],[Carrillo, Gómez-Castro & Vázquez, 2022c]
Self-similar analysis \(\qquad \frac{\partial \rho}{\partial t} = \diver (\rho^\alpha \nabla V_t ), \quad-\Delta V_t = \rho\)
We can look for solutions of the form \[ B(t,x)= t^{-\gamma} F(|x|\,t^{-\beta}). \]
Plugging this in the equation we recover the self-similar solution of mass 1 for \(\alpha \in (0,1)\) \[ \scriptsize B(t,x) = t^{-\frac 1 \alpha}\left( \alpha + \left( \frac{ \omega_d |x|^dt^{-\frac 1 {\alpha}} } { \alpha} \right)^{\frac {\alpha} {1-\alpha }} \right)^{-1/\alpha }. \]
The case \(\alpha = 1\) was already studied in [Serfaty & Vázquez, 2014].
The same algebra works for \(\alpha > 1\) but gives no finite-mass solutions
Mass variable \(\qquad \frac{\partial \rho}{\partial t} = \diver (\rho^\alpha \nabla V_t ), \quad-\Delta V_t = \rho\)
Let us define \[ m(t,x) = \int_{-\infty}^x \rho(t,y) \diff y \]
Integrating the equation for \(\rho\) in \(x\): \(\displaystyle\qquad \frac{\partial m}{\partial t} = \rho^\alpha \frac{\partial V_t}{\partial x}.\)
The equation for \(V_t\) is \(-\frac{\partial^2 V}{\partial x^2} = \rho = \frac{\partial m}{\partial x}\). Setting \(\frac{\partial V}{\partial x} (-\infty) = 0\) we get \(-\frac{\partial V}{\partial x} = m\), so
\[ \frac{\partial m}{\partial t} = -\left( \frac{\partial m}{\partial x} \right)^\alpha m \]
If \(d > 1\) and \(\rho\) is radially symmetric, we get to the same equation for \[ m(t,v) = \int_{A_v} \rho(t,x) \diff x , \quad A_v = B(0,r) \text{ such that }|A_v| = v. \]
Mass conservation implies \(m(t,0) = 0, m(t,\infty) = m_0(\infty)\).
Time asymptotics for \(m\)
Assume that \(0 \le \rho_0 \in L^\infty_c(\Rd)\) is radially symmetric, and let \(M = \| \rho_0 \|_{L^1}\).
Let \(\alpha \in (0,1)\) and denote \[ G_M(\kappa)=\int \limits _{\omega_d |y|^d \le \kappa } F_M (|y|) \diff y \] Then, \[\begin{equation} \label{eq:asymp mass} \sup_{\kappa \ge \ee } \left| \dfrac{ m (t, t^{\frac 1 \alpha} \kappa ) } {G_M (\kappa)} - 1 \right| \longrightarrow 0 , \end{equation}\] as \(t \to \infty\) for any \(\ee > 0\).
Let \(\alpha > 1\) and define \[\begin{equation} G( y ) = \begin{dcases} y & y\le 1 \\ 1 & y > 1. \end{dcases} \end{equation}\] Then, \[\begin{equation} \sup_ { y \ge \varepsilon } \left| \frac{m\left( t , M (\alpha t)^{\frac 1 \alpha} y \right) }{M G (y)} - 1\right| \to 0, \end{equation}\] as \(t \to +\infty\) for any \(\varepsilon > 0\).