---
title: Inviscid Incompressible Porous Medium Equations
url: https://www.emergentmind.com/topics/inviscid-incompressible-porous-medium-equations-ipm
type: topic
---

# Inviscid Incompressible Porous Medium Equations

The inviscid incompressible porous medium equations (IPM) are a class of active-scalar systems in which a transported scalar—typically denoted \(\rho\) or \(\theta\)—determines an incompressible Darcy velocity through a nonlocal elliptic law. In two dimensions, a standard normalization is
\[
\partial_t \rho + u\cdot \nabla \rho = 0,\qquad u = -\nabla p -(0,\rho),\qquad \nabla\cdot u=0,
\]
while equivalent sign conventions and rotated formulations, such as \(u+\nabla p=(0,\rho)\) or \(u+\nabla p=(-\rho,0)\), also occur in the literature. The subject includes Eulerian and Lagrangian formulations, stratified steady states, asymptotic stability and instability mechanisms, convex-integration weak solutions, and several singularity scenarios ranging from reduced one-dimensional models to finite-time blow-up in wedge geometries [1403.5749], [1411.6958], [1102.2597], [2410.01297], [2511.01827].

## 1. Governing equations and equivalent formulations

In the whole-plane formulation used for the 2D IPM equation, the scalar satisfies a pure transport law and the velocity is determined by Darcy’s law plus incompressibility:
\[
\partial_t \rho + (u\cdot\nabla)\rho = 0,\qquad u = -\nabla \pi -(0,\rho),\qquad \nabla\cdot u = 0.
\]
Taking divergence gives
\[
-\Delta \pi = \partial_{x_2}\rho,
\]
so the velocity becomes a zero-order singular integral of \(\rho\). One convenient representation is
\[
u=\nabla^\perp(-\Delta)^{-1}R_1\rho,
\]
with \(\nabla^\perp=(-\partial_{x_2},\partial_{x_1})\) and \(R_1\) the first Riesz transform [2410.23727].

The same equation is often written with a transported scalar \(\theta\) and constitutive relation
\[
\partial_t \theta + (u\cdot\nabla)\theta = 0,\qquad u=-\nabla p-(0,\theta),
\]
posed on \(\mathbb R^2\). In that normalization, the vorticity satisfies
\[
\omega=\nabla^\perp\cdot u = -\partial_{x_1}\theta,
\]
which makes IPM structurally comparable to 2D Euler and SQG, but with a scalar-to-velocity law tied to Darcy forcing rather than a stream-function law [1403.5749].

The same active-scalar structure is studied on \(\mathbb T^2\), on the strip \(\mathbb T\times[-1,1]\) with no-penetration, and on wedge domains. In the wedge setting one frequently uses the rotated normalization
\[
u+\nabla p = (-\rho,0),\qquad \partial_t\rho + u\cdot\nabla \rho =0,\qquad \nabla\cdot u=0,
\]
which is equivalent up to orientation of gravity. The stream-function form then becomes
\[
u=\nabla^\perp \Psi,\qquad \Delta \Psi = \partial_{x_2}\rho.
\]
This suggests that much of the IPM literature is invariant under simple sign and rotation conventions, while the analytic questions depend more strongly on geometry, regularity class, and the structure of the scalar field than on a particular normalization [2511.01827].

## 2. Lagrangian formulation and analytic trajectories

For classical IPM solutions, the Lagrangian flow map \(X(a,t)\) is defined by
\[
\frac{d}{dt}X(a,t)=u(X(a,t),t),\qquad X(a,0)=a.
\]
Because \(\theta\) is transported, one has
\[
\theta(X(a,t),t)=\theta_0(a).
\]
Using the vorticity identity \(\omega=-\partial_{x_1}\theta\), the 2D Biot–Savart law, and incompressibility of the flow map, IPM admits a closed Lagrangian singular-integral system:
\[
\frac{dX(a,t)}{dt}
=
\frac{1}{2\pi}\int_{\mathbb R^2}
\frac{(X(a,t)-X(b,t))^\perp}{|X(a,t)-X(b,t)|^2}\,
\{\theta_0(b),X_2(b,t)\}\,db,
\]
where \(\{\cdot,\cdot\}\) is the Poisson bracket. In this representation the source term is a conserved Lagrangian quantity built from the initial scalar and the second component of the flow map [1403.5749].

The same work derives an evolution equation for \(\nabla_a X\) involving a Calderón–Zygmund matrix kernel
\[
K(y)=\frac{1}{2\pi |y|^4}
\begin{pmatrix}
2y_1y_2 & y_2^2-y_1^2\\
y_2^2-y_1^2 & -2y_1y_2
\end{pmatrix},
\]
analytic away from the origin and of mean zero on spheres. Those structural properties are the key input in the time-analyticity argument.

