COSAC: COmpact and Scalable Arbitrary-Centered Discrete Gaussian Sampling over Integers

Abstract

The arbitrary-centered discrete Gaussian sampler is a fundamental subroutine in implementing lattice trapdoor sampling algorithms. However, existing approaches typically rely on either a fast implementation of another discrete Gaussian sampler or pre-computations with regards to some specific discrete Gaussian distributions with fixed centers and standard deviations. These approaches may only support sampling from standard deviations within a limited range, or cannot efficiently sample from arbitrary standard deviations determined on-the-fly at run-time.

In this paper, we propose a compact and scalable rejection sampling algorithm by sampling from a continuous normal distribution and performing rejection sampling on rounded samples. Our scheme does not require pre-computations related to any specific discrete Gaussian distributions. Our scheme can sample from both arbitrary centers and arbitrary standard deviations determined on-the-fly at run-time. In addition, we show that our scheme only requires a low number of trials close to 2 per sample on average, and our scheme maintains good performance when scaling up the standard deviation. We also provide a concrete error analysis of our scheme based on the Rényi divergence. We implement our sampler and analyse its performance in terms of storage and speed compared to previous results. Our sampler’s running time is center-independent and is therefore applicable to implementation of convolution-style lattice trapdoor sampling and identity-based encryption resistant against timing side-channel attacks.

Appendices

A Precision Analysis

To avoid sampling a uniformly random real r with high absolute precisions at rejection steps 11 and 23 in Algorithm 2, and step 4 in Algorithm 3, we adapt the comparison approach similar to [26]. Assume an IEEE-754 floating-point value $f\in (0,1)$ with $(\delta _{f}+1)$-bit precision is represented by $f=\left( 1+mantissa\cdot 2^{-\delta _{f}}\right) \cdot 2^{exponent}$, where integer mantissa has $\delta _{f}$ bits and $exponent\in \mathbb {Z}^{-}$. To check $r<f$, one can sample $r_{m}\hookleftarrow \mathcal {U}\left( \{0,1\}^{\delta _{f}+1}\right)$, $r_{e}\hookleftarrow \mathcal {U}\left( \{0,1\}^{l}\right)$, and check $r_{m}<mantissa+2^{\delta _{f}}$ and $r_{e}<2^{l+exponent+1}$ instead for some l such that $l+exponent+1\ge 0$.

Here we analyse the precision requirement of $r_e$. We have the following theorem for the worst-case acceptance rate in Algorithm 2:

Theorem 5

Assume $x\in [-\tau ,\tau ]$ and $y\in [-\tau \sigma -1,\tau \sigma +1]$. In worst case, step 11 in Algorithm 2 has the acceptance rate:

$$\begin{aligned} p_1\ge \exp \left( -\frac{\left( -2\tau \sigma +c_{F}-3/2\right) \left( c_{F}-3/2\right) }{2\sigma ^2}\right) , \end{aligned}$$

and step 23 in Algorithm 2 has the acceptance rate:

$$\begin{aligned} p_2\ge \exp \left( -\frac{\left( 2\tau \sigma +c_{F}+3/2\right) \left( c_{F}+3/2\right) }{2\sigma ^2}\right) . \end{aligned}$$

Proof

For $b=0$ and $y\le -1/2$, we have the acceptance rate $p_1=\exp \left( -Y_1/\left( 2\sigma ^2\right) \right)$ at step 11 in Algorithm 2 where:

$$\begin{aligned} Y_1&=\left( \lfloor y\rceil +c_{F}\right) ^2-\left( y+1\right) ^2 \\&= \left( y+\delta +c_{F}\right) ^2-\left( y+1\right) ^2 \quad (\lfloor y\rceil =y+\delta \text { where } \delta \in [-1/2,1/2]) \\&= \left( 2y+\delta +c_{F}+1\right) \left( \delta +c_{F}-1\right) \\&\le \left( -2\tau \sigma +c_{F}-3/2\right) \left( c_{F}-3/2\right) . \quad \text {(when }\delta =-1/2 \text { and }y=-\tau \sigma -1) \end{aligned}$$

Similarly, for $b=1$ and $y\ge 1/2$, we have the acceptance rate $p_2=\exp \left( -Y_2/\left( 2\sigma ^2\right) \right)$ at step 23 in Algorithm 2 where:

$$\begin{aligned} Y_2&=\left( \lfloor y\rceil +c_{F}\right) ^2-\left( y-1\right) ^2 \\&= \left( y+\delta +c_{F}\right) ^2-\left( y-1\right) ^2 \quad (\lfloor y\rceil =y+\delta \text { where }\delta \in [-1/2,1/2]) \\&= \left( 2y+\delta +c_{F}-1\right) \left( \delta +c_{F}+1\right) \\&\le \left( 2\tau \sigma +c_{F}+3/2\right) \left( c_{F}+3/2\right) . \quad \text {(when }\delta =1/2\text { and }y=\tau \sigma +1) \end{aligned}$$

   $\square$

