We present a novel approach to the study of di-Higgs production via gluon-gluon fusion at the LHC. The relevant Feynman diagrams involving two Standard Model-like Higgs bosons $hh$ are computed within a simplified model approach that enables one to interpret possible signals of new physics in a model-independent way as well as to map these onto specific theories. This is possible thanks to a decomposition of such a signal process into all its squared amplitudes and their relative interferences, each of which has a well-defined coupling structure. We illustrate the power of this procedure for the case of both a minimal and next-to-minimal representation of supersymmetry, for which the new physics effects are due to top squarks entering the loops of $gg\ensuremath{\rightarrow}hh$. The squarks yield both a change of the integrated cross section and peculiar kinematic features in its differential distributions with respect to the Standard Model. These effects can in turn be traced back to the relevant diagrammatic and coupling structures and allow for a detailed analysis of the process. In order to do so, we perform systematic scans of the parameter spaces of such new physics scenarios and identify benchmark points which exhibit potentially observable features during the current and upcoming runs of the LHC.