We have developed an improved LevenburgMarquart technique to rapidly invert Bouguer gravity data for a 3-D density distribution as a source of the observed field. This technique is designed to replace tedious forward modeling with an automatic solver that determines density models constrained by geologic information supplied by the user. Where such information is not available, objective models are generated. The technique estimates the density distribution within the source volume using a least-squares inverse solution that is obtained iteratively by singular value decomposition using orthogonal decomposition of matrices with sequential Householder transformations. The source volume is subdivided into a series of right rectangular prisms of specified size but of unknown density. This discretization allows the construction of a system of linear equations relating the observed gravity field to the unknown density distribution.Convergence of the solution to the system is tightly controlled by a damping parameter which may be varied at each iteration. The associated algorithm generates statistical measures of solution quality not available with most forward methods. Along with the ability to handle large data sets within reasonable time constraints, the advantages of this approach are: (1) the ease with which pre-existing geological information can be included to constrain the solution, (2) its minimization of subjective user input, (3) the avoidance of difficulties encountered during wavenumber domain transformations, and (4) the objective nature of the solution.Application to a gravity data set from Hamilton County, Indiana, has yielded a geologically reasonable result that agrees with published models derived from interpretation of gravity, magnetic, seismic, and drilling data.