Simplified spherical harmonics approximation (SPN) equations are widely used in modeling light propagation in biological tissues. However, with the increase of order N, its computational burden will severely aggravate. We propose a graphics processing unit (GPU) accelerated framework for SPN equations. Compared with the conventional central processing unit implementation, an increased performance of the GPU framework is obtained with an increase in mesh size, with the best speed-up ratio of 25 among the studied cases. The influence of thread distribution on the performance of the GPU framework is also investigated.