  Extending our previous studies, we perform high-resolution simulations of
inspiraling binary neutron stars in numerical relativity. We thoroughly carry
through a convergence study in our currently available computational resources
with the smallest grid spacing of $\approx 63$--86~meter for the neutron-star
radius 10.9--13.7\,km. The estimated total error in the gravitational-wave
phase is of order 0.1~rad for the total phase of $\gtrsim 210$\,rad in the last
$\sim 15$--16 inspiral orbits. We then compare the waveforms (without
resolution extrapolation) with those calculated by the latest
effective-one-body formalism (tidal SEOBv2 model referred to as TEOB model). We
find that for any of our models of binary neutron stars, the waveforms
calculated by the TEOB formalism agree with the numerical-relativity waveforms
up to $\approx 3$\,ms before the peak of the gravitational-wave amplitude is
reached: For this late inspiral stage, the total phase error is $\lesssim
0.1$\,rad. Although the gravitational waveforms have an inspiral-type feature
for the last $\sim 3$\,ms, this stage cannot be well reproduced by the current
TEOB formalism, in particular, for neutron stars with large tidal deformability
(i.e., lager radius). The reason for this is described.
