We investigate the occurrence of multiple jet zonal flows in the 2D rotating annulus model, extended to include the possibility of boundary friction. We consider Rayleigh numbers up to 10 times critical. Without boundary friction the majority of our solutions are single-jet zonal flows, but when boundary friction is present, persistent multiple jet solutions are found much more easily. Compared to the stress-free case, the number of jets increases, though the strength of the zonal flow decreases. The dependence of these features on Ekman and Rayleigh number is discussed, suggesting that at values well beyond the reach of present 3D simulations, solutions resembling the observed jovian zonal flow may exist. The boundary condition at the metallic/insulating hydrogen interface will be an important part of any explanation of the occurrence of multiple jets on the planet Jupiter.