This is my code for the Pareto/NBD: Problem with initialization - #6 by sonicking
The problem with your Pareto/NBD code above is that you are using pretty informative priors for r, alpha, s, and beta. But they are parameters of interest themselves and should come from less-informative priors. When you do, you will see Stan is unable to get convergence.
To answer your question, this means I would consider using lognormal distribution for lambda, mu and perhaps k, too.
My biggest issue with your approach for the Pareto/GGG is that it does not appear that one needs to evaluate that integral in each iteration. IIt should be performed only once when you are trying to compute P(alive).