Over the past two decades, the availability of high-performance computing and the development of neural networks of larger capacity have conspired to fuel a revolution in the way we think at the optimisation of complex systems. When the dimensionality of the space of relevant design parameters exceeds a few units, and brute-force scans cease be a viable option for its exploration. We nowadays, have the option of letting automated systems find their way to configurations that correspond to advantageous extrema of carefully specified objective functions. The engine under the hood of these optimisation searches is automatic differentiation, which allows computer programs to keep track of the gradient of the objective function, through the chain rule of differential calculus, as computer code performs arbitrarily complex successions of operations to model the behaviour of the system. Crucial to a successful optimisation of the system is the inclusion in the model of all relevant effects that have an impact on the precision of the inference that the data generated by the system may produce. An incomplete description of the inference itself, or a mock up of the reconstruction techniques performing the dimensionality reduction step which translates raw data into high-level features informing the objective function, are likely to prevent the identification of designs that maximise the true objective, as they introduce a misalignment.

Cosmological shock waves are essential to understanding the formation of cosmological structures. To study them, scientists run computationally expensive high-resolution 3D hydrodynamic simulations. Interpreting the simulation results is challenging because the resulting data sets are enormous, and the shock wave surfaces are hard to separate and classify due to their complex morphologies and multiple shock fronts intersecting. We introduce a novel pipeline, Virgo, combining physical motivation, scalability, and probabilistic robustness to tackle this unsolved unsupervised classification problem. To this end, we employ kernel principal component analysis with low-rank matrix approximations to denoise data sets of shocked particles and create labeled subsets. We perform supervised classification to recover full data resolution with stochastic variational deep kernel learning. We evaluate on three state-of-the-art data sets with varying complexity and achieve good results. The proposed pipeline runs automatically, has only a few hyperparameters, and performs well on all tested data sets. Our results are promising for large-scale applications, and we highlight now enabled future scientific work.