In constant pH MD simulations, titratable groups can dynamically change their protonation state. These changes are driven by interactions between the group and the chemical environment, modeled with a force field, and the aqueous proton concentration, modeled with a pH potential. Because at the force field level, a number of contributions to the free energy of (de)protonation are not included explicitly, i.e. quantum mechanical interactions associated with bond breakage and formation, as well as the actual proton particle, corrections to the force field are needed in λ-dynamics based constant pH MD. In GROMACS these corrections are implemented as analytical functions, V MM (λ j ), fitted to the free energy profile associated with the deprotonation of a titratable residue j at the force field level. The accuracy of such free energy profiles not only depends on how closely the force field model represents the true potential energy surface but also on the convergence of sampling of all other degrees of freedom in the system. Therefore, whereas in normal MD, the accuracy of the dynamics depends solely on the quality of the force field, the accuracy of λ-dynamics based constant pH MD depends additionally on whether all relevant degrees of freedom are sampled sufficiently in the simulations required for parameterizing the correction potentials.