ABSTRACT Mixed integer linear programming (MILP) has been gaining traction as a method for solving optimal pump scheduling problems in water distribution networks (WDNs). However, inclusion of variable speed pumps (VSPs) in MILP pump scheduling frameworks has not been given adequate treatment. This article addresses this gap by describing a methodology for formulating and solving optimal pump scheduling problems with VSPs using MILP and piece-linear approximations of network components. The methodology proceeds in four steps: (a) WDN simulation with initial pump schedule(s), (b) approximation of network components, including VSP, using linear and piece-linear functions around the chosen operating points, (c) formulation of a fully parameterised mixed integer linear programme, and (d) solution of the optimisation problem and WDN simulation with optimal pump schedule(s). The methodology is coded in MATLAB/OCTAVE and Python and is publicly available on GitHub. It was applied to solve a pump scheduling problem on a two variable speed pump single-tank network that allows the reader to easily understand how the methodology works and how it is applied in practice. The results show that the formulation is robust and the optimiser is able to return a globally optimal solution for a range of operating points.