  Out-of-equilibrium reactions between different particle species are the main
processes contributing to bulk viscosity in neutron stars. In this work, we
numerically compare three different approaches to the modeling of bulk
viscosity: the multi-component fluid with reacting particle species and two
bulk stress formalism based on the M\"uller-Israel-Stewart theory, namely the
Hiscock-Lindblom and the Maxwell-Cattaneo models, whose flux-conservative
formulation in radial gauge-polar slicing coordinates and spherical symmetry is
derived in a companion paper. To our knowledge, this is the first time that a
neutron star is simulated with the complete Hiscock-Lindblom model of bulk
viscosity. We find that the Hiscock-Lindblom and Maxwell-Cattaneo models are
good approximations of the multi-component fluid for small perturbations and
when the non-equilibrium equation of state of the fluid depends on only one
independent particle fraction. For more than one independent particle fraction
and for large perturbations, the bulk stress approximation is still valid but
less accurate. In addition, we include the energy loss due to the luminosity of
the reactions in the bulk stress formulation. We find that the energy loss due
to bulk viscosity has a larger effect on the dynamics than the bulk stress or
the variation in particle composition per se. The new one-dimensional,
general-relativistic hydrodynamic code developed for this work, hydro-bulk-1D,
is publicly available.
