Merging binary neutron stars are among the strongest known sources of gravitational waves, have features compatible with the events producing short-hard gamma-ray bursts and represent formation sites of high-mass number r-process elements. Their fully nonlinear numerical modelling, however, constitutes a major challenge. Since the pioneering simulations of the early 2000s numerical relativity has reached a stage where a complete description of the inspiral, merger, and post-merger phases of the late evolution of close binary neutron star systems is possible. This lecture will review the basic framework that is required to model the dynamics of these systems with a focus on the hydrodynamics and MHD equations along with the methods best suited for obtaining accurate numerical solutions. Key results from the growing list of available simulations will also be discussed.