We present a numerical framework for Fluorescence Diffuse Optical Tomography (fDOT) that combines a forward model together with an iterative reconstruction procedure. Using rapid linear solvers, we derived an efficient reconstruction strategy for quadratic regularizers. The method outperforms traditional reconstruction approaches. Starting from quadradic regularization, we then extend the framework to more general L-p constraints. We present reconstruction experiments that confirm the superiority of non-quadratic sparsity promoting regularization.