Exact simulations of huge systems of charged particles are impossible in practice, because their cost is proportional to N 2 , where N is the number of particles. We present an approximate and simple O(N log N) algorithm based upon the idea of recursive bisection of the set of particles. A small fraction of the particle-particle interactions is evaluated exactly by direct summation, while the rest is approximated via multipole expansions. Our algorithm is easy to program (in particular, in parallel, which is trivial), and it is well suited for systems with a non-uniform distribution of particles. For uniform systems, its accuracy and execution time are comparable to other fast methods such as tree codes and the fast multiple method.