Techniques an developed here for evaluating the r-modes of rotating neutron
stars through second order in the angular velocity of the star. Second-ord
er corrections to the frequencies and eigenfunctions for these modes are ev
aluated for neutron star models. The second-order eigenfunctions for these
modes are determined here by solving an unusual inhomogeneous hyperbolic bo
undary-value problem. The numerical techniques developed to solve this unus
ual problem are somewhat non-standard and may well be of interest beyond th
e particular application here. The bulk-viscosity coupling to the r-modes,
which appears first at second order, is evaluated. The bulk-viscosity time
scales are found here to be longer than previous estimates for normal neutr
on stars, but shorter than previous estimates for strange stars. These new
time scales do not substantially affect the current picture of the gravitat
ional radiation driven instability of the r-modes either for neutron stars
or for strange stars. [S0556-2821(99)01518-0].