Molecular dynamics approaches have been applied to study the transitions between polymorphs and, in the context of CSP, to identify structures that interconvert at non-zero temperatures . In this study, we present the implementation of the threshold algorithm, which is based on Monte Carlo sampling of the energy landscape , to molecular crystals, using an accurate force field with an atomic multipole description of electrostatics. The aim of the threshold algorithm is to find the lowest energy at which transitions are possible between local energy minima. By combining trajectories from multiple starting structures, a global picture of the connectivity of minima can be constructed. The threshold algorithm has previously been applied to fairly simple inorganic crystal structures, to investigate the energy landscape of MgF 2 and to study the entropic stabilization of high energy phases of CaF2 . Here, we investigate its application to more complex molecular organic crystals, which are characterized by a balance of several types of intermolecular interactions, and discuss the insight that this method provides to help understand their polymorphism. Four crystal systems (Fig. ), including single-component crystals and a co-crystal, were chosen for study, each of which has known polymorphism. All molecules have reasonably rigid molecular structures, so that rigid-molecule simulations should provide a realistic picture of the crystal energy landscape.