We first characterize the increasing eigenfunctions associated to the following family of integro-differential operators, for any alpha, x > 0, gamma > 0 and f a smooth function on R+, L-(gamma) f(x) = x(-alpha) (sigma/2x(2)f ''(x) + (sigma gamma + b)xf'(x) + integral(infinity)(0) (f(e(-r) x) - f(x))e(-r gamma) + xf'(x)rI({r <= 1})nu(dr)), (0,1) where the coefficients b is an element of R, sigma >= 0 and the measure nu, which satisfies the integrability condition integral(infinity)(0) (1 boolean AND r(2))nu(dr) < +infinity, are uniquely determined by the distribution of a spectrally negative, infinitely divisible random variable, with characteristic exponent psi. L-(gamma) is known to be the infinitesimal generator of a positive alpha-self-similar Feller process, which has been introduced by Lamperti [Z. Wahrsch. Verw. Gebiete 22 (1972) 205-225]. The eigenfunctions are expressed in terms of a new family of power series which includes, for instance, the modified Bessel functions of the first kind and some generalizations of the Mittag-Leffler function. Then, we show that some specific combinations of these functions are Laplace transforms of self-decomposable or infinitely divisible distributions concentrated on the positive line with respect to the main argument, and, more surprisingly, with respect to the parameter psi(gamma). In particular, this generalizes a result of Hartman [Ann. Sc. Norm. Super Pisa Cl. Sci. IV-III (1976) 267-287] for the increasing solution of the Bessel differential equation. Finally, we compute, for some cases, the associated decreasing eigenfunctions and derive the Laplace transform of the exponential functionals of some spectrally negative Levy processes with a negative first moment.