---
title: 'Leibenson Equation: Nonlinear Diffusion'
url: https://www.emergentmind.com/topics/leibenson-equation
type: topic
---

# Leibenson Equation: Nonlinear Diffusion

The Leibenson equation, also spelled Leibenzon equation, is the doubly nonlinear parabolic equation
$$
\partial_t u=\Delta_p(u^q)=\operatorname{div}\big(|\nabla(u^q)|^{p-2}\nabla(u^q)\big),
$$
posed on Euclidean spaces or, more generally, on Riemannian manifolds, with parameters \(p>1\) and \(q>0\). It simultaneously generalizes the heat equation, the porous medium and fast diffusion equations, and the parabolic \(p\)-Laplacian. In the recent geometric PDE literature it is studied as a nonlinear diffusion model for filtration in porous media, turbulent flow, and non-Newtonian transport, with existence, propagation, decay, gradient, and probabilistic results now available on broad classes of manifolds [2601.20640] [2604.23227].

## 1. Definition, operator structure, and special cases

On a Riemannian manifold \((M,g)\), the \(p\)-Laplacian of a function \(v\) is
$$
\Delta_p v=\operatorname{div}\big(|\nabla v|^{p-2}\nabla v\big),
$$
where \(\nabla\) and \(\operatorname{div}\) are defined with respect to the metric \(g\). The Leibenson equation therefore diffuses the nonlinear quantity \(u^q\) through a nonlinear spatial operator. In abstract evolution form it can be written as
$$
\partial_t u + A(\Phi(u))=0,\qquad \Phi(u)=u^q,\quad A(w)=-\Delta_p w,
$$
which is why the equation is termed doubly nonlinear: the state dependence enters through \(u^q\), while the flux is nonlinear in the gradient through \(\Delta_p\) [2601.20640].

Several standard equations appear as parameter reductions of the same model [2604.23227].

| Parameters | Equation | Standard name |
|---|---|---|
| \(p=2,\ q=1\) | \(\partial_t u=\Delta u\) | Heat equation |
| \(p=2,\ q>1\) | \(\partial_t u=\Delta(u^q)\) | Porous medium equation |
| \(p=2,\ 0<q<1\) | \(\partial_t u=\Delta(u^q)\) | Fast diffusion equation |
| \(q=1\) | \(\partial_t u=\Delta_p u\) | Parabolic \(p\)-Laplace equation |

The PDE literature also places the equation in the context of filtration of turbulent compressible fluids through porous media, and the name “Leibenson” is attached to that hydrodynamic lineage [2603.27791]. A plausible implication is that the equation occupies the same structural role for doubly nonlinear diffusion that the heat equation occupies for linear diffusion.

## 2. Weak formulation and analytical framework

Because the operator may be degenerate or singular, the natural solution concept is weak. For an interval \((0,T)\) and domain \(\Omega\subset M\), a typical weak solution framework requires
$$
u\in C\big([0,T];L^{q+1}(\Omega)\big),\qquad
u^q\in L^p\big((0,T);W^{1,p}_0(\Omega)\big),
$$
or, in the global manifold setting, the corresponding spaces on \(M\). The weak formulation is
$$
\Big[\int_M u\,\varphi\Big]_{t_1}^{t_2}
+\int_{t_1}^{t_2}\!\!\int_M
\Big(-u\,\partial_t\varphi
+|\nabla u^q|^{p-2}\langle\nabla u^q,\nabla\varphi\rangle\Big)\,d\mu\,dt=0,
$$
for test functions \(\varphi\) in the natural Bochner–Sobolev class; subsolutions and supersolutions are defined by replacing equality with the corresponding inequality [2601.20640].

This formulation is compatible with local Dirichlet problems on precompact domains and with the Cauchy problem on the whole manifold. It also aligns with the regularity class used in the sharp upper-bound theory, where one imposes
$$
u\in C\big(I;L^{1+q}(M)\big),\qquad
u^q\in L_{\mathrm{loc}}^p\big(I;W^{1,p}(M)\big),
$$
together with admissible time-dependent test functions in
$$
W_{\mathrm{loc}}^{1,\,1+\frac1q}\big(I;L^{1+\frac1q}(M)\big)\cap
L_{\mathrm{loc}}^p\big(I;W^{1,p}(M)\big)
$$
[2604.23227]. For weak subsolutions on open subsets, an equivalent formulation is used in the sharp long-distance estimates, again based on \(u\in C(I;L^{1+q})\) and \(u^q\in L^p_{\mathrm{loc}}(I;W^{1,p})\) [2603.27791].

A recurrent feature of this framework is compatibility with comparison arguments. In the manifold existence theory, if \(u\) is a weak supersolution and \(v\) a weak subsolution with ordered initial data, then \(u\ge v\) almost everywhere, and the positive part \((v-u)_+\) has decaying \(L^2\)-norm [2601.20640].

