Known formulas for variational bounds on Darcy's constant for slow flow through porous media depend on two-point and three-point spatial correlation functions. Certain bounds due to Prager and Doi depending only a two-point correlation functions have been calculated for the first time for random aggregates of spheres with packing fractions ( η) up to η = 0.64. Three radial distribution functions for hard spheres were tested for η up to 0.49: (1) the uniform distribution or "well-stirred approximation," (2) the Percus-Yevick approximation, and (3) the semi-empirical distribution of Verlet and Weis. The empirical radial distribution functions of Bennett and Finney were used for packing fractions near the random-close-packing limit ( ηRCP ∼ 0.64). An accurate multidimensional Monte Carlo integration method (VEGAS) developed by Lepage was used to compute the required two-point correlation functions. The results show that Doi's bounds are preferred for η < 0.10 while Prager's bounds are preferred for η > 0.10. The "upper bounds" computed using the well-stirred approximation actually become negative (which is physically impossible) as η increases, indicating the very limited value of this approximation. The other two choices of radial distribution function give reasonable results for η up to 0.49. However, these bounds do not decrease with η as fast as expected for large η. It is concluded that variational bounds dependent on three-point correlation functions are required to obtain more accurate bounds on Darcy's constant for large η.