the unit ball. The algorithm now repeatedly samples a random point $\vx$ from the unit ball. With probability $\exp(\lambda(M-f(\vx)))$, this point is accepted and we output it. Otherwise we continue.

Since $\exp(\lambda(M-f(\vx))) \in [\exp(-2\lambda), 1]$, this is a rejection sampling algorithm, and it outputs an exact sample from ${{\mathcal{D}}_f^{\lambda}}$. Each step accepts with probability at least $\exp(-\lambda)$. Thus the algorithm terminates in an expected $O(\exp(2\lambda))$ many steps, each of which requires one evaluation of $f$.

The above algorithm assumes access to an oracle to sample from $B_d(\mathbf{0}, 1)$ to arbitrary precision. This can be replaced by sampling from a grid of finite precision points in the ball. This creates two sources of error. Firstly, the function is not constant in the grid cell. This error is easily bounded since $f$ is Lipschitz. Secondly, some grid cells may cross the boundary of $B_d(\mathbf{0}, 1)$. This is a probability $d 2^{-b}$ event when sampling a grid point with $b$ bits of precision. Taking these errors into account gives us a sample within Wasserstein distance at most $O((d+\lambda) 2^{-b})$.

The above is a Las Vegas algorithm. One can similarly derive a Monte Carlo algorithm by aborting and outputting a random $\vx$ after $\exp(2\lambda) 