## 3. Existence theory on Riemannian manifolds

A central recent result establishes global existence on arbitrary geodesically complete Riemannian manifolds under the conditions
$$
p>1,\qquad q>0,\qquad pq\ge 1,
$$
for nonnegative initial data
$$
u_0\in L^1(M)\cap L^\infty(M).
$$
Under these assumptions, the Cauchy problem
$$
\partial_t u=\Delta_p(u^q),\qquad u(\cdot,0)=u_0,
$$
admits a nonnegative bounded weak solution on \(M\times(0,\infty)\) [2601.20640].

The construction is local-to-global. First, one solves a truncated Dirichlet problem on a bounded domain \(\Omega\subset M\) and finite time interval, replacing \(u\) by
$$
\tilde u=\min\Big(N,\max\Big(u,\frac1N\Big)\Big)
$$
and using the regularized flux
$$
A(u,\nabla u)=q^{p-1}\tilde u^{(q-1)(p-1)}|\nabla u|^{p-2}\nabla u.
$$
A Galerkin approximation is then built in a finite-dimensional basis of \(L^2(\Omega)\cap W_0^{1,p}(\Omega)\). The key analytic input is the standard monotonicity estimate for the \(p\)-Laplacian, together with energy bounds that control the approximate solutions in \(L_t^\infty L_x^2\) and \(L_t^pW_x^{1,p}\). Uniform lower and upper bounds show that the truncation becomes inactive, after which compactness and strong convergence identify a weak solution of the original equation on \(\Omega\times(0,T)\). Finally, an exhaustion by increasing geodesic balls \(B_k\nearrow M\) yields a global solution on the whole manifold [2601.20640].

One notable aspect of this existence theorem is geometric minimality. The hypothesis is only geodesic completeness: no curvature bounds, no global volume-growth assumptions, and no global Sobolev or Poincaré inequalities are required beyond the local ones available on precompact sets [2601.20640]. This sharply contrasts with earlier manifold results that required Cartan–Hadamard geometry, nonnegative Ricci curvature, or relative Faber–Krahn inequalities.

The same work records several qualitative consequences of the construction. For nonnegative initial data, solutions remain nonnegative, and for every \(r\ge1\),
$$
\|u(\cdot,t)\|_{L^r(M)}\le \|u_0\|_{L^r(M)}.
$$
Thus the mass and all \(L^r\)-norms are non-increasing in time [2601.20640].

## 4. Diffusion regimes and sharp quantitative behavior

The basic regime parameter is
$$
D:=1-q(p-1).
$$
Its sign separates slow and fast diffusion phenomena. When \(q(p-1)>1\), equivalently \(D<0\), the equation is in the slow diffusion range. When \(q(p-1)<1\), equivalently \(D>0\), it is in the fast diffusion range [2506.07221].

In the slow diffusion case, compactly supported data exhibit finite propagation speed. Recent manifold work states that under
$$
q(p-1)>1,
$$
solutions with compactly supported initial data have support expanding at a finite rate, and the associated mean value inequality is a main tool in the propagation analysis [2601.20640]. Gradient estimates in this regime are obtained after introducing the pressure-type variable
$$
v=\frac{q(p-1)}{\delta}u^{\frac{\delta}{p-1}},\qquad \delta=q(p-1)-1>0,
$$
and controlling the Li–Yau-type quantity
$$
\frac{|\nabla v|^p}{v}-\frac{\partial_t v}{v}
$$
under a Ricci lower bound and two-sided bounds on \(|\nabla v|^{p-2}v\) [2506.07221].

In the fast diffusion case, \(D>0\) governs both decay and extinction. For the weighted equation
$$
\rho\,\partial_t u=\Delta_p u^q,
$$
posed on a weighted manifold \((M,\mu)\), finite extinction is proved when \(D>0\), a weighted Sobolev inequality holds, and
$$
\left\|\frac{\rho}{\omega}\right\|_{L^\theta(M,d\mu)}<\infty
$$
for the exponent
$$
\theta=\frac{\kappa}{\kappa-1-\frac{D}{\sigma}}.
$$
Under these hypotheses, bounded nonnegative weak subsolutions vanish identically after a finite time [2412.06496].

Sharp global upper bounds are also available in the fast diffusion regime. On manifolds with non-negative Ricci curvature, or more generally under a relative Faber–Krahn inequality plus polynomial volume growth, one has sharp long-distance estimates for weak subsolutions in the whole range
$$
p>1,\qquad q>0,\qquad q(p-1)<1.
$$
The resulting off-diagonal bounds involve the intrinsic scaling variable \(|x|/t^{1/\beta}\), where
$$
\beta=p-D,
$$
and the spatial decay exponent \(p/D\); the work explicitly states that this sharpens earlier bounds by removing an extra logarithmic factor and proves a previous conjecture [2603.27791].

