Owing to their increased electron mobility compared to conventional semiconductors, threedimensional (3D) Dirac semimetals are considered to be promising candidates for integration into next-generation electronic devices. In these materials, the low-energy dynamics of the charge carriers are governed by an effective tilted Dirac equation, in which a mass term appears when strain is applied to the crystal lattice. In this work, we present a novel finite-difference scheme capable of numerically solving the 3D tilted Dirac equation in the time domain. The method employs fourth-order accurate finite differences to discretize the spatial derivatives and a symplectic partitioned Runge-Kutta (PRK) integrator to propagate the Dirac spinor in time. To this end, a careful separation of the complex-valued spinor into two real-valued parts is performed to ensure compatibility with the PRK technique. Moreover, to account for the additional term in the Hamiltonian arising from the tilt of the Dirac cones, fourth-order accurate averaging operators are incorporated into the spatial discretization, without compromising the key properties of the scheme. The resulting numerical method is explicit and is shown to conserve the norm, energy, and momentum of the system. Its stability condition is derived, and the numerical dispersion thoroughly investigated. Illustrative numerical experiments are performed, demonstrating the excellent properties of the proposed method and its applicability to a more realistic scenario in which retroreflection is predicted to occur in the Dirac semimetal Cd3As2, as a result of its tilted dispersion relation.
