Bayesian inference for stochastic models is often challenging because evaluating the likelihood function typically requires integrating over a large number of latent variables. However, if only few parameters need to be inferred, it can be more efficient to perform the inference based on a comparison of the observations with (a large number of) model simulations, in terms of only few summary statistics. In Machine Learning (ML), Simulation Based Inference (SBI) using neural density estimation is often considered superior to the traditional sampling-based approach known as Approximate Bayesian Computation (ABC). Here, we present a new set of ABC algorithms based on Simulated Annealing and demonstrate that they are competitive with ML approaches, whilst requiring much less hyper-parameter tuning. For the design of these sampling algorithms we draw intuition from non-equilibrium thermodynamics, where we associate each summary statistic with a state variable (energy) quantifying the distance to the observed value as well as a temperature that controls the degree to which the associated statistic contributes to the posterior. We derive an optimal annealing schedule on a Riemannian manifold of state variables based on a minimal entropy production principle. Our new algorithms generalize the established Simulated Annealing based ABC to multiple state variables and temperatures. In situations where the information-content is unevenly distributed among the summary statistics, this can greatly improve performance of the algorithm. Our method also allows monitoring the convergence of individual statistics, which is a great diagnostic tool in out-of-sample situations. We validate our approach on standard benchmark tasks from the SBI literature and a hard inference problem from solar physics and demonstrate that it is highly competitive with the state-of-the-art.