Long-time regularization is tied to Sobolev geometry. If \(1<p<n\), the Euclidean-type Sobolev inequality
$$
\left(\int_M |v|^{\frac{pn}{n-p}}\,d\mu\right)^{\frac{n-p}{n}}
\le C\int_M |\nabla v|^p\,d\mu
$$
is equivalent to a sharp \(L^\infty\) upper bound for bounded weak solutions. In particular, when
$$
p>nD=n\big(1-q(p-1)\big),
$$
one has
$$
\|u(t)\|_{L^\infty(M)}
\le C\,\|u_0\|_{L^1(M)}^{\frac{p}{p-nD}}\,t^{-\frac{n}{p-nD}},
$$
and conversely this type of upper bound forces the Sobolev inequality [2604.23227]. This establishes a nonlinear analogue of the classical heat-kernel/Sobolev equivalence.

## 5. Geometry, weighted settings, and reactive extensions

The manifold theory of the Leibenson equation is now spread across several geometric regimes. For bare existence, geodesic completeness suffices [2601.20640]. For sharp long-distance decay in the fast regime, the relevant analytic inputs are a relative Faber–Krahn inequality and a polynomial lower volume bound, with non-negative Ricci curvature singled out as a standard sufficient condition [2603.27791]. For long-time \(L^\infty\)-regularization, the decisive structure is a Euclidean-type Sobolev inequality [2604.23227]. For gradient estimates, one assumes a Ricci lower bound by a non-positive constant [2506.07221]. The literature therefore separates existence, regularization, and pointwise gradient control according to different geometric hypotheses.

Weighted manifolds introduce an additional layer of structure. In the finite-extinction theory, one works on a weighted manifold \((M,\mu)\) with
$$
d\mu=h(x)\,d\mathrm{vol}_g(x),
$$
and the ambient geometry is encoded through a weighted Sobolev inequality involving a weight \(\omega\). The integrability of \(\rho/\omega\) then determines whether the weighted Leibenson equation extinguishes in finite time [2412.06496]. This suggests that, in the fast diffusion range, geometry and inhomogeneity can be traded against each other through functional inequalities.

A distinct extension is the reaction–diffusion equation
$$
u_t=\Delta_p u^m+u^q
$$
on complete non-compact manifolds of infinite volume. Under
$$
1<p<N,\qquad m(p-1)\ge1,\qquad m>1,\qquad q>m(p-1),
$$
global weak solutions exist for sufficiently small initial data. If the manifold supports a Sobolev inequality, the result holds provided
$$
q>m(p-1)+\frac{p}{N},
$$
while under an additional Poincaré-type inequality it extends to the full interval \(q>m(p-1)\) [2505.08304]. The latter is described as having no Euclidean counterpart, because the non-compact infinite-volume geometry alters the balance between diffusion and reaction.

Several open directions are explicitly identified in this body of work. The condition \(pq\ge1\) in the general existence theorem is regarded as a minor technical restriction that might be removed; finer regularity and gradient bounds remain active questions; and long-time asymptotics on general manifolds, including Barenblatt-type behavior beyond Euclidean settings, are still largely open [2601.20640].

## 6. Probabilistic interpretation and nomenclature issues

A recent probabilistic development identifies the Leibenson equation on \(\mathbb{R}^d\) as a nonlinear Fokker–Planck equation and constructs its McKean–Vlasov counterpart. In that formulation the coefficients depend pointwise on the time-marginal density and on its first and second order derivatives. The Barenblatt solutions are realized as the one-dimensional marginal density curve of unique solutions to the associated McKean–Vlasov SDE, and these solutions form a nonlinear Markov process in the sense of McKean called the Leibenson process [2508.12979].

The same work emphasizes that the resulting SDE is highly singular: the diffusion is strongly degenerate and the drift is merely of bounded variation. Nevertheless, the solutions are shown to be probabilistically strong, meaning measurable functionals of the driving Brownian motion and the initial condition [2508.12979]. This provides a stochastic counterpart to the PDE that is structurally analogous, in a nonlinear setting, to the relation between Brownian motion and the heat equation. A plausible implication is that the Leibenson equation now belongs simultaneously to nonlinear diffusion theory and to nonlinear Markov-process theory.

A recurrent nomenclature issue is the confusion between Leibenson and Levinson. One arXiv paper with “Leibenson” in discussion is in fact about Levinson’s log–log theorem in complex analysis and linear elliptic PDE; it explicitly states that “Leibenson” there is a misspelled reference to Levinson’s theorem and not an equation [2012.07169]. In PDE and geometric-analysis usage, by contrast, the Leibenson equation denotes the doubly nonlinear evolution
$$
\partial_t u=\Delta_p(u^q).
$$

Source: https://www.emergentmind.com/topics/leibenson-equation