Under the local well-posedness hypothesis recalled there—specifically \(\theta_0\in C^{1,\gamma}(\mathbb R^2)\) with \(\gamma\in(0,1)\), so that \(u(\cdot,t)\in C^{1,\gamma}\)—the Lagrangian trajectories are real analytic in time as long as the chord-arc parameter
\[
\Lambda(t)=\exp\Big(\int_0^t \|\nabla u(s)\|_{L^\infty}\,ds\Big)
\]
remains finite. Thus, for 2D IPM in the classical regime, Eulerian regularity at the level of Hölder-continuous velocity gradient implies real-analytic temporal regularity of particle paths [1403.5749].

## 3. Stratified equilibria and asymptotic stability

A basic structural fact is that sufficiently regular stationary IPM solutions in a bounded simply connected domain with \(u\cdot n=0\) are necessarily stratified:
\[
\rho(x,y)=f(y),\qquad u=0.
\]
The proof uses the stationary transport equation \(u\cdot\nabla\rho=0\), incompressibility, and the Darcy law to show first that
\[
\int_\Omega u_2\rho = 0,
\]
and then that the same quantity equals \(\int_\Omega |u|^2\), forcing \(u\equiv 0\). This identifies vertical stratification as the canonical steady-state class for IPM [1411.6958].

Around the linear stratification \(\bar\rho(y)=y\), the perturbation \(\theta=\rho-\bar\rho\) satisfies a partially damped dynamics. On \(\mathbb R^2\), one obtains the relaxation identity
\[
\frac{d}{dt}\int_{\mathbb R^2} (\rho-y)^2\,dx\,dy
=
-2\int_{\mathbb R^2}|u|^2\,dx\,dy,
\]
which shows that the deviation from the background profile is a Lyapunov functional and that the velocity is dissipated by stratification rather than by viscosity. The same paper proves global asymptotic stability for sufficiently small perturbations of \(\bar\rho(y)=y\): on \(\mathbb R^2\), small data in \(W^{4,1}\cap H^s\) with \(s>20\) generate global solutions converging back to \(\bar\rho\), while on \(\mathbb T^2\) small data in \(H^s\), \(s\ge 20\), converge to a stratified stationary limit, because purely vertical modes are undamped and persist [1411.6958].

A quantitative refinement is available for the three fundamental domains \(\mathbb R^2\), \(\mathbb T^2\), and \(\mathbb T\times[-1,1]\). For the linear stratification \(Nx_2\), and on \(\mathbb R^2\) also for quasi-linearly stratified backgrounds
\[
\rho_s(x_2)=N x_2 + \int_0^{x_2}\sigma(\eta)\,d\eta,
\]
global asymptotic stability holds provided the buoyancy frequency \(N\) dominates the perturbation size and, in the quasi-linear case, the weighted \(L^1\) size of \(\hat\sigma\). For \(m>3\), the whole-space result gives
\[
\|u(t)\|_{H^m}\le \frac{C}{(1+Nt)^2}\|\theta_0\|_{H^m},\qquad
\|\nabla u_2(t)\|_{H^{m-1}}\le \frac{C}{1+Nt}\|\theta_0\|_{H^m},
\]
and \(\|\theta(t)\|_{L^2}\to 0\) as \(Nt\to\infty\). On \(\mathbb T^2\) and \(\mathbb T\times[-1,1]\), the solution converges not to \(Nx_2\) itself but to an asymptotic vertical profile \(Nx_2+\theta_\infty(x_2)\), with sharp decay rates for the decaying \(x_1\)-dependent modes. The argument uses an anisotropic commutator estimate that lowers the regularity threshold to any real \(m>3\) in the purely linear case [2210.11437].

## 4. Weak solutions, relaxation, and Muskat mixing

For weak-solution theory, IPM is often rewritten as a differential inclusion. In the constant-viscosity case one introduces
\[
u=2v+(0,\rho)
\]
and a flux variable \(m\), so that the system becomes a linear conservation law
\[
\partial_t\rho+\operatorname{div}m=0,\qquad
\operatorname{div}(u-(0,\rho))=0,\qquad
\operatorname{curl}(u+(0,\rho))=0,
\]
together with the pointwise constitutive constraint
\[
|\rho|=1,\qquad m=\tfrac12 \rho u.
\]
The corresponding wave cone is
\[
\Lambda=\{(\rho,u,m): |\rho|^2=|u|^2\},
\]
and the exact relaxation of the constitutive set is
\[
K^\Lambda=
\left\{
(\rho,u,m):
|\rho|\le 1,\ 
\left|m-\tfrac12\rho u\right|\le \tfrac12(1-\rho^2)
\right\}.
\]
This explicit relaxed set replaces earlier \(T4\)-configuration arguments and is the basis of a convex-integration construction of infinitely many weak solutions, including nontrivial weak solutions with compact support in time and mixing solutions for the unstable Muskat problem with flat interface. In that setting, the coarse-grained density selected by Otto’s relaxed formulation corresponds to the maximally mixing subsolution [1102.2597].

