SIDM-Hydro Hybrid (SHH) Simulations

• SIDM-Hydro Hybrid (SHH) is a simulation framework that couples self-interacting dark matter dynamics with hydrodynamic galaxy formation physics using both Monte Carlo and fluid methods.
• The approach isolates the distinct influences of dark matter self-scattering and baryonic processes, enabling controlled comparisons across diverse halo mass scales.
• SHH techniques improve predictions of core formation, halo shapes, and thermalization by interpolating between kinetic and fluid regimes in cosmological settings.

Searching arXiv for the specified SHH-related papers and closely related SIDM+baryon simulation work. SIDM–Hydro Hybrid (SHH) denotes a class of numerical strategies that combine self-interacting dark matter (SIDM) dynamics with hydrodynamical galaxy-formation physics in a single simulation framework. In the literature summarized here, the term covers two closely related usages: first, cosmological hydrodynamical simulations in which CDM and SIDM are evolved with the same baryonic subgrid model so that the relative roles of self-scattering and baryon-driven structure formation can be isolated; second, a more formal hybrid method that interpolates between kinetic SIDM scattering and a fluid description through a continuous function of the local Knudsen number (Fry et al., 2015, Schon et al., 14 Aug 2025). Across dwarf, Milky Way-like, early-type-galaxy, massive-galaxy, and cluster scales, SHH methods are primarily used to determine when baryons erase, preserve, or invert the structural signatures that SIDM would produce in dark-matter-only calculations (Robles et al., 2017, Rose et al., 2022, Despali et al., 2022).

1. Definition and conceptual scope

In its broadest usage, SHH is an organizational label for simulations that re-run the same halos in CDM and SIDM, with and without baryons, under otherwise matched numerical conditions. This structure appears explicitly in several suites: early-type-galaxy zooms with CDM-h, CDM-d, SIDM-h, and SIDM-d realizations (Despali et al., 2022); dwarf-galaxy zooms comparing CDM-DMO, SIDM-DMO, CDM-Hydro, and SIDM-Hydro (Robles et al., 2017); and Milky Way-like or massive-galaxy simulations that contrast dark-matter-only and hydrodynamical variants under common galaxy-formation models (Rose et al., 2022, Mastromarino et al., 2022).

The central rationale is that SIDM modifies halo thermalization, core formation, and shape, while baryons alter the potential through cooling, star formation, feedback, and black-hole growth. SHH calculations therefore aim to separate three effects that are degenerate in total mass profiles: dark-matter self-scattering, baryonic contraction or expansion, and the numerical consequences of mock-observable pipelines. This is why several studies compare not only density profiles, but also projected ellipticities, X-ray isophotes, Einstein radii, and dark-matter fractions (Despali et al., 2022, Robertson et al., 2018, Mastromarino et al., 2022).

A narrower and more formal usage appears in the later “SHH method for SIDM,” where SHH is a continuous interpolation between two limits of the collisional Boltzmann equation: a stochastic, pairwise SIDM N-body treatment in dilute regions and an ideal-fluid Euler solver in frequent-scattering regions (Schon et al., 14 Aug 2025). This suggests that “SHH” has evolved from a descriptive label for SIDM-plus-hydrodynamics cosmological runs into a specific multiregime solver architecture.

2. Numerical formulations and coupling strategies

Most SHH implementations treat SIDM scatterings with a Monte Carlo pairwise algorithm. In the early-type-galaxy zooms of Despali et al., elastic collisions with constant transfer cross section $\sigma_T/m_\chi = 1\ \mathrm{cm}^2 \mathrm{g}^{-1}$ are implemented explicitly in the N-body solver via a nearest-neighbour Monte Carlo method, with interaction probability proportional to $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$ (Despali et al., 2022). In the dwarf-galaxy ChaNGa runs, the local scattering rate is written as $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$ and the corresponding characteristic time scale is $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$ (Fry et al., 2015). In the FIRE-2 dwarf simulations, the pairwise acceptance probability is written as $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$, after which accepted pairs scatter elastically and isotropically in the center-of-mass frame (Robles et al., 2017). Related Monte Carlo schemes also appear in SWIFT–EAGLE/TangoSIDM, modified GADGET-3 with EAGLE, and BAHAMAS-SIDM (Correa et al., 2024, Mastromarino et al., 2022, Robertson et al., 2018).

