Hi, u can use gaussian calculations to predict mulliken charge density distribution. Once you set up the calculations in gaussian, try for the OPT+Frequency method and save NBO's. Also please dont forget to check the box for compute polarisation. After the calculations gets over, u can right click and go to charge distribution, where you can see the mulliken charge distribution and even colour the atoms based on the charge density.

http://www.chemistry4.me/Gaussian/G09W/help/results.htm

Mulliken population analysis is not good, because it lacks any physical or mathematical limits as the basis set is improved towards completeness.

Here are some things you should keep in mind when choosing a method to compute the partial atomic charges:

Net atomic charges, also called partial charges or partial atomic charges, are commonly used for two different purposes: (1) to quantify the transfer of electrons between atoms in a material as computed by quantum chemistry and (2) to construct force-fields (for classical molecular dynamics or Monte Carlo simulations) that approximately reproduce the electrostatic potential surrounding the material. While quantum chemical topology (QCT), also called quantum theory of atoms in molecules (QTAIM), often describes electron transfer in buried atoms with reasonable accuracy, the QCT atomic charges are not suitable to construct point-charge models for classical force-fields, because the QCT atomic multipoles are large in magnitude. (The QCT method can be used to construct accurate force-fields using multi-centered polyatomic multipole expansions, but these are more complicated than the simple point-charge based force-fields often used.) On the other hand, electrostatic potential derived point charges (e.g., CHELP, CHELPG, ESP, Merz-Singh-Kollman, REPEAT) have reasonable accuracy for reproducing the electrostatic potential surrounding the material, but these do not accurately describe electron transfer for buried atoms.

The density –derived electrostatic and chemical (DDEC) methods solve this problem by partitioning the electron density to simultaneously reproduce atomic chemical states and the electrostatic potential surrounding the material with excellent accuracy. The latest generation (called DDEC6) is described in these papers: (a) T. A. Manz and N. Gabaldon Limas, “Introducing DDEC6 atomic population analysis: part 1. Charge partitioning theory and methodology,” RSC Advances, 6 (2016) 47771-47801 (http://dx.doi.org/10.1039/C6RA04656H) and (b) N. Gabaldon Limas and T. A. Manz, “Introducing DDEC6 atomic population analysis: part 2. Computed results for a wide range of periodic and nonperiodic materials,” RSC Advances, 6 (2016) 45727-45747 (http://dx.doi.org/10.1039/C6RA05507A). It gives net atomic charges that are closely correlated to many experimental properties as described in those articles. The DDEC6 method is implemented in the free Chargemol program (http://ddec.sourceforge.net) which can analyze results of GAUSSIAN, VASP, and other quantum chemistry programs.