The availability of accurate numerical waveforms is an important requirement
for the creation and calibration of reliable waveform models for gravitational
wave astrophysics. For black hole-neutron star binaries, very few accurate
waveforms are however publicly available. Most recent models are calibrated to
a large number of older simulations with good parameter space coverage for
low-spin non-precessing binaries but limited accuracy, and a much smaller
number of longer, more recent simulations limited to non-spinning black holes.
In this paper, we present long, accurate numerical waveforms for three new
systems that include rapidly spinning black holes, and one precessing
configuration. We study in detail the accuracy of the simulations, and in
particular perform for the first time in the context of BHNS binaries a
detailed comparison of waveform extrapolation methods to the results of Cauchy
Characteristic Extraction. The new waveforms have $<0.1\,{\rm rad}$ phase
errors during inspiral, rising to $\sim (0.2-0.4)\,{\rm rad}$ errors at merger,
and $\lesssim 1\%$ error in their amplitude. We compute the faithfulness of
recent analytical models to these numerical results, and find that models
specifically designed for BHNS binaries perform well ($F>0.99$) for binaries
seen face-on. For edge-on observations, particularly for precessing systems,
disagreements between models and simulations increase, and models that include
precession and/or higher-order modes start to perform better than BHNS models
that currently lack these features.