The hydrodynamical side varies by code family. AREPO with the IllustrisTNG model is used for early-type galaxies and MW-like systems, including gas cooling, star formation, stellar feedback, black-hole growth, AGN feedback, and ideal-MHD (Despali et al., 2022, Rose et al., 2022). ChaNGa with Gasoline-inherited modules is used in one dwarf SHH suite, with metal-line cooling, a uniform UV background, self-shielding approximation, star formation above $n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$, “blastwave” SN feedback, early thermal feedback, and turbulent diffusion (Fry et al., 2015). GIZMO’s mesh-free finite-mass solver with FIRE-2 physics is used in another dwarf suite, while SWIFT with SPHENIX and SWIFT-EAGLE modules underlies the TangoSIDM project (Robles et al., 2017, Correa et al., 2024).

The formal SHH interpolation method differs from these cosmological setups by hybridizing the equations of motion themselves. Starting from the Boltzmann equation, it introduces a local Knudsen number,

$Kn \equiv \lambda_{\rm mfp}/\ell_c,$

with $\lambda_{\rm mfp}=1/(\sigma_m \rho)$, and a bridging function

$$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$

The hybrid updates then become

$$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$

and

$\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$0

As $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$1, the fluid solver dominates; as $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$2, the stochastic SIDM solver dominates (Schon et al., 14 Aug 2025).

3. Cosmological implementations across halo-mass regimes

At dwarf scales, SHH simulations have been used to compare SIDM core formation with baryon-driven potential fluctuations. In the Bastidas-Fry et al. suite, a constant $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$3 was chosen to maximize SIDM effects, and the halos span from $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$4 field dwarfs up to $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$5 systems (Fry et al., 2015). The principal finding is strongly mass dependent: in very low-mass halos with $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$6, low central densities and low velocity dispersions imply $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$7 across the central $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$8, so fixed-cross-section SIDM does not significantly lower the central density. In more massive dwarfs, both CDM + hydro and SIDM + hydro form cores once bursty feedback is efficient (Fry et al., 2015).

The FIRE-2 SHH suite focuses on four isolated halos with $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$9–$\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$0, re-run in CDM and SIDM with and without baryons (Robles et al., 2017). Here the key contrast with the ChaNGa results is that every SIDM hydrodynamic simulation with $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$1 yields a flatter central profile, with $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$2 over $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$3–$\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$4, while CDM-Hydro remains steeper than $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$5 at low stellar mass and only becomes core-like in the most massive dwarf (Robles et al., 2017). This suggests that the sensitivity of SHH results to baryons depends not only on halo mass but also on the adopted feedback model.

At Milky Way scale, SHH takes on a different role. A suite of six realizations of a single MW-mass halo, run with IllustrisTNG baryonic physics and an ETHOS-3 velocity-dependent SIDM model, finds that baryons alter the inner thermal structure more strongly than SIDM scatterings for galaxies above $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$6 (Rose et al., 2022). In these runs, SIDM does not form a central core once baryons are included; rather, baryon contraction steepens the density profile and raises the central velocity dispersion, producing outward heat flow and a mild core-collapse tendency (Rose et al., 2022).

At early-type-galaxy and massive-galaxy scales, SHH has been used to assess halo shapes, density profiles, and lensing observables. Despali et al. simulate five halos hosting early-type galaxies with $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$7 and $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$8 in four variants, while Mastromarino et al. analyze $\Gamma(r)=\rho(r)\,\sigma\,\bar v(r)$9 cosmological boxes of massive galaxies under CDM, SIDM1, and velocity-dependent SIDM, each with dark-matter-only and hydro versions (Despali et al., 2022, Mastromarino et al., 2022). At cluster scale, BAHAMAS-SIDM extends the SHH logic to a $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$0 volume and focuses on X-ray and strong-lensing tests (Robertson et al., 2018).

4. Structural predictions: cores, cusps, shapes, and thermalization

A recurring SHH result is that baryons generally reduce the contrast between CDM and SIDM inferred from dark-matter-only calculations. In the early-type-galaxy zooms, the dark-matter-only inner halo at $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$1 is significantly rounder in SIDM, with axis ratios changing from $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$2 to $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$3, but once baryons are included the inner shapes become $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$4 for CDM-h and $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$5 for SIDM-h, and the projected ellipticity difference at $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$6 is reduced from $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$7 versus $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$8 in dark-matter-only runs to $$\tau_{\rm SIDM}(r)\equiv[\rho(r)\,\sigma\,\bar v(r)]^{-1}$$9 versus $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$0 in hydro runs (Despali et al., 2022).

Massive-galaxy SHH simulations show a similarly strong baryonic masking effect in density profiles. In hydro runs, low-mass systems with $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$1–$$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$2 are denser in SIDM1b than in CDMb inside the inner $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$3, intermediate-mass systems show nearly coincident profiles, and only the highest masses develop shallow mild cores in SIDM1b (Mastromarino et al., 2022). The ratio $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$4 reaches $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$5–$$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$6 for SIDM1 but only $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$7–$$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$8 for CDM at $$P_{ij}=(\sigma/m)\,|v_i-v_j|\,\Delta t/V_{\mathrm{kernel}}$$9–$n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$0, so the net CDM–SIDM difference in hydro profiles is under $n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$1 per cent over $n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$2–$n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$3 (Mastromarino et al., 2022).

