A grand canonical ensemble Monte Carlo simulation is performed to investigate the adsorption, heat of adsorption, and distributions of thiophene in all-silica Y and Na-Y zeolites. Biased particle insertions and deletions were implemented to allow the computation of equilibrium adsorption isotherms of such molecules. The calculated number of absorbed thiophene molecules in these zeolites is in good agreement with the experimental data. The calculated results show that the number absorbed of thiophene molecules in Na-Y is much greater than that in all-silica Y over the range of pressure. The calculated heat of adsorption is in good agreement with experimental results. The Na-Y zeolite, rather than all-silica Y, preferentially adsorbs the thiophene. A distribution analysis of the adsorbed phase structure reveals a different adsorption site in the zeolites.