POLYRATE is a computer program for the calculation of chemical reaction rates of polyatomic species (and also atoms and diatoms as special cases) using variational transition state theory (VTST) with multidimensional tunneling. The original POLYRATE program was especially well suited for simple barrier reactions. But numerous reactions, especially in atmospheric and combustion chemistry, have no potential energy barrier between reactants and products. We have implemented a flexible VTST treatment for barrierless association reaction by employing a variable reaction coordinate and multifaceted dividing surfaces. The most time-consuming part of the calculation is an integral on eight- dimensional multifaceted hypersurface in the coordinate space of the interacting reagents; this hypersurface is conventionally called the dividing surface. The integral on the dividing surface is evaluated by the Monte Carlo method. The Monte Carlo integration requires thousands of samples of the energy to achieve a reasonably accurate result, and each sample requires a full electronic structure calculation. The electronic structure energies are calculated “on the fly” by interfacing the dynamics program with an electronic structure program; this procedure is called direct dynamics. The Monte Carlo integration is massively parallelized using the MPI message passing interface. By using an efficient strategy for achieving load balance and by using parallel I/O, we developed a new parallel version of POLYRATE that can be applied to large-scale problems on distributed-memory computers with almost linear scaling up to 150 processors and with 85% efficiency for 250 processors. The program has been well documented and is made available on the Internet for general use.