The problem with that is it will not be differentiable if you just set gammln(1) and especially gammln(2). There is more accurate code here. http://www.johndcook.com/cpp_gamma.html