The same h-principle extends to viscosity jump \(|A_\mu|<1\). In that case weak solutions in \(C_tL^\infty_{w^*}\) are again recovered by convex integration once a subsolution is available; nontrivial weak solutions with compact support in time and mixing solutions to the unstable Muskat problem with initial flat interface persist. The viscosity contrast changes the geometry of the relaxation set: the paper identifies a pinch singularity preventing the two fluids from mixing wherever there is neither Rayleigh–Taylor nor vorticity at the interface. It also verifies that the connection between subsolutions and Otto’s Lagrangian relaxed solution, previously established for \(A_\mu=0\), remains valid for \(|A_\mu|<1\) [2004.03307].

This weak-solution theory is conceptually distinct from homogenization of ideal flow through a perforated solid region. For 2D Euler in a perforated domain, the macroscopic limit is either the full-plane Euler flow or Euler in the exterior of an effective impermeable obstacle, depending on the scaling of hole size and distance; no Darcy-type or classical IPM law appears in the limit. That distinction separates active-scalar IPM from inviscid homogenization through porous geometries [1407.2792].

## 5. Critical regularity, instability, and small-scale formation

A different branch of the literature studies how much of the stratified stability picture survives at critical or near-critical regularity. One mechanism is long-time small-scale formation. For smooth IPM solutions on \(\mathbb R^2\), \(\mathbb T^2\), and the strip, the potential energy
\[
E(t)=\int_\Omega x_2\rho(x,t)\,dx
\]
satisfies
\[
E'(t)=-\|\partial_{x_1}\rho(t)\|_{\dot H^{-1}(\Omega)}^2.
\]
This monotonicity implies \(\int_0^\infty \|\partial_{x_1}\rho(t)\|_{\dot H^{-1}}^2\,dt<\infty\). For several geometric classes of initial data—odd-in-\(x_2\) data in \(\mathbb R^2\), symmetric torus data, “bubble” data, and “layered” rearrangements of stratified states—the same paper shows that any global smooth solution must develop unbounded Sobolev norms as \(t\to\infty\). These constructions yield nonlinear instability of a large class of stratified steady states in the sense of infinite-time derivative growth, even when global smoothness is assumed [2102.05213].

At the critical Lipschitz level, IPM is mildly ill-posed in \(W^{1,\infty}\) near arbitrary vertically stratified backgrounds \(g(x_2)\). Writing \(\rho=g(x_2)+\eta\), the perturbation equation contains the variable-coefficient term \(-g'(x_2)R_2\eta\). For every sufficiently small \(\varepsilon>0\), and for any \(g\in \dot B_{p,1}^{2+2/p}\) with \(2<p<\infty\) and \(g'(x_2)\neq 0\), there exists initial data \(\rho_0\) with
\[
\|\rho_0-g(x_2)\|_{W^{1,\infty}}\le \varepsilon
\]
such that the local solution satisfies
\[
\sup_{0\le t\le c\varepsilon}\|\rho(\cdot,t)-g(x_2)\|_{W^{1,\infty}}\ge c,
\]
for universal \(c>0\). The construction uses smooth compactly supported perturbations \(f_N\) for which the Riesz term \(g'R_2\nabla^\perp f_N\) is large while the \(W^{1,\infty}\) norm is small, producing norm inflation in arbitrarily short time. Notably, this applies even when \(g'(x_2)<0\), a regime usually regarded as physically stable [2410.23727].

