We have investigated and compared the ability of numerical and Gaussian-type basis sets to accurately describe the geometries and binding energies of a selection of hydrogen bonded systems that are well studied theoretically and experimentally. The numerical basis sets produced accurate results for geometric parameters but tended to overestimate binding energies. However, a comparison of the time taken to optimize phosphinic acid dimer, the largest complex considered in this study, shows that calculations using numerical basis sets offer a definitive advantage where geometry optimization of large systems is required.