Abstract: The ability to observe small celestial bodies has grown drastically over the last decade. The increase in interest for these bodies has increased demand for higher fidelity trajectory simulations in order to assure mission success. Most methods that are available for simulating trajectories about asymmetric bodies assume they are of uniform density. Here we propose a modification to two well-known methods: the mascon model and the spherical harmonic series approximation, for use in simulating trajectories about variable density bodies. In particular, we will look at contact binaries which are bodies consisting of two different densities.