Considering the potent cytotoxic activities of hybrid benzofuran–imidazolium and quinazolinone derivatives on the breast cancer cell line MCF-7, novel hybrid derivatives incorporating benzofuran, imidazole and quinazolinone pharmacophores were designed using a molecular hybridization approach. Since aromatase is highly expressed in the MCF-7 cell line, we tried to put these pharmacophores together in such a way that they would arrange themselves in a symmetrical shape, similar to aromatase inhibitors. Subsequently, the binding of these novel hybrid compounds to aromatase have been investigated using a docking procedure applying a combined quantum mechanical/molecular mechanical (QM/MM) method. The QM/MM calculation was performed on the reference structures to obtain atomic charges on the ligand atoms. The results indicated that the hybrid compounds were adopted properly within the aromatase binding site, suggesting that they could be potential inhibitors of aromatase. These novel designed compounds engage in hydrophobic and H-bond interactions with the aromatase binding site, which are in agreement with the basic physicochemical features of known aromatase inhibitors. To obtain more accurate results for the binding energies of the ligands, the structures of the ligands with the best interaction energies, obtained from the docking study, were re-optimized using a three-layer ONIOM method (QM:QM:MM) in which the binding pocket of the enzyme was considered as medium-level. The results demonstrated that when the optimized geometrical structures were subjected to re-docking, a better interaction energy was obtained, which strengthens the ability of these compounds to act as potential inhibitors of aromatase.