  To date, essentially all simulation codes that solve the full set of
Einstein's equations are performed in the framework of Eulerian hydrodynamics.
The exception is our recently developed Numerical Relativity code SPHINCS_BSSN
which solves the commonly used BSSN formulation of the Einstein equations on a
structured mesh and the matter equations via Lagrangian particles. We show
here, for the first time, SPHINCS_BSSN neutron star merger simulations with
piecewise polytropic approximations to four nuclear matter equations of state.
In this set of neutron star merger simulations we focus on perfectly symmetric
binary systems that are irrotational and have 1.3 $M_\odot$ masses. We
introduce some further methodological refinements (a new way of steering
dissipation, an improved particle-mesh mapping) and we explore the impact of
the exponent that enters in the calculation of the thermal pressure
contribution. We find that it leaves a noticeable imprint on the gravitational
wave amplitude (calculated via both quadrupole approximation and the
$\Psi_4$-formalism) and has a noticeable impact on the amount of dynamic
ejecta. Consistent with earlier findings, we only find a few times $10^{-3}$
\Msun as dynamic ejecta in the studied equal mass binary systems, with softer
equations of state (which are more prone to shock formation) ejecting larger
amounts of matter. In all of the cases, we see a credible high-velocity
($\sim0.5 .. 0.7c$) ejecta component of $\sim 10^{-4}$ \Msun that is launched
at contact from the interface between the two neutron stars. Such a
high-velocity component has been suggested to produce an early, blue precursor
to the main kilonova emission and it could also potentially cause a kilonova
afterglow.
