Our goal in this paper is a computational method which can reproduce equilibrium binding mode pop-ulations with much less computational time than treating binding modes separately. Specifically, each alchemical binding free energy calculation requires simulation at N different alchemical intermediate states (where N is typically at least 12-20 ), where each alchemical state is associated with a vector of alchemical parameters (which we will refer to as values). If we consider M different binding modes, the total cost of a binding free energy calculation that covers all binding modes separately is MNx where x is the cost of a single simulation (Figure (A) and (C)). This becomes impractical as the number of potential binding modes grows. Instead, the approach we envision is one where we calculate an absolute binding free energy for a single, reasonably populated binding mode and, in a separate calculation that can be run concurrently, efficiently determine the relative free energies (or equilibrium populations) of all M potential binding modes (Figure ). Then, we can combine the populations of the individual binding modes and the free energy estimate for a single binding mode into a binding free energy that includes all of the possible binding modes. Thus this approach would have a computational cost of Nx + y, where y is the simulation time to determine the binding mode populations. Such an approach could work as outlined in Ref., 1 providing a way to compute interconversion free energies between different metastable binding modes. This would have implications for both absolute and relative binding free energy calculations. For absolute calculations, depending on how many binding modes are being considered, such an approach could drastically reduce the total amount of simulation time, as long as y ⌧ Nx (see Figure )) (and the wallclock time as long as y ⌧ x since the N independent alchemical simulations could be run in parallel). This is not currently feasible because we have no suitable, general-purpose method for efficiently sampling binding modes, and thus the cost of these calculations (y) is far too expensive in terms of both human effort and computational time. Our focus here is on developing a method for obtaining binding mode populations which has a cost y which is relatively favorable compared to Nx, or ideally even x so that a parallel calculation could complete in at most x wall-clock time.