The accurate representation of nitrogen-containing heterocycles is essential for modeling biological systems. In this study, the general effective fragment potential (EFP2) method is used to model dimers of benzene and pyridine, complexes for which high-level theoretical data including large basis spin-component-scaled second-order perturbation theory (SCS-MP2), symmetry-adapted perturbation theory (SAPT), and coupled cluster with singles, doubles, and perturbative triples (CCSD(T))-are available. An extensive comparison of potential energy curves and components of the interaction energy is presented for sandwich, T-shaped, parallel displaced, and hydrogen-bonded structures of these dimers. EFP2 and CCSD(T) potential energy curves for the sandwich, T-shaped, and hydrogen-bonded dimers have an average root-mean-square deviation (RMSD) of 0.49 kcal/mol; EFP2 and SCS-MP2 curves for the parallel displaced dimers have an average RMSD of 0.52 kcal/mol. Additionally, results are presented from an EFP2 Monte Carlo/ simulated annealing (MC/SA) computation to sample the potential energy surface of the benzene-pyridine and pyridine dimers.