Heat conduction in dielectric crystals originates from the propagation of atomic vibrations, whose microscopic dynamics is well described by the linearized phonon Boltzmann transport equation. Recently, it was shown that thermal conductivity can be resolved exactly and in a closed form as a sum over relaxons, $\mathit{i.e.}$ collective phonon excitations that are the eigenvectors of Boltzmann equation's scattering matrix [Cepellotti and Marzari, PRX $\mathbf{6}$ (2016)]. Relaxons have a well-defined parity, and only odd relaxons contribute to the thermal conductivity. Here, we show that the complementary set of even relaxons determines another quantity --- the thermal viscosity --- that enters into the description of heat transport, and is especially relevant in the hydrodynamic regime, where dissipation of crystal momentum by Umklapp scattering phases out. We also show how the thermal conductivity and viscosity parametrize two novel viscous heat equations --- two coupled equations for the temperature and drift-velocity fields --- which represent the thermal counterpart of the Navier-Stokes equations of hydrodynamics in the linear, laminar regime. These viscous heat equations are derived from a coarse-graining of the linearized Boltzmann transport equation for phonons, and encompass both limits of Fourier's law and of second sound, taking place, respectively, in the regimes of strong or weak momentum dissipation. Last, we introduce the Fourier deviation number as a descriptor that captures the deviations from Fourier's law due to hydrodynamic effects. We showcase these findings in a test case of a complex-shaped device made of graphite, obtaining a remarkable agreement with the very recent experimental demonstration of hydrodynamic transport in this material. The present findings also suggest that hydrodynamic behavior can appear at room temperature in micrometer-sized diamond crystals.