We present a numerical method for the study of correlated quantum impurity problems out of equilibrium, which is particularly suited to address steady state properties within Dynamical Mean Field Theory. The approach, recently introduced in [Arrigoni et al., Phys. Rev. Lett. 110, 086403 (2013)], is based upon a mapping of the original impurity problem onto an auxiliary open quantum system, consisting of the interacting impurity coupled to bath sites as well as to a Markovian environment. The dynamics of the auxiliary system is governed by a Lindblad master equation whose parameters are used to optimize the mapping. The accuracy of the results can be readily estimated and systematically improved by increasing the number of auxiliary bath sites, or by introducing a linear correction. Here, we focus on a detailed discussion of the proposed approach including technical remarks. To solve for the Green's functions of the auxiliary impurity problem, a non-hermitian Lanczos diagonalization is applied. As a benchmark, results for the steady state current-voltage characteristics of the single impurity Anderson model are presented. Furthermore, the bias dependence of the single particle spectral function and the splitting of the Kondo resonance are discussed. In its present form the method is fast, efficient and features a controlled accuracy.

17 pages, 6 figures