The paper presents the development of a novel approach to the solution of detailed chemistry in internal combustion engine simulations, which relies on the analytical computation of the ordinary differential equations (ODE) system Jacobian matrix in sparse form. Arbitrary reaction behaviors in either Arrhenius, third-body or fall-off formulations can be considered, and thermodynamic gas-phase mixture properties are evaluated according to the well-established 7-coefficient JANAF polynomial form. The current work presents a full validation of the new chemistry solver when coupled to the KIVA-4 code, through modeling of a single cylinder Caterpillar 3401 heavy-duty engine, running in two-stage combustion mode. The code has been tested on a wide range of simulations, at different injection timings, intake pressures, and EGR mass fractions, and considering two reaction mechanisms: a skeletal one with 29 species and 52 reactions, and a comprehensive, semi-detailed one with 160 species and 1540 reactions. The results show that the developed approach allows computational time savings of more than one order of magnitude in comparison to a reference chemistry solver, even with no reduction of the combustion mechanism size.