We develop a combined machine learning (ML) and quantum mechanics approach that enables data-efficient reconstruction of flexible molecular force fields from high-level ab initio calculations, through the consideration of fundamental physical constraints. We discuss how such constraints are recovered and incorporated into ML models. Specifically, we use conservation of energy—a fundamental property of closed classical and quantum mechanical systems—to derive an efficient gradient-domain machine learning (GDML) model. The challenge of constructing conservative force fields is accomplished by learning in a Hilbert space of vector-valued functions that obey the law of energy conservation. We proceed with the development of a multi-partite matching algorithm that enables a fully automated recovery of physically relevant point group and fluxional symmetries from the training dataset into a symmetric variant of our model. The symmetric GDML (sGDML) approach is able to faithfully reproduce global force fields at the accuracy high-level ab initio methods, thus enabling sample intensive tasks like molecular dynamics simulations at that level of accuracy. (This chapter is adapted with permission from Chmiela (Towards exact molecular dynamics simulations with invariant machine-learned models, PhD thesis. Technische Universität, Berlin, 2019).).