Most of simulations often require the calculation of all pairwise interaction in large ensembles of particles, such as N-body problem of gravitation, electrostatic interaction and magnetic dipolar interaction, etc. The main difficulty in the calculation of long-range interaction is how to accelerate the slow convergence of the occurring sums. In this work, we are interested in the dipolar interaction in the two dimensional (2D) magnetic dipolar nanoparticle systems, which have attracted much attention due to both their important technological applications such as high-density patterned recording media and their rich and often unusual experimental behaviours. We develop a high efficiency algorithm based on the Lekner method to evaluate the magnetic dipolar energy for such systems, where the simulation cell is periodically replicated in the plane. Taking advantage of the symmetry of the systems, the dipolar interaction energy is expressed by rapidly converging series of modified Bessel functions in our algorithm. We found that our algorithm is better than the traditional Ewald summation method in efficiency for the regular arrays. Moreover, two simple formulas are obtained to evaluate the self-energy, which is important in the simulation of the dipolar systems.