Relativistic sub-photospheric shocks are a possible mechanism for producing 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 follows electron-positron pair production in photon-photon collisions. We use the code to simulate shocks with properties relevant to GRBs. We focus on plane-parallel solutions, which are accurate deep below the photosphere. The shock generates a power-law photon spectrum through the first-order Fermi mechanism, extending upward from the typical upstream photon energy. Strong (high Mach number) shocks produce rising nu F-nu spectra. We observe that in non-relativistic shocks the spectrum extends to E-max similar to m(e)v(2), where v is the speed difference between the upstream and downstream. In relativistic shocks the spectrum extends to energies E > 0.1 m(e)c(2) where its slope softens due to Klein-Nishina effects. Shocks with Lorentz factors gamma > 1.5 are prolific producers of electron-positron pairs, yielding hundreds of pairs per proton. The main effect of pairs is to reduce the shock width by a factor of similar to Z(+/-)(-1). Most pairs annihilate far downstream of the shock, and the radiation spectrum relaxes to a Wien distribution, reaching equilibrium with the plasma at a temperature determined by the shock jump conditions and the photon number per proton. We discuss the implications of our results for observations of radiation generated by sub-photospheric shocks.