Sampling from a Poisson distribution (infinite support) I would like to sample from an infinite discrete distribution, i.e. Poisson distribution
\begin{align}
\Bbb P(X=k)=e^{-\lambda}\frac{\lambda^k}{k!}
\end{align}
where $\lambda$ is a fixed parameter.
In case of a finite distribution I could use ALIAS method, but it doesn't work in this case. I've found that I could use a transformation
\begin{align}
X=\min\{ k: U \leq \sum\limits_{i=1}^n p_i \}
\end{align}
but is that the way to go?
 A: There is a quite surprising way to sample from the Poisson distribution with any $\lambda$ parameter.  The basis for this is the following algorithm to sample a Poisson variate with $\lambda = 1$, which involves only integer arithmetic and no floating-point operations (Duchon and Duvignau, "Preserving the number of cycles of length k in a growing uniform permutation", Electronic Journal of Combinatorics 23(4), 2016):

*

*Set ret to 1, a to 1, and b to 0.

*Generate $j$, a uniform random integer in [0, a].

*If $j<a$ and $j<b$, return ret.

*If $j=a$, add 1 to ret.  Otherwise, subtract 1 from ret and set b to $1+a$.

*Add 1 to a and go to step 1.

Poisson variates with other $\lambda$ parameters follow by splitting $\lambda$ into pieces of size no greater than 1.  Then, for each piece, a Poisson variate with $\lambda$ equal to that piece is generated.  Then all these variates are added together.
Finally, let $\mu < 1$, then to generate a Poisson variate with parameter $\mu$, it's enough to:

*

*Generate $N$, a Poisson variate with parameter 1.

*Flip $N$ coins that each show heads with probability $\mu$, then count the number of heads.

A: Yes, exactly. This is known as inverse transform sampling, and this earlier thread may be helpful.
Actually, you have a typo in your formula (note that there is only a single $k$ and a single $n$ in your formula). What you would do is to generate a value $p$ that is uniformly distributed in $[0,1]$, and then the resulting Poisson variable $k$ is the smallest $k$ such that $\sum_{i=0}^kp_i\geq p$. You would of course pre-tabulate your $p_i$ up to some sufficiently large number, so you don't need to recalculate them multiple times.
