If you look at http://users.wfu.edu/natalie/s11phy752/lecturenote/Ewald.pdf, you note that

   * the 1=erf + erfc trick works also for a 2-dimension lattice since the arguments depend only on the distance of the lattice points.

   * the representation of the three-dimensional delta function by a sum of plane waves over the three-dimensional reciprocal lattice can easily be modified to the 2D case

   * instead of 1/|\vec r_i -\vec r_j| one can in principle use any pair potential v(\vec r_i -\vec r_j)

So, I do not see why the Ewald method should not work for your problem. Of course, you have to work out the details.

http://users.wfu.edu/natalie/s11phy752/lecturenote/Ewald.pdf

@Herbert H H Homeier, i am currently not experted in this domain. But i believe that the Ewald summation which works with the case of 3D systems with 3D periodicity don't work with those with 2D periodicity. I believe in that because i found when searching some papers which aims to generalize  the ES to systems with p dimensional periodicity. There is much theoretical background there which i currently don't master so i need someone who is experted in this domain to explain for me why the standard ES don't work in 2D and what to do. You can take a look at this two papers below

See also RHHR method (third file)

@Rami Zouari Actually the first of the .pdf is _based_ on Ewald partitioning, and it explicitly states that the method is exponentially converging. That it can be improved by FFT is nothing that says that Ewald summation is not working in 2d conditions.  Perhaps you should contact the authors of the preprint directly  :)

The explicit expression for the dipole-dipole interaction energy in 2D in terms of the Riemann zeta function is presented in Eq. (5.8), on p. 155 (Sec. 5.1.1), of the book Two-Dimensional Crystals, by I Lyuksyutov, AG Naumovets, and V Pokrovsky (Academic Press, New York, 1992). The double sum in this expression is fast-converging. In practice, one can use a variety of transformation methods, such as that by Daniel Shanks (J. Math. Phys. 34, 1 (1955)), for generating a fast-converging sequence of partial sums for this sum.*)

Two remarks are in order. Firstly, the Ewald summation technique is used for sums that are conditionally convergent, such as that for the electrostatic energy of a periodic lattice of charged ions. The sum associated with the dipole-dipole interaction energy is absolutely convergent. Secondly, the dipole-dipole interaction potential in 3D decays like 1/r3, where r denotes the distance between two dipoles (for this, see Eqs. (5.2) and (5.3) in the book Theory of Magnetism I, 2nd ed, by DC Mattis (Springer, Berlin, 1988)). The power 3 is often used also in 2D, as evidenced by the value s = 3/2 in the above-mentioned book by Lyuksyutov et al. (see remark following Eq. (5.6) herein). Strictly, since the Green function of the Laplace operator in 2D is proportional to ln(r), and not to 1/r, which is specific to 3D, the dipole-dipole interaction potential in 2D decays like 1/r2, implying s = 1. Use of s = 3/2 is based on the consideration that whereas dipoles are confined to a two-dimensional plane, their interaction follows the electrostatic laws in 3D. Be it as it may, in performing your calculations you will have to ascertain the appropriate value for s before using the above-mentioned expression for the dipole-dipole interaction energy. □

__

*) The routine C06BAF of the NAG Library makes use of this transformation.

Thank you Dr Farid for your response. I will read the references that you gave to me and i will work in all the details inside.