Diffusive Instabilities in Dusty Disks: Linear Growth and Nonlinear Breakdown

Abstract

We revisit the diffusive instability in dusty disks that arises when the dust mass diffusivity and/or viscosity decreases sufficiently steeply with increasing dust density. Our updated model includes an incompressible, viscous gas that responds azimuthally and couples to the dust through drag. We show that the basic criterion for diffusion-slope-driven instability remains approximately $β_\mathrm{diff}\lesssim -2$ for small dust stopping times, with gas feedback providing only modest quantitative changes for parameters motivated by streaming-instability turbulence. We perform nonlinear numerical calculations and confirm linear growth and mode selection toward the fastest-growing wavenumber. However, for power-law closures $D\proptoΣ^{β_\mathrm{diff}}$ with $β_\mathrm{diff}<0$, the nonlinear evolution does not saturate. Instead, steepening gradients amplify the nonlinear dust-pressure term and drive finite-time collapse into increasingly sharp spikes. Motivated by the absence of multidimensional saturation channels in our 1D framework, we test a simple piecewise closure in which the negative diffusion slope operates only over a finite density interval. This modification eliminates blowup and produces peak densities controlled by the imposed saturation scale. Our results support diffusive instabilities as a linear organizing mechanism in dusty turbulence, while highlighting that realistic nonlinear saturation requires additional physics beyond the present closure.

Figures (7)

• Figure 1: Growth rates $\gamma_\mathrm{max}$ (top row) for the fastest growing modes $K_\mathrm{max}$ (bottom row) of the diffusive instability in $\beta_\mathrm{diff}$ - $\tau_\mathrm{s}$ - space for different values of $\mathrm{Sc}$ and $\mathrm{Sc}_\mathrm{g}$ in different columns. Parameters kept constant between the panels are $\delta = 10^{-5}$, $\alpha_\mathrm{g} = 10^{-3}$, $\varepsilon = 1$. All three 'types' of instabilities discussed in Gerbig_Lin_Lehmann_2024 are present in the depicted parameter space: the diffusive instability driven by the diffusion slope is permitted for sufficiently negative $\beta_\mathrm{diff}$ and sufficiently small $\tau_\mathrm{s}$ (bottom left region of all panels), the diffusive instability driven by the viscosity slope is present for a sufficiently negative viscosity slope $\beta_\mathrm{visc}$ and large $\tau_\mathrm{s}$ (bottom right of second column, and right side of fourth column), and the diffusive overstability driven by the viscosity slope is present for sufficiently positive $\beta_\mathrm{visc}$ and large $\tau_\mathrm{s}$ (top right of second column, and right side of third column). Letters in the second column correspond to setups in the numerical experiments in Sect. \ref{['sect:numericaltests']}.
• Figure 2: Linear-nonlinear consistency test for different branches of the diffusive instability. An exact linear eigenmode solution for $K=100$ is evolved forward in time using the pseudo-spectral IMEX solver. Panels show the time evolution of the amplitude of the dominant Fourier component of the normalized surface density perturbation (red), compared to the prediction from the linear eigenvalue problem (black). Panels (a) - (d) correspond, respectively, to: (a) a stable, exponentially damped mode ($\tau_\mathrm{s} = 0.1, \beta_\mathrm{visc} = \beta_\mathrm{diff} = -0.5$); (b) the diffusion-slope-driven instability ($\tau_\mathrm{s} = 0.1, \beta_\mathrm{visc} = \beta_\mathrm{diff} = -2.5$); (c) the (dust) viscosity-slope-driven instability ($\tau_\mathrm{s} = 10, \beta_\mathrm{visc} = \beta_\mathrm{diff} = -1.3$); and (d) the overstable mode driven by the (dust) viscosity-slope ($\tau_\mathrm{s} = 50, \beta_\mathrm{visc} = \beta_\mathrm{diff} = 0.7$). The setups are marked in the instability map in Fig. \ref{['fig:instability_map']}. The insets in panels (a)-(c) show how the instantaneous numerical growth rate converge to the analytic eigenvalue. The inset in panel (d) shows the spiralling trajectory of the overstable mode in the $(\Re(\hat{v}_x), \Re(\hat{v}_y))$ phase space.
• Figure 3: Comparison between dust-only models, testing different values for dust-to-gas ratio $\varepsilon$ and gas viscosity $\alpha_\mathrm{g}$. All setups share $\delta = \alpha = 10^{-5}$, $\tau_\mathrm{s} = 0.1, \beta_\mathrm{visc} = \beta_\mathrm{diff} = -2.5$ such that a $K=100$ mode is subject to diffusive instability driven by diffusion slope. Left: Linear growth rates from the eigenvalue problem. Right: Evolution of a $K=100$ eigenmode with amplitudes that remain in the linear regime.
• Figure 4: Nonlinear evolution and blow-up for the diffusive instability. Evolution of the maximum dust surface density for resolutions 32, 64, and 128 using both the native IMEX solver and Dedalus (left), and snapshots of $\Sigma(x,t)$ for the $N_x = 128$ runs at selected times for the native solver (upper right), and Dedalus (lower right). The initial $k=100$ eigenmode grows slightly below the linear rate (black dotted line), then transitions into a regime of super-linear growth once nonlinear gradients dominate. The system does not saturate as steepening gradients reinforce the nonlinear dust-pressure term, while the viscosity is unable to regularize. Both solvers show the same qualitative evolution, where the initially sinusoidal mode sharpens into narrow spikes with width approaching the grid scale. Minor differences in peak sharpness reflect the different high-$k$ regularization properties of the two methods.
• Figure 5: Evolution of a noise initial condition. Maximum dust density as a function of time (left) and snapshots of $\Sigma(x,t)$ for various times. Independently of resolution, the initial noise field organizes into a close-to pure Fourier mode with a wavenumber approaching the fastest-growing $k$ from linear theory.