The free energy profiles of the metal-ligand exchange reactions were evaluated by umbrella sampling (US) simulations using the i-PI driver combined with the Plumed v2.9.0 library. The PMFs were constructed using the Weighted Histogram Analysis Method (WHAM) code v2.0.11. For Mg 2+ in water, the Mg•••O distance was chosen as a reaction coordinate (RC). The starting structures for each window were generated by a steered MD, pulling a water molecule from the first solvation shell to a distance ranging from 1.5 to 7 Å. The US covered the distance from 1.5 to 6.0 Å, split into 48 windows with a spacing of 0.075 Å. In each window, the trajectory was propagated for 1 ns, with the first 50 ps used for equilibration and skipped from the analysis, corresponding to 45.6 ns of sampling time. The US windows were propagated independently at 300 K. A harmonic umbrella restraint with force constant 500 kcal/(mol• Å2 ) was employed in the windows from 1.5 to 3.00 Å, while the force constant was lowered to 250 kcal/(mol• Å2 ) in the region from 3.00 to 6.00 Å. The PMF was computed as an average from 2 US runs with random seeds for generating the initial velocity distributions, with an uncertainty estimated as a standard deviation. The final PMF was corrected by the entropy term 2k B T ln(r/r 0 ) that accounts for the increasing volume of configuration space with increasing distance r. For Pd 2+ in acetonitrile, the coordination number was chosen as a collective variable for the US (further details in section SI §S4.2). The starting structures for the US runs were generated by steered MD, during which one of the two MeCN axial ligands was pushed towards the Pd centre, guiding the system towards the ligand exchange event through the formation of the pentacoordinated TS. In the US, the coordination number was varied from