DFT/B3LYP calculations were carried out on complexes formed by NH4+ with aromatics, viz. benzene, phenol, pyrrole, imidazole, pyridine, indole, furane, and thiophene, to characterize the forces involved in such interactions and to gain further insight into the nature and diversity of cation-aromatic interactions. Such calculations may provide valuable information for understanding molecular recognition in biological systems and for force-field development. B3LYP/6-31G** optimization on 35 initial structures resulted in 11 different finally optimized geometries, which could be divided into three types: NH4+-pi complexes, protonated heterocyclic-NH3 hydrogen bond complexes, and heterocyclic-NH4+ hydrogen bond complexes. For NH4+-pi complexes, NH4+ always tilts toward the carbon-carbon bond rather than toward the heteroatom or the carbon-heteroatom bond. The calculated CHelpG charges suggest that the charge distribution of a free heterocyclic may be used to predict the geometry of its complex. Charge population and electrostatic interaction estimations show that the NH4+-pi interaction has the largest nonelectrostatic interaction fraction (similar to 47%) of the total binding energy, while the NH4+-aromatic hydrogen bond interaction has the largest electrostatic fraction (similar to 90%). A good correlation between binding energy and electrostatic interaction in the NH4+-pi complexes is found, which shows that nonelectrostatic interaction is important for cation-pi binding. The results calculated with basis sets from 6-31G to 6-311++G(2df, 2dp) show that DeltaE(corr) and DeltaH(corr) do not require a basis-set superposition error (BSSE) correction, in view of experimental error, if a larger basis set is used in the calculation. The calculated DeltaH(corr) values for the NH4+-C6H6 complex with different basis sets suggest that the experimental DeltaH may be overestimated.