Nonlocal Degenerate Diffusion-Aggregation PDE

• The topic is a nonlinear PDE system that models the interplay of vanishing diffusion and nonlocal aggregation, capturing complex dynamics in collective phenomena.
• It employs advanced mathematical tools such as Wasserstein gradient flow, convexity analysis, and bifurcation theory to study existence, uniqueness, and phase transitions.
• Robust numerical schemes that ensure positivity and energy dissipation enable detailed simulations and precise analysis of steady states and dynamic behavior.

A multi-dimensional nonlocal degenerate diffusion-aggregation equation is a nonlinear PDE system featuring a degenerate diffusion operator—vanishing for small values of the unknown—and a nonlocal aggregation or interaction effect, typically modeled via convolution with an interaction kernel. Such systems are canonical models for collective phenomena in physics, biology, and social dynamics where both dispersive (diffusive) and attractive (aggregative) forces act at different length scales and their interplay determines pattern formation, equilibrium, and the overall dynamics. The mathematical analysis of these PDEs combines several advanced tools: gradient-flow theory in Wasserstein space, convexity analysis of free-energy functionals, bifurcation and phase-transition theory, sharp existence and blow-up dichotomies, as well as the regularity and uniqueness of compactly supported steady states. The field has advanced through rigorous results on equilibrium selection, long-time asymptotics, moment estimates, sharp blow-up criteria, and computational schemes that respect the underlying energetic structure.

1. Structural Form and Mathematical Framework

A prototypical multi-dimensional nonlocal degenerate diffusion-aggregation equation for a scalar density $\rho(x,t)\ge0$ in $\mathbb{R}^d$ (or a domain $\Omega\subset\mathbb{R}^d$) has the general structure: $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$ where $m>1$ is the degeneracy exponent, and $K$ is a symmetric, typically radial, interaction kernel that models nonlocal pairwise attraction (or repulsion). The operator $\nabla\rho^m$ produces nonlinear degenerate diffusion (vanishing where $\rho=0$). The convolution $K*\rho$ is the nonlocal potential generated by the density profile. The equation admits an interpretation as a $2$-Wasserstein gradient flow of the free-energy

$\mathbb{R}^d$0

with the first term promoting spreading and the second promoting aggregation. This structure underlies both the variational (energy) and dynamic (PDE) properties of the system (Matthes et al., 2024, Carrillo et al., 2019, Kaib, 2016).

For systems (multiple species), cross-diffusion and cross-interactions are included: $\mathbb{R}^d$1 with additional energy terms representing cross-interactions, as studied in (Matthes et al., 2024).

The interaction kernel $\mathbb{R}^d$2 may be bounded, integrable, or singular (e.g., $\mathbb{R}^d$3 or a Riesz kernel $\mathbb{R}^d$4), and its properties (regularity, convexity, decay) crucially affect the qualitative analysis (Zhou et al., 21 Dec 2025, Zhang, 2018).

2. Existence, Uniqueness, and Regularity Theory

Existence and uniqueness of (weak) solutions rely primarily on energy-dissipation properties and compactness arguments. Suppose the free energy is displacement convex (or suitably split into convex and small nonconvex parts), then JKO (minimizing-movement) schemes in Wasserstein space yield global weak solutions that dissipate energy in time (Matthes et al., 2024, Bailo et al., 2018). In the scalar case, the existence of unique compactly supported global minimizers for $\mathbb{R}^d$5 or for $\mathbb{R}^d$6 with subcritical interaction parameters is well established (Kaib, 2016, Delgadino et al., 2019).

Uniqueness of steady states in the class of radially decreasing densities fundamentally depends on the degeneracy exponent. Specifically:

• For $\mathbb{R}^d$7 and attractive $\mathbb{R}^d$8, there exists at most one radially decreasing compactly supported steady state (up to translation) for a fixed mass (Delgadino et al., 2019, Kaib, 2016).
• For $\mathbb{R}^d$9, nonuniqueness can occur, and multiple compactly supported equilibria may exist, depending on the potential and parameter regime (Delgadino et al., 2019).

Regularity theory yields $\Omega\subset\mathbb{R}^d$0 bounds for solutions in the subcritical and critical regimes, as well as Hölder continuity for $\Omega\subset\mathbb{R}^d$1, leveraging degenerate parabolic regularity and De Giorgi–Nash–Moser iteration (Zhang, 2018, Chayes et al., 2012). When the diffusion degenerates at zero, the analysis employs restriction to the positivity set, technical cut-offs, and energy arguments to prevent singularities (Matthes et al., 2024, Kaib, 2016, Eckardt et al., 2024).

3. Equilibrium Structure, Bifurcation, and Phase Transitions

The free-energy landscape for these systems exhibits a variety of behaviors as parameters vary:

• Subcritical regime: For sufficiently strong diffusion ($\Omega\subset\mathbb{R}^d$2 large, or weak aggregation), the only minimizer is the spatially homogeneous state, corresponding to the uniform density (Carrillo et al., 2019, Chayes et al., 2012).
• Bifurcation: Local bifurcations of nonhomogeneous minimizers can occur from the uniform state as the interaction parameter increases, dictated by the instability of Fourier modes of the convolution kernel (Carrillo et al., 2019).
• Phase transitions: Both continuous (supercritical) and discontinuous (subcritical) transitions are classified, with the order determined by the sign of higher-order expansions (e.g., cubic term in Crandall–Rabinowitz theory) and the diffusion exponent $\Omega\subset\mathbb{R}^d$3 (Carrillo et al., 2019).
• Mesa limit: As $\Omega\subset\mathbb{R}^d$4, solutions converge towards hard height-constrained states (patches of density bounded by one), and the free-energy functionals $\Omega\subset\mathbb{R}^d$5-converge accordingly (Carrillo et al., 2019).

For singular kernels (e.g., Riesz, Newtonian), precise thresholds for global existence versus finite-time blow-up have been derived, involving sharp constants from the Hardy–Littlewood–Sobolev inequality and optimal interpolation arguments (Zhou et al., 21 Dec 2025). For $\Omega\subset\mathbb{R}^d$6 in the intermediate (critical) regime, equivalence of $\Omega\subset\mathbb{R}^d$7 and $\Omega\subset\mathbb{R}^d$8 criteria for blow-up/global existence is established (Zhou et al., 21 Dec 2025).

4. Gradient-Flow Structure and Long-Time Asymptotics

The degenerate diffusion-aggregation equation is a 2-Wasserstein gradient flow of the associated free energy on the space of probability measures with finite second moment (Carrillo et al., 2018, Matthes et al., 2024). This structure underpins:

• Energy dissipation identities: $\Omega\subset\mathbb{R}^d$9 for all solutions,
• Contractivity and exponential convergence towards equilibrium under $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$0-convexity or approximate convexity of the energy functional,
• Sharp decay rates in $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$1 and Wasserstein distance, controlled by the convex part of the energy (Matthes et al., 2024, Chayes et al., 2012).

For systems with small cross-diffusion, exponential convergence to a unique (deformed) steady state is preserved for small coupling, at a rate close to, but less than, that without cross-diffusion (Matthes et al., 2024). In particular, compact support of equilibrium is maintained, and degeneracy at zero is handled by technical cut-offs, moment control, and boundary-vanishing chemical potentials.

Self-similar asymptotics and the emergence of compactly supported steady states are rigorously proved in models with both classical and fractional diffusion, repulsive-attractive or Morse-type kernels, conditional on the regimes ($$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$2, $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$3, spatial dimension) (Carrillo et al., 2015, Carrillo et al., 2018).

5. Numerical Schemes and Computational Considerations

Robust fully discrete schemes for these equations must maintain both positivity of the density and monotonic decay of a discrete free energy. Finite volume methods, implicit in time, are constructed using dimensional splitting and flux reconstructions that preserve positivity (via M-matrix structure or suitable CFL conditions) and guarantee energy dissipation at each step (Bailo et al., 2018). These schemes are developed for general inhomogeneous diffusion terms ($$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$4 possibly degenerating at $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$5) as well as nonlocal interactions (with possible singularity at the origin).

For high-dimensional problems, sweeping splitting methods decouple the updates in each coordinate, enabling efficient and scalable computation without sacrificing the variational structure. These methods accurately capture sharp features, including free boundaries, metastability, merging, and phase transitions.

6. Extensions: Anisotropic Diffusion, Systems, and Open Problems

Anisotropic and degenerate diffusion tensors introduce additional complexities: when the degeneracy set is "thin" in fractal dimension (box-dimension less than $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$6 for $$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$7 data), global existence of very weak solutions persists, even if diffusion vanishes on a fractal set (Eckardt et al., 2024). Nonlocal adhesion operators can control the evolution even in regions where diffusion is absent, provided the aggregation term is sufficiently regular and the reaction term is present.

Cross-diffusion systems with energy non-convexities (e.g., two-component aggregation–diffusion–cross-diffusion PDEs) are analyzed by splitting the energy into convex and small non-convex pieces, with control estimates ensuring dominance of convexity and exponential decay of energy (Matthes et al., 2024).

Among the open challenges and ongoing directions:

• Fine regularity and free boundary regularity for weak solutions, especially in the presence of degenerate diffusion and singular kernels,
• Quantitative rates and uniform regularity in large domains/$$\partial_t\rho = \nabla\cdot\left(\nabla\rho^m - \rho\,\nabla (K*\rho)\right), \qquad \rho(\cdot,0) = \rho_0 \ge 0,$$8,
• Characterization of global minimizers and phase diagrams in regimes with strong aggregation or singular interactions,
• Stability and attractor properties in multi-species and anisotropic variants,
• Generalization of sharp blow-up criteria and threshold phenomena to more complex nonlinearities and coupling structures (Zhou et al., 21 Dec 2025, Eckardt et al., 2024, Carrillo et al., 2019).

These results jointly establish a comprehensive theoretical and computational framework for understanding the multi-dimensional nonlocal degenerate diffusion-aggregation equation and its equilibrium and dynamical phenomenology.