A new technique is presented for the solution of Poisson's equation in spherical coordinates. The method employs an expansion of the solution in a new set of functions defined herein for the first time, called 'spectral forms'. The spectral forms have spherical harmonics as their angular part, but use a new set of radial functions that automatically statisfy the boundary conditions, up to a multiplicative constant, on the Poisson solution. The resultant problem reduces to a set of simultaneous equations AC C ¼B B for the expansion coefficientsC C. The matrix A is block diagonal in the spherical harmonic indices l; m and is independent of any parameters. The simultaneous equations may be solved by LU decomposition. The LU decomposition only needs to be done once and multiple right hand sides (B-vectors) can be treated by a matrix-vector multiply. For a parallel computing platform, each such B-vector may be dealt with on a separate processor. Thus the algorithm is highly parallel. This technique may be used to calculate Coulomb energy integrals efficiently on a parallel computer.