At the \(L^2\)-critical Sobolev threshold \(H^2\), the stable equation near \(-x_2\) is strongly ill-posed. For arbitrarily small \(\varepsilon>0\), there exist perturbations \(\rho_{\mathrm{in}}\) with
\[
\|\rho_{\mathrm{in}}\|_{H^2(\mathbb R^2)}\le \varepsilon
\]
such that the corresponding classical solution of IPM with initial data \(-x_2+\rho_{\mathrm{in}}\) satisfies
\[
\|\rho(t)+x_2\|_{H^2(\mathbb R^2)}=\infty
\quad\text{for every }t\in(0,T_\varepsilon],
\]
while remaining bounded in \(H^{2-\varepsilon'}\) for every \(\varepsilon'>0\). The mechanism combines a small \(H^2\) perturbation that cancels the background profile near the origin with a multiscale construction generating strong hyperbolic deformation and explosive growth of second derivatives [2410.01297].

## 6. Reduced models, blow-up mechanisms, and related interpretations

Several recent results isolate singularity mechanisms in symmetry classes or reduced models. One infinite-energy class of Castro–Córdoba–Gancedo–Orive type solutions uses the ansatz
\[
\psi(\tau,x,y)=y f(\tau,x),\qquad \rho(\tau,x,y)=y f(\tau,x)+yG(\tau),
\]
which reduces 2D IPM to a one-dimensional nonlocal equation for \(b=\partial_x f\):
\[
\partial_\tau b+\left(\int_{-\pi}^x b(\tau,\bar x)\,d\bar x\right)\partial_x b-b^2+\frac1\pi\int_{-\pi}^{\pi}b^2\,dx
=
\frac1\pi\left(\int_0^\tau\int_{-\pi}^{\pi} b^2\,dx\,d\bar\tau\right)b.
\]
This reduced equation has the explicit self-similar blow-up solution
\[
b(\tau,x)=\frac{\mu}{\cos(\mu\tau)}\cos x,
\]
blowing up at \(\tau^*=\pi/(2\mu)\). The stability theory identifies a sharp regularity threshold: small \(C^3\) perturbations preserve finite-time blow-up with profile \((\tau^*-\tau)^{-1}\cos x\), whereas arbitrarily small \(C^{2-\epsilon}\) perturbations can destroy asymptotic self-similarity. Via a nonlinear change of variables, this becomes asymptotic stability of the steady family \(\mu\cos x\) for the Proudman–Johnson equation [2507.17381].

Finite-time singularity formation has also been established for the full 2D IPM in wedge domains, without boundary mass. Using 1-homogeneous solutions
\[
\rho(t,r,\theta)=rP(t,\theta),\qquad \Psi(t,r,\theta)=r^2G(t,\theta),
\]
the equation reduces to a one-dimensional angular system
\[
\partial_t P + 2G P' = G'P,\qquad
G''+4G = P\sin\theta + P'\cos\theta,
\]
with \(G(t,0)=G(t,L)=0\). A self-similar ansatz \(P(t,\theta)=(1-t)^{-1}P_*(\theta)\), \(G(t,\theta)=(1-t)^{-1}G_*(\theta)\) yields finite-time gradient blow-up at \(t=1\). The resulting 2D solutions are Lipschitz initially, can be chosen compactly supported, vanish on the boundary, are smooth away from the origin, and satisfy
\[
\limsup_{t\to 1}\int_0^t \|\nabla \rho(s)\|_{L^\infty}\,ds = +\infty.
\]
The vanishing of density on the boundary is a defining feature: the blow-up mechanism overcomes the full regularizing effect of transport rather than relying on boundary-pinned mass [2511.01827].

A different reduced scenario starts from the 2D periodic half-plane and extracts a boundary-layer model for the boundary trace. The resulting one-dimensional periodic equation is
\[
\partial_t p + u\,\partial_x p = 0,\qquad
u=gH_a\partial_x p,
\]
where
\[
H_af(x)=\operatorname{P.V.}\int_{\mathbb T} f(x-y)\,\frac{1}{\pi}\frac{y}{y^2+a^2}\,dy.
\]
This model is locally well-posed for smooth periodic data and blows up in finite time for smooth, bounded, nonnegative, even data with \(p_0(0)=0\) and \(p_0'\ge 0\) on \([0,\pi)\). The proof uses a Beale–Kato–Majda-type criterion together with a weighted integral inequality for \(H_a\), and positions the model as a boundary-layer analogue of the Córdoba–Córdoba–Fontelos equation [2412.16376].

Broader porous-medium terminology should be handled carefully. A geometric variational theory for an incompressible fluid moving through an elastic porous medium derives coupled fluid–solid equations and recovers Biot’s wave equations in suitable parameter regimes; that framework is a poromechanical generalization rather than the active-scalar IPM equation. Conversely, homogenization of ideal Euler flow through perforated media yields transparency or impermeability, not an effective Darcy/IPM law. This suggests that “porous medium” in PDE usage covers several mathematically distinct limits, of which active-scalar IPM is only one [2007.02605], [1407.2792].

Source: https://www.emergentmind.com/topics/inviscid-incompressible-porous-medium-equations-ipm