In boundary methods, piecewise particular solutions are employed to solve a given elliptic equation within subdomains of some region of interest. A boundary approximation is then obtained by satisfying the interior and exterior boundary conditions in a least squares sense. In this paper, we examine convergence, derive error norm bounds for approximate solutions and conduct a stability analysis of the associated algebraic problem. The aim of this analysis is to help choosing good partitions of subdomains. Finally, numerical experiments are carried out for a typical interface problem, demonstrating that very accurate solutions can be obtained while at the same time keeping small the condition numbers of the associated coefficient matrices.