Probability density approximation (PDA) is a nonparametric method of calculating probability densities. When integrated into Bayesian estimation, it allows researchers to fit psychological processes for which analytic probability functions are unavailable, significantly expanding the scope of theories that can be quantitatively tested. PDA is, however, computationally intensive, requiring large numbers of Monte Carlo simulations in order to attain good precision. We introduce Parallel PDA (pPDA), a highly efficient implementation of this method utilizing the Armadillo C++ and CUDA C libraries to conduct millions of model simulations simultaneously in graphics processing units (GPUs). This approach provides a practical solution for rapidly approximating probability densities with high precision. In addition to demonstrating this method, we fit a piecewise linear ballistic accumulator model (Holmes, Trueblood, & Heathcote, 2016) to empirical data. Finally, we conducted simulation studies to investigate various issues associated with PDA and provide guidelines for pPDA applications to other complex cognitive models.