At dwarf scale, the structural picture is not uniform across SHH implementations. In the ChaNGa runs, once hydrodynamics and bursty feedback are included, both CDM + hydro and SIDM + hydro develop cores in systems with $n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$4–$n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$5, and enclosed-mass differences within $n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$6–$n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$7 become essentially zero once feedback-generated cores are present (Fry et al., 2015). In the FIRE-2 runs, by contrast, SIDM-Hydro remains distinct from CDM-Hydro: SIDM runs, with or without baryons, retain $n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$8–$n_{\rm thresh}=100\ \mathrm{amu\ cm^{-3}}$9 cores and inner slopes between roughly $Kn \equiv \lambda_{\rm mfp}/\ell_c,$0 and $Kn \equiv \lambda_{\rm mfp}/\ell_c,$1, while CDM-Hydro develops comparable cores only at the highest stellar mass in the sample (Robles et al., 2017). This indicates that SHH does not imply a single baryonic outcome; it is a framework for quantifying that outcome under a specified galaxy-formation model.

The thermal interpretation becomes explicit in MW-like SHH calculations. In dark-matter-only SIDM, the halo develops a $Kn \equiv \lambda_{\rm mfp}/\ell_c,$2 core and a flat $Kn \equiv \lambda_{\rm mfp}/\ell_c,$3 plateau inside $Kn \equiv \lambda_{\rm mfp}/\ell_c,$4. In SHH runs with baryons, the central region instead exhibits a large isothermal plateau from $Kn \equiv \lambda_{\rm mfp}/\ell_c,$5 down to the softening scale, while the density profile remains cusped inward of $Kn \equiv \lambda_{\rm mfp}/\ell_c,$6 (Rose et al., 2022). The stated mechanism is that baryon contraction makes the center the hottest region, reversing the temperature gradient so that SIDM transports heat outward and removes thermal support from the center (Rose et al., 2022).

5. Observable consequences and empirical constraints

SHH methods are especially important when predictions are compared directly with lensing, X-ray, and kinematic data. In early-type galaxies, projected ellipticities within $Kn \equiv \lambda_{\rm mfp}/\ell_c,$7 were compared to three HST strong-lens samples: SLACS, BELLS, and SL2S (Despali et al., 2022). A Kolmogorov–Smirnov test shows that all three samples strongly reject CDM dark-matter-only predictions with $Kn \equiv \lambda_{\rm mfp}/\ell_c,$8; SLACS is most compatible with CDM-h at $Kn \equiv \lambda_{\rm mfp}/\ell_c,$9 but is also marginally compatible with SIDM-h; BELLS and SL2S show a slight preference for SIDM-h; and the combined sample of 152 lenses mildly favors CDM-h with $\lambda_{\rm mfp}=1/(\sigma_m \rho)$0 over SIDM-h with $\lambda_{\rm mfp}=1/(\sigma_m \rho)$1 (Despali et al., 2022). The same study finds that mock Chandra and XMM-Newton X-ray isophotes agree with observational measurements of nine relaxed ellipticals for both CDM-h and SIDM-h over $\lambda_{\rm mfp}=1/(\sigma_m \rho)$2 (Despali et al., 2022).

At cluster scale, SHH simulations enable strong-lensing constraints that are not accessible in dark-matter-only runs. BAHAMAS-SIDM finds that including baryons restores strong-lensing critical curves in clusters that would produce none in SIDM-only calculations, but the effective Einstein radii remain systematically smaller than in CDM (Robertson et al., 2018). When compared with 19 X-ray-selected CLASH clusters, the result is that an SIDM cross section $\lambda_{\rm mfp}=1/(\sigma_m \rho)$3 is in tension with observed cluster lensing strengths at velocities around $\lambda_{\rm mfp}=1/(\sigma_m \rho)$4 (Robertson et al., 2018).

