We present a survey of a direct method of boundary integral equations for the numerical solution of the Cauchy problem for the Laplace equation in doubly connected domains. The domain of solution is located between two closed boundary surfaces (curves in the case of two-dimensional domains). This Cauchy problem is reduced to finding the values of a harmonic function and its normal derivative on one of the two closed parts of the boundary according to the information about these quantities on the other boundary surface. This is an ill-posed problem in which the presence of noise in the input data may completely destroy the procedure of finding the approximate solution. We describe and present the results for a procedure of regularization aimed at the stable determination of the required quantities based on the representation of the solution to the Cauchy problem in the form a single-layer potential. For given data, this representation yields a system of boundary integral equations with two unknown densities. We establish the existence and uniqueness of these densities and propose a method for the numerical discretization in two- and three-dimensional domains. We also consider the cases of simply connected domains of the solution and unbounded domains. Numerical examples are presented both for two- and three-dimensional domains. These numerical results demonstrate that the proposed method gives good accuracy with relatively small amount of computations.