Let $\varDelta \le 1/2$ be the maximum relative error of the right hand side computations at rejection steps 11 and 23 in Algorithm 2, and step 4 in Algorithm 3. For $\exp \left( -Y_1/\left( 2\sigma ^2\right) \right)$ at step 11 in Algorithm 2, we have:

$$\begin{aligned}&exponent_{1} \ge \left\lfloor \log _2\left( (1-\varDelta )\cdot \exp \left( -\frac{Y_1}{2\sigma ^2}\right) \right) \right\rfloor \\&\ge \left\lfloor -1-\frac{\left( -2\tau \sigma +c_{F}-3/2\right) \left( c_{F}-3/2\right) }{2\sigma ^2}\cdot \log _{2}e\right\rfloor \quad \text {(by Thm. 5 and }\varDelta \le 1/2) \\&\ge \left\lfloor -1-\frac{2\tau \sigma +2}{\sigma ^2}\cdot \log _{2}e\right\rfloor . \quad \text {(when }c_{F}=-1/2) \end{aligned}$$

Similarly, for $\exp \left( -Y_2/\left( 2\sigma ^2\right) \right)$ at step 23 in Algorithm 2, we have:

$$\begin{aligned}&exponent_{2} \ge \left\lfloor \log _2\left( (1-\varDelta )\cdot \exp \left( -\frac{Y_2}{2\sigma ^2}\right) \right) \right\rfloor \\&\ge \left\lfloor -1-\frac{\left( 2\tau \sigma +c_{F}+3/2\right) \left( c_{F}+3/2\right) }{2\sigma ^2}\cdot \log _{2}e\right\rfloor \quad \text {(by Thm. 5 and }\varDelta \le 1/2) \\&\ge \left\lfloor -1-\frac{2\tau \sigma +2}{\sigma ^2}\cdot \log _{2}e\right\rfloor . \quad \text {(when }c_{F}=1/2) \end{aligned}$$

For $\exp \left( -c_{F}^2/\left( 2\sigma ^2\right) \right) /S$ at step 4 in Algorithm 3, we have:

$$\begin{aligned}&exponent_{3} \ge \left\lfloor \log _2\left( (1-\varDelta )\cdot \exp \left( -\frac{c_{F}^2}{2\sigma ^2}\right) /S\right) \right\rfloor \\&\ge \left\lfloor -1-\frac{1}{8\sigma ^2}\cdot \log _{2}e-\log _{2}\left( \sigma \sqrt{2\pi }\right) \right\rfloor . \quad \text {(when }c_{F}=\pm 1/2 \text { and }\varDelta \le 1/2) \end{aligned}$$

Therefore, we have:

$$\begin{aligned} exponent\ge \min \left\{ \left\lfloor -1-\frac{2\tau \sigma +2}{\sigma ^2}\cdot \log _{2}e\right\rfloor , \left\lfloor -1-\frac{1}{8\sigma ^2}\cdot \log _{2}e-\log _{2}\left( \sigma \sqrt{2\pi }\right) \right\rfloor \right\} . \end{aligned}$$

Since the probability $\mathrm {Pr}\left[ -\tau \le x\le \tau \right] =\text {erf}\left( \tau /\sqrt{2}\right)$ for $x\hookleftarrow \mathcal {N}(0,1)$, to ensure $1-\mathrm {Pr}\left[ -\tau \le x\le \tau \right] \le 2^{-\lambda }$, we need $\tau \ge \sqrt{2}\cdot \text {erf}^{-1}\left( 1-2^{-\lambda }\right)$. Therefore, for $\lambda =128$ and $\sigma \in \left[ 2,2^{20}\right]$, we have $\tau \ge 13.11$, $exponent\ge -23$, and thus $l\ge 22$, i.e. $r_e$ needs to have at least 22 bits.

B Proof of Algorithm 1

Since Algorithm 1 was an exercise in [4] without solutions, here we provide a brief proof of Algorithm 1.

Normalisation Factor. By definition, we have the normalisation factor:

