Point multipole expansions are widely used to gain physical insight into complex distributions of charges and to reduce the cost of computing interactions between such distributions. However, practical applications that typically retain only a few leading terms may suffer from unacceptable loss of accuracy in the near-field. We propose an alternative approach for approximating electrostatic charge distributions, Optimal Physical Multipoles (OPMs), which optimally represent the original charge distribution with a set of point charges. By construction, approximation of electrostatic potential based on OPMs retains many of the useful properties of the corresponding point multipole expansion, including the same asymptotic behavior of the approximate potential for a given multipole order. At the same time, OPMs can be significantly more accurate in the near field: up to 5 times more accurate for some of the charge distributions tested here which are relevant to biomolecular modeling. Unlike point multipoles, for point charge distributions the OPM always converges to the original point charge distribution at finite order. Furthermore, OPMs may be more computationally efficient and easier to implement into existing molecular simulations software packages than approximation schemes based on point multipoles. In addition to providing a general framework for computing OPMs to any order, closed-form expressions for the lowest order OPMs (monopole and dipole) are derived. Thus, for some practical applications Optimal Physical Multipoles may represent a preferable alternative to point multipoles.