  Sub-photospheric shock dissipation is one of the main proposed mechanisms for
producing the prompt gamma-ray burst (GRB) emission. Such shocks are mediated
by scattering of radiation. We introduce a time dependent, special relativistic
code which dynamically couples Monte Carlo radiative transfer to the flow
hydrodynamics. The code also self-consistently implements electron-positron
pair production and annihilation. We simulate shocks with properties relevant
for GRBs and study the steady-state solutions, which are accurate deep below
the jet photosphere. The shock generates a power-law photon spectrum through
the first-order Fermi mechanism, extending upwards from the typical upstream
photon energy. Strong shocks (for which the downstream pressure is much larger
than the upstream pressure) have rising $\nu F_\nu$ shock spectra. The spectrum
extends up to $\epsilon_{max} \equiv E_{max}/m_e c^2 \sim v^2$ for
non-relativistic shocks, where $m_e$ is the electron rest mass and $v$ is the
relative speed between the upstream and downstream in units of the speed of
light $c$. For mildly relativistic shocks the power law softens at $\epsilon
\gtrsim 10^{-1}$ due to Klein-Nishina effects, and shocks with $v\gamma \gtrsim
1$, where $\gamma \equiv (1-v^2)^{-1/2}$, produce electron-positron pairs. As
an example, a strong shock with $v\gamma = 3$ and a photon-to-proton ratio of
$n_\gamma/n_p = 2 \times 10^5$ has a peak pair-to-proton ratio of $Z_\pm
\approx 225$. The main effect of pairs in a steady-state shock is to decrease
its spatial width by a factor of $\sim Z_\pm$. The post-shock spectrum
thermalizes in the downstream. In absence of emission and absorption processes,
kinetic equilibrium at temperature $\theta_d \equiv kT_d/m_e c^2 \approx
\epsilon_d/3$ is reached at an optical depth of $\tau \gg \theta_d^{-1}$ behind
the shock, where $\epsilon_d$ is the average downstream photon energy.
