Continuum modeling using the Lennard-Jones potential has been shown to provide a good estimation for the interaction energy between regular-shaped homogeneous molecules comprising the same type of atoms. However, this method may not be accurate for heterogeneous molecules, which are made up of more than one chemical element. The traditional method to deal with this involves approximating the molecule via multiple surfaces in a piecemeal fashion. While this approach works well for small sized molecules, calculations become intensive for large sized molecules as a large number of sums from multiple surface interactions are involved. To address this issue, we propose a new model that approximates a heterogeneous molecule with a single surface or volume, where attractive and repulsive constants (A and B) in the Lennard-Jones potential are replaced by functions A (r) and B (r), which depend on the parameterization of the surface r. We comment that this technique is suitable for regular-shaped nanostructures where their heterogeneity can be modeled by surface (or volume) parameterization. Validation of the new approach is carried out via two problems, namely, carbon nanotube-methane and carbon nanotube-coronene interactions. For coronene and methane, which are assumed to be radially symmetric, we propose A (r) and B (r) to be sigmoidal functions for which the interaction strength decreases from the inner region of the carbon atoms toward the outer region of the hydrogen atoms. Our results for both cases show that using the sigmoidal profiles for A (r) and B (r) gives rise to interaction energies that are in better agreement with those obtained from molecular dynamics studies compared to results using constant A and B. The new approach provides a significant improvement to the current continuum modeling using the Lennard-Jones potential.