  We present axisymmetric two-temperature general relativistic radiation
magnetohydrodynamic (GRRMHD) simulations of the inner region of the accretion
flow onto the supermassive black hole M87. We address uncertainties from
previous modeling efforts through inclusion of models for (1) self-consistent
dissipative and Coulomb electron heating (2) radiation transport (3)
frequency-dependent synchrotron emission, self-absorption, and Compton
scattering. We adopt a distance $D=16.7$ Mpc, an observer angle $\theta =
20^{\circ}$, and consider black hole masses $M/M_{\odot} = (3.3\times10^{9},
6.2\times10^{9})$ and spins $a_{\star} = (0.5, 0.9375)$ in a four-simulation
suite. For each $(M, a_{\star})$, we identify the accretion rate that recovers
the 230 GHz flux from VLBI measurements. We report on disk thermodynamics at
these accretion rates ($\dot{M}/\dot{M}_{\mathrm{Edd}} \sim 10^{-5}$). The disk
remains geometrically thick; cooling does not lead to a thin disk component.
While electron heating is dominated by Coulomb rather than dissipation for $r
\gtrsim 10 GM/c^2$, the accretion disk remains two-temperature. Radiative
cooling of electrons is not negligible, especially for $r \lesssim 10 GM/c^2$.
The Compton $y$ parameter is of order unity. We then compare derived and
observed or inferred spectra, mm images, and jet powers. Simulations with
$M/M_{\odot} = 3.3\times10^{9}$ are in conflict with observations. These
simulations produce mm images that are too small, while the low-spin simulation
also overproduces X-rays. For $M/M_{\odot} = 6.2\times10^{9}$, both simulations
agree with constraints on radio/IR/X-ray fluxes and mm image sizes. Simulation
jet power is a factor $10^2-10^3$ below inferred values, a possible consequence
of the modest net magnetic flux in our models.
