The Laplace-Cauchy problem of propagating Dirichlet and Neumann data from a portion to the rest of the boundary is an ill-posed inverse problem. Many regularizing algorithms have been recently proposed, in order to stabilize the solution with respect to noisy or incomplete data. Our main application is in electro-encephalography (EEG) where potential measurements available at part of the scalp are used to reconstruct the potential and the current on the inner skull surface. This problem, known as cortical mapping, and other applications --- in fields such as nondestructive testing, or biomedical engineering --- require to solve the problem in realistic, three-dimensional geometry. The goal of this article is to present a new boundary element based method for solving the Laplace-Cauchy problem in three dimensions, in a multilayer geometry. We validate the method experimentally on simulated data.