  The method of harmonic balance (HB) is a spectrally accurate method used to
obtain periodic steady state solutions to dynamical systems subjected to
periodic perturbations. We adapt HB to solve for the stress response of the
Giesekus model under large amplitude oscillatory shear (LAOS) deformation. HB
transforms the system of differential equations to a set of nonlinear algebraic
equations in the Fourier coefficients. Convergence studies find that the
difference between the HB and true solutions decays exponentially with the
number of harmonics ($H$) included in the ansatz as $e^{-m H}$. The decay
coefficient $m$ decreases with increasing strain amplitude, and exhibits a "U"
shaped dependence on applied frequency. The computational cost of HB increases
slightly faster than linearly with $H$. The net result of rapid convergence and
modest increase in computational cost with increasing $H$ implies that HB
outperforms the conventional method of using numerical integration to solve
differential constitutive equations under oscillatory shear. Numerical
experiments find that HB is simultaneously about three orders of magnitude
cheaper, and several orders of magnitude more accurate than numerical
integration. Thus, it offers a compelling value proposition for parameter
estimation or model selection.
