Numerical tests were made, using least squares fitting of a spherical harmonic model, to a selection of Magsat data to determine the practical limits of this technique with modern computers. The resulting (M102189) model, whose coefficients were adjusted up to n=50, was compared with M07AV6, a previous model which used least squares (on vector data) for coefficients up to n=29, and Gauss-Legendre quadrature (on Z residuals) to adjust the coefficients up to n=63. For the new least squares adjustment to n=50 a condition number of 115 was obtained for the solution matrix, with a resulting precision of 11 significant figures. The M102189 model shows a lower and more Gaussian residual distribution than did M07AV6, though the Gaussian envelope fits to the residual distributions, even for the scalar field, gives "standard deviations" never lower than 6 nT, a factor of three higher than the estimated Magsat observational errors. Ionospheric currents are noted to have a significant effect on the coefficients of the internal potential functions.