A monolithic first--order BSSNOK formulation of the Einstein--Euler equations and its solution with path-conservative finite difference CWENO schemes

We present a new, monolithic first--order (both in time and space) BSSNOK formulation of the coupled Einstein--Euler equations. The entire system of hyperbolic PDEs is solved in a completely unified manner via one single numerical scheme applied to both the conservative sector of the matter part and to the first--order strictly non--conservative sector of the spacetime evolution. The coupling between matter and space-time is achieved via algebraic source terms. The numerical scheme used for the solution of the new monolithic first order formulation is a path-conservative central WENO (CWENO) finite difference scheme, with suitable insertions to account for the presence of the non--conservative terms. By solving several crucial tests of numerical general relativity, including a stable neutron star, Riemann problems in relativistic matter with shock waves and the stable long-time evolution of single and binary puncture black holes up and beyond the binary merger, we show that our new CWENO scheme, introduced two decades ago for the compressible Euler equations of gas dynamics, can be successfully applied also to numerical general relativity, solving all equations at the same time with one single numerical method. In the future the new monolithic approach proposed in this paper may become an attractive alternative to traditional methods that couple central finite difference schemes with Kreiss-Oliger dissipation for the space-time part with totally different TVD schemes for the matter evolution and which are currently the state of the art in the field.