$$\begin{aligned} \frac{1}{c}&=\sum _{k\in \mathbb {Z}}\exp \left( -\frac{\left( |k|+1/2\right) ^2}{2\sigma ^2}\right) \\&=\sum _{k\in \mathbb {Z}^-}\exp \left( -\frac{\left( k-1/2\right) ^2}{2\sigma ^2}\right) +\exp \left( -\frac{1}{8\sigma ^2}\right) +\sum _{k\in \mathbb {Z}^+}\exp \left( -\frac{\left( k+1/2\right) ^2}{2\sigma ^2}\right) \\&=\rho _{1/2,\sigma }\left( \mathbb {Z}^-\right) +\rho _{-1/2,\sigma }\left( \mathbb {Z}^+\right) +\exp \left( -\frac{1}{8\sigma ^2}\right) \\&\ge \frac{1-\epsilon }{1+\epsilon }\cdot \rho _{\sigma }\left( \mathbb {Z}\right) +\exp \left( -\frac{1}{8\sigma ^2}\right) -2. \quad \text {(By Lemma 1)} \end{aligned}$$

Correctness. Let $z_0=\lfloor x\rceil$. We have $Y=\left( \left| z_0\right| +1/2\right) ^2-x^2\ge 0$ for any $x\in \mathbb {R}$. Therefore, the rejection condition $\exp \left( -Y/\left( 2\sigma ^2\right) \right) \in \left( 0,1\right]$. We have the output distribution:

$$\begin{aligned} \mathrm {Pr}\left[ z=z_0\right]&\propto \int _{z_0-1/2}^{z_0+1/2} \exp \left( -\frac{x^2}{2\sigma ^2}\right) \cdot \exp \left( -\frac{\left( \left| z_0\right| +1/2\right) ^2-x^2}{2\sigma ^2}\right) \,\mathrm {d}x \\&=\int _{z_0-1/2}^{z_0+1/2} \exp \left( -\frac{\left( \left| z_0\right| +1/2\right) ^2}{2\sigma ^2}\right) \,\mathrm {d}x=\exp \left( -\frac{\left( \left| z_0\right| +1/2\right) ^2}{2\sigma ^2}\right) . \end{aligned}$$

Rejection Rate. By definition, we have the output probability density function $f(x)=c\cdot \exp \left( -\frac{\left( |\lfloor x\rceil |+1/2\right) ^2}{2\sigma ^2}\right)$ and the input probability density function $g(x)=\rho _{\sigma }\left( x\right) /\left( \sigma \sqrt{2\pi }\right)$. The expected number of trials can be written as:

$$\begin{aligned} M=\max \frac{f(x)}{g(x)}=\max \left( \frac{\exp \left( -\frac{\left( |\lfloor x\rceil |+1/2\right) ^2}{2\sigma ^2}\right) }{\rho _{\sigma }\left( x\right) }\cdot \frac{\sigma \sqrt{2\pi }}{1/c}\right) . \end{aligned}$$

We have:

$$\begin{aligned} \frac{\exp \left( -\frac{\left( |\lfloor x\rceil |+1/2\right) ^2}{2\sigma ^2}\right) }{\rho _{\sigma }\left( x\right) }=\frac{\exp \left( -\frac{\left( |\lfloor x\rceil |+1/2\right) ^2}{2\sigma ^2}\right) }{\exp \left( -\frac{x^2}{2\sigma ^2}\right) }=\exp \left( -\frac{\left( |\lfloor x\rceil |+1/2\right) ^2-x^2}{2\sigma ^2}\right) \le 1. \end{aligned}$$

Therefore,

$$\begin{aligned} M\le \frac{\sigma \sqrt{2\pi }}{1/c}&\le \frac{\sigma \sqrt{2\pi }}{\frac{1-\epsilon }{1+\epsilon }\cdot \rho _{\sigma }\left( \mathbb {Z}\right) +\exp \left( -\frac{1}{8\sigma ^2}\right) -2} \\&\le \frac{\sigma \sqrt{2\pi }}{\frac{1-\epsilon }{1+\epsilon }\cdot \left( \sigma \sqrt{2\pi }-1\right) +\exp \left( -\frac{1}{8\sigma ^2}\right) -2}, \end{aligned}$$

where the second inequality follows from the inequality of 1/c and the third inequality follows from the fact that $\rho _{\sigma }\left( \mathbb {Z}\right) \ge \sigma \sqrt{2\pi }-1$. Thus, we have the rejection probability:

$$\begin{aligned} 1-\frac{1}{M}\le \frac{\left( 1-\frac{1-\epsilon }{1+\epsilon }\right) \cdot \sigma \sqrt{2\pi }+\frac{1-\epsilon }{1+\epsilon }-\exp \left( -\frac{1}{8\sigma ^2}\right) +2}{\sigma \sqrt{2\pi }} \approx \frac{3-\exp \left( -\frac{1}{8\sigma ^2}\right) }{\sigma \sqrt{2\pi }}\le \frac{2}{\sigma }\sqrt{\frac{2}{\pi }}, \end{aligned}$$

when $\epsilon$ is small.