We have developed a parallel CUDA implementation of the GLST algorithm for accelerated electrostatic interactions in MD simulations. We tested the algorithm on multiple salt water boxes of varying dimensions, and examined the dependence of its performance and accuracy on the target error threshold and the direct space cutoff. Linear scaling of the performance versus the system size and defined control of the RMSE have been demonstrated. The present implementation does not include periodic boundary condition, and further improvements are also needed. Nevertheless, given its highly parallelizable nature, we believe it will be a competitive algorithm to adopt especially for simulations of large systems running on computing clusters where multiple GPUs or a large numbers of CPU cores or threads can be utilized.