The simulation of the light flux density distribution on a receiver plays an important role in energy estimation, design and optimization of a heliostat field and the focusing strategy for a central receiver system (CRS). However, this simulation is a time-consuming procedure. In this paper, we propose a fast simulation method that fully exploits the tremendous rendering and parallel computing capacities of contemporary graphics processing units (GPUs). First, an auxiliary spatial data structure is employed to organize the heliostats in the field, and a parallel light beam traversal algorithm is designed and performed on the GPU to determine the shadowing and blocking heliostats for each reference heliostat. Then, the flux spot reflected by each heliostat on the receiver is computed using the HFLCAL model and accumulated for the final flux density distribution. Both the computing stage and accumulation stage are accomplished via GPU rendering pipeline. The proposed method is verified by taking the PS10 power plant as an example. Because this method considers both shadowing and blocking effects, the simulation results are consistent with those in the official report. Due to its high efficiency, the proposed method has potential applications in CRS design and optimization. (C) 2017 Elsevier Ltd. All rights reserved.