  Cosmic reionization was driven by the imbalance between early sources and
sinks of ionizing radiation, both of which were dominated by small-scale
structure and are thus usually treated in cosmological reionization simulations
by subgrid modelling. The recombination rate of intergalactic hydrogen is
customarily boosted by a subgrid clumping factor,
${\left<n^2\right>/\left<n\right>^2}$, which corrects for unresolved
fluctuations in gas density ${n}$ on scales below the grid-spacing of
coarse-grained simulations. We investigate in detail the impact of this
inhomogeneous subgrid clumping on reionization and its observables, as follows:
(1) Previous attempts generally underestimated the clumping factor because of
insufficient mass resolution. We perform a high-resolution $N$-body simulation
that resolves haloes down to the pre-reionization Jeans mass to derive the
time-dependent, spatially-varying local clumping factor and a fitting formula
for its correlation with local overdensity. (2) We then perform a large-scale
$N$-body and radiative transfer simulation that accounts for this inhomogeneous
subgrid clumping by applying this clumping factor-overdensity correlation.
Boosting recombination significantly slows the expansion of ionized regions,
which delays completion of reionization and suppresses 21 cm power spectra on
large scales in the later stages of reionization. (3) We also consider a
simplified prescription in which the globally-averaged, time-evolving clumping
factor from the same high-resolution $N$-body simulation is applied uniformly
to all cells in the reionization simulation, instead. Observables computed with
this model agree fairly well with those from the inhomogeneous clumping model,
e.g. predicting 21 cm power spectra to within 20% error, suggesting it may be a
useful approximation.