For massive galaxies, the observables most sensitive to SHH effects are not total density profiles but decomposed inner mass budgets. In the EAGLE-based SHH runs, the distribution of strong-lensing Einstein radii from CDM and SIDM hydro simulations is comparable to observed samples of strong lenses with $\lambda_{\rm mfp}=1/(\sigma_m \rho)$5, whereas the dark-matter fraction inside the stellar half-mass radius differs more distinctly: for $\lambda_{\rm mfp}=1/(\sigma_m \rho)$6, SIDM1b yields central dark-matter fractions typically about $\lambda_{\rm mfp}=1/(\sigma_m \rho)$7 per cent lower than CDMb in 3D and lower by $\lambda_{\rm mfp}=1/(\sigma_m \rho)$8–$\lambda_{\rm mfp}=1/(\sigma_m \rho)$9 per cent in projection (Mastromarino et al., 2022). The paper therefore emphasizes techniques that separate luminous and dark mass in the inner regions (Mastromarino et al., 2022).

A different SHH observable appears in the TangoSIDM project, where the stellar-mass Tully–Fisher relation is used to constrain velocity-dependent SIDM (Correa et al., 2024). In Milky Way-like galaxies, where baryons dominate the central potential, SIDM thermalizes and dark matter accumulates in the central regions, producing density profiles steeper than CDM with adiabatic contraction (Correa et al., 2024). The resulting extra central mass shifts $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$0 upward at fixed stellar mass and steepens the fitted relation: the observed compiled sample gives $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$1, whereas SigmaVel30 yields $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$2 with $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$3, and SigmaVel60 yields $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$4 with $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$5 for $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$6 and $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$7 when restricted to $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$8 (Correa et al., 2024).

6. Degeneracies, misconceptions, and open directions

A common misconception is that halo shapes or inner density slopes by themselves provide robust, baryon-independent SIDM constraints. The SHH literature repeatedly rejects that view. For early-type galaxies, previous shape-based constraints suggesting $$h(\rho)=\frac{1}{1+(\rho_{\rm SHH}/\rho)^{n_{\rm SHH}}}.$$9 are described as overestimating discrimination power because they relied on dark-matter-only predictions (Despali et al., 2022). For MW-like halos, the statement is stronger: any simulations used to constrain SIDM cross sections for galaxies with stellar masses between $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$0 and at least $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$1 require baryons to make accurate predictions (Rose et al., 2022).

Another misconception is that baryons always wash out SIDM signatures. SHH results are more conditional. In some dwarf regimes, baryons and SIDM generate nearly indistinguishable cores, as in the ChaNGa suite above $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$2–$$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$3 (Fry et al., 2015). In the FIRE-2 dwarf runs, however, SIDM remains comparatively robust to baryonic uncertainty, and galaxies with $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$4 are identified as potentially ideal targets because SIDM predicts substantial cores while CDM predicts cusps (Robles et al., 2017). At Milky Way scale and above, several SHH studies instead find the opposite trend: baryon-dominated centers can make SIDM halos cuspier than their CDM counterparts through contraction and outward heat transfer (Rose et al., 2022, Correa et al., 2024).

The formal SHH method of 2025 identifies a separate limitation of standard Monte Carlo SIDM solvers in the frequent-scattering regime. In isolated NFW halos with $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$5 and $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$6, pure-SIDM runs stall near $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$7 as timesteps become extremely small, whereas SHH runs with a fiducial bridging function proceed to $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$8 and achieve two orders of magnitude higher central densities in about $$\frac{d u_i}{dt}=h(\rho_i)\left[-\frac{1}{\rho_i}\nabla P\right]_i+(1-h(\rho_i))\,{\rm Kick}_{\rm SIDM,i}-\nabla\Phi_i,$$9 less wall-clock time (Schon et al., 14 Aug 2025). This suggests that, beyond galaxy formation, SHH may also become a method for following gravothermal collapse in dense SIDM cores.

Open directions are stated explicitly across the literature. Early-type-galaxy studies call for larger simulated samples, unified end-to-end mock-observation pipelines for lensing, X-ray, and stellar kinematics, and exploration of velocity-dependent or inelastic SIDM microphysics (Despali et al., 2022). MW-like SHH work introduces a modified-DMO tracer prescription that reproduces full-hydro inner profiles to better than $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$00 per cent inside $\rho\,\sigma_T\,v_{\rm rel}\,\Delta t$01 and recommends its use as a lower-cost proxy for baryonic effects in parameter surveys (Rose et al., 2022). The formal SHH method points to the inclusion of non-ideal fluid terms, conductive and viscous corrections, and dark-sector cooling as the next step beyond ideal-fluid interpolation (Schon et al., 14 Aug 2025).

Taken together, these studies define SHH less as a single code or single subgrid model than as a methodological principle: SIDM predictions must be evaluated in tandem with baryonic structure formation, and in some regimes with an explicit interpolation between kinetic and fluid limits. The principal implication is that reliable SIDM inference depends on the mass scale, the cross-section model, the feedback prescription, and the observable under consideration, rather than on dark-matter-only intuition alone.