arnaucube
arnaucube/poulpyPublic Self-hosted repository
poulpy/poulpy-schemes/src/tfhe/blind_rotation/cggi_algo.rs
poulpy-schemes/src/tfhe/blind_rotation/cggi_algo.rsRaw
  1use itertools::izip;
  2use poulpy_hal::{
  3    api::{
  4        ScratchAvailable, SvpApplyDftToDft, SvpPPolAllocBytes, TakeVecZnx, TakeVecZnxBig, TakeVecZnxDft, TakeVecZnxDftSlice,
  5        TakeVecZnxSlice, VecZnxAddInplace, VecZnxBigAddSmallInplace, VecZnxBigAllocBytes, VecZnxBigNormalize,
  6        VecZnxBigNormalizeTmpBytes, VecZnxCopy, VecZnxDftAdd, VecZnxDftAddInplace, VecZnxDftAllocBytes, VecZnxDftApply,
  7        VecZnxDftSubInplace, VecZnxDftZero, VecZnxIdftApply, VecZnxIdftApplyConsume, VecZnxIdftApplyTmpBytes,
  8        VecZnxMulXpMinusOneInplace, VecZnxNormalize, VecZnxNormalizeInplace, VecZnxNormalizeTmpBytes, VecZnxRotate,
  9        VecZnxSubInplace, VmpApplyDftToDft, VmpApplyDftToDftAdd, VmpApplyDftToDftTmpBytes,
 10    },
 11    layouts::{Backend, DataMut, DataRef, Module, Scratch, SvpPPol, VecZnx, ZnxView, ZnxZero},
 12};
 13
 14use poulpy_core::{
 15    Distribution, GLWEOperations, TakeGLWECt,
 16    layouts::{GGSWInfos, GLWE, GLWEInfos, GLWEToMut, LWE, LWECiphertextToRef, LWEInfos},
 17};
 18
 19use crate::tfhe::blind_rotation::{
 20    BlincRotationExecute, BlindRotationKeyInfos, BlindRotationKeyPrepared, CGGI, LookUpTable, LookUpTableRotationDirection,
 21};
 22
 23#[allow(clippy::too_many_arguments)]
 24pub fn cggi_blind_rotate_scratch_space<B: Backend, OUT, GGSW>(
 25    module: &Module<B>,
 26    block_size: usize,
 27    extension_factor: usize,
 28    glwe_infos: &OUT,
 29    brk_infos: &GGSW,
 30) -> usize
 31where
 32    OUT: GLWEInfos,
 33    GGSW: GGSWInfos,
 34    Module<B>: VecZnxDftAllocBytes
 35        + VmpApplyDftToDftTmpBytes
 36        + VecZnxNormalizeTmpBytes
 37        + VecZnxBigAllocBytes
 38        + VecZnxIdftApplyTmpBytes
 39        + VecZnxBigNormalizeTmpBytes,
 40{
 41    let brk_size: usize = brk_infos.size();
 42
 43    if block_size > 1 {
 44        let cols: usize = (brk_infos.rank() + 1).into();
 45        let dnum: usize = brk_infos.dnum().into();
 46        let acc_dft: usize = module.vec_znx_dft_bytes_of(cols, dnum) * extension_factor;
 47        let acc_big: usize = module.vec_znx_big_bytes_of(1, brk_size);
 48        let vmp_res: usize = module.vec_znx_dft_bytes_of(cols, brk_size) * extension_factor;
 49        let vmp_xai: usize = module.vec_znx_dft_bytes_of(1, brk_size);
 50        let acc_dft_add: usize = vmp_res;
 51        let vmp: usize = module.vmp_apply_dft_to_dft_tmp_bytes(brk_size, dnum, dnum, 2, 2, brk_size); // GGSW product: (1 x 2) x (2 x 2)
 52        let acc: usize = if extension_factor > 1 {
 53            VecZnx::bytes_of(module.n(), cols, glwe_infos.size()) * extension_factor
 54        } else {
 55            0
 56        };
 57
 58        acc + acc_dft
 59            + acc_dft_add
 60            + vmp_res
 61            + vmp_xai
 62            + (vmp | (acc_big + (module.vec_znx_big_normalize_tmp_bytes() | module.vec_znx_idft_apply_tmp_bytes())))
 63    } else {
 64        GLWE::bytes_of(glwe_infos) + GLWE::external_product_inplace_scratch_space(module, glwe_infos, brk_infos)
 65    }
 66}
 67
 68impl<D: DataRef, B: Backend> BlincRotationExecute<B> for BlindRotationKeyPrepared<D, CGGI, B>
 69where
 70    Module<B>: VecZnxBigAllocBytes
 71        + VecZnxDftAllocBytes
 72        + SvpPPolAllocBytes
 73        + VmpApplyDftToDftTmpBytes
 74        + VecZnxBigNormalizeTmpBytes
 75        + VecZnxIdftApplyTmpBytes
 76        + VecZnxIdftApply<B>
 77        + VecZnxDftAdd<B>
 78        + VecZnxDftAddInplace<B>
 79        + VecZnxDftApply<B>
 80        + VecZnxDftZero<B>
 81        + SvpApplyDftToDft<B>
 82        + VecZnxDftSubInplace<B>
 83        + VecZnxBigAddSmallInplace<B>
 84        + VecZnxRotate
 85        + VecZnxAddInplace
 86        + VecZnxSubInplace
 87        + VecZnxNormalize<B>
 88        + VecZnxNormalizeInplace<B>
 89        + VecZnxCopy
 90        + VecZnxMulXpMinusOneInplace<B>
 91        + VmpApplyDftToDft<B>
 92        + VmpApplyDftToDftAdd<B>
 93        + VecZnxIdftApplyConsume<B>
 94        + VecZnxBigNormalize<B>
 95        + VecZnxNormalizeTmpBytes,
 96    Scratch<B>: TakeVecZnxDftSlice<B> + TakeVecZnxDft<B> + TakeVecZnxBig<B> + TakeVecZnxSlice + TakeVecZnx + ScratchAvailable,
 97{
 98    fn execute<DR: DataMut, DI: DataRef>(
 99        &self,
100        module: &Module<B>,
101        res: &mut GLWE<DR>,
102        lwe: &LWE<DI>,
103        lut: &LookUpTable,
104        scratch: &mut Scratch<B>,
105    ) {
106        match self.dist {
107            Distribution::BinaryBlock(_) | Distribution::BinaryFixed(_) | Distribution::BinaryProb(_) | Distribution::ZERO => {
108                if lut.extension_factor() > 1 {
109                    execute_block_binary_extended(module, res, lwe, lut, self, scratch)
110                } else if self.block_size() > 1 {
111                    execute_block_binary(module, res, lwe, lut, self, scratch);
112                } else {
113                    execute_standard(module, res, lwe, lut, self, scratch);
114                }
115            }
116            _ => panic!("invalid CGGI distribution"),
117        }
118    }
119}
120
121fn execute_block_binary_extended<DataRes, DataIn, DataBrk, B: Backend>(
122    module: &Module<B>,
123    res: &mut GLWE<DataRes>,
124    lwe: &LWE<DataIn>,
125    lut: &LookUpTable,
126    brk: &BlindRotationKeyPrepared<DataBrk, CGGI, B>,
127    scratch: &mut Scratch<B>,
128) where
129    DataRes: DataMut,
130    DataIn: DataRef,
131    DataBrk: DataRef,
132    Module<B>: VecZnxBigAllocBytes
133        + VecZnxDftAllocBytes
134        + SvpPPolAllocBytes
135        + VmpApplyDftToDftTmpBytes
136        + VecZnxBigNormalizeTmpBytes
137        + VecZnxIdftApplyTmpBytes
138        + VecZnxIdftApply<B>
139        + VecZnxDftAdd<B>
140        + VecZnxDftAddInplace<B>
141        + VecZnxDftApply<B>
142        + VecZnxDftZero<B>
143        + SvpApplyDftToDft<B>
144        + VecZnxDftSubInplace<B>
145        + VecZnxBigAddSmallInplace<B>
146        + VecZnxRotate
147        + VecZnxAddInplace
148        + VecZnxSubInplace
149        + VecZnxNormalize<B>
150        + VecZnxNormalizeInplace<B>
151        + VecZnxCopy
152        + VecZnxMulXpMinusOneInplace<B>
153        + VecZnxBigNormalize<B>
154        + VmpApplyDftToDft<B>,
155    Scratch<B>: TakeVecZnxDftSlice<B> + TakeVecZnxDft<B> + TakeVecZnxBig<B> + TakeVecZnxSlice + ScratchAvailable + TakeVecZnx,
156{
157    let n_glwe: usize = brk.n_glwe().into();
158    let extension_factor: usize = lut.extension_factor();
159    let base2k: usize = res.base2k().into();
160    let dnum: usize = brk.dnum().into();
161    let cols: usize = (res.rank() + 1).into();
162
163    let (mut acc, scratch_1) = scratch.take_vec_znx_slice(extension_factor, n_glwe, cols, res.size());
164    let (mut acc_dft, scratch_2) = scratch_1.take_vec_znx_dft_slice(extension_factor, n_glwe, cols, dnum);
165    let (mut vmp_res, scratch_3) = scratch_2.take_vec_znx_dft_slice(extension_factor, n_glwe, cols, brk.size());
166    let (mut acc_add_dft, scratch_4) = scratch_3.take_vec_znx_dft_slice(extension_factor, n_glwe, cols, brk.size());
167    let (mut vmp_xai, scratch_5) = scratch_4.take_vec_znx_dft(n_glwe, 1, brk.size());
168
169    (0..extension_factor).for_each(|i| {
170        acc[i].zero();
171    });
172
173    let x_pow_a: &Vec<SvpPPol<Vec<u8>, B>>;
174    if let Some(b) = &brk.x_pow_a {
175        x_pow_a = b
176    } else {
177        panic!("invalid key: x_pow_a has not been initialized")
178    }
179
180    let mut lwe_2n: Vec<i64> = vec![0i64; (lwe.n() + 1).as_usize()]; // TODO: from scratch space
181    let lwe_ref: LWE<&[u8]> = lwe.to_ref();
182
183    let two_n: usize = 2 * n_glwe;
184    let two_n_ext: usize = 2 * lut.domain_size();
185
186    mod_switch_2n(two_n_ext, &mut lwe_2n, &lwe_ref, lut.rotation_direction());
187
188    let a: &[i64] = &lwe_2n[1..];
189    let b_pos: usize = ((lwe_2n[0] + two_n_ext as i64) & (two_n_ext - 1) as i64) as usize;
190
191    let b_hi: usize = b_pos / extension_factor;
192    let b_lo: usize = b_pos & (extension_factor - 1);
193
194    for (i, j) in (0..b_lo).zip(extension_factor - b_lo..extension_factor) {
195        module.vec_znx_rotate(b_hi as i64 + 1, &mut acc[i], 0, &lut.data[j], 0);
196    }
197    for (i, j) in (b_lo..extension_factor).zip(0..extension_factor - b_lo) {
198        module.vec_znx_rotate(b_hi as i64, &mut acc[i], 0, &lut.data[j], 0);
199    }
200
201    let block_size: usize = brk.block_size();
202
203    izip!(
204        a.chunks_exact(block_size),
205        brk.data.chunks_exact(block_size)
206    )
207    .for_each(|(ai, ski)| {
208        (0..extension_factor).for_each(|i| {
209            (0..cols).for_each(|j| {
210                module.vec_znx_dft_apply(1, 0, &mut acc_dft[i], j, &acc[i], j);
211            });
212            module.vec_znx_dft_zero(&mut acc_add_dft[i])
213        });
214
215        // TODO: first & last iterations can be optimized
216        izip!(ai.iter(), ski.iter()).for_each(|(aii, skii)| {
217            let ai_pos: usize = ((aii + two_n_ext as i64) & (two_n_ext - 1) as i64) as usize;
218            let ai_hi: usize = ai_pos / extension_factor;
219            let ai_lo: usize = ai_pos & (extension_factor - 1);
220
221            // vmp_res = DFT(acc) * BRK[i]
222            (0..extension_factor).for_each(|i| {
223                module.vmp_apply_dft_to_dft(&mut vmp_res[i], &acc_dft[i], skii.data(), scratch_5);
224            });
225
226            // Trivial case: no rotation between polynomials, we can directly multiply with (X^{-ai} - 1)
227            if ai_lo == 0 {
228                // Sets acc_add_dft[i] = (acc[i] * sk) * X^{-ai} - (acc[i] * sk)
229                if ai_hi != 0 {
230                    // DFT X^{-ai}
231                    (0..extension_factor).for_each(|j| {
232                        (0..cols).for_each(|i| {
233                            module.svp_apply_dft_to_dft(&mut vmp_xai, 0, &x_pow_a[ai_hi], 0, &vmp_res[j], i);
234                            module.vec_znx_dft_add_inplace(&mut acc_add_dft[j], i, &vmp_xai, 0);
235                            module.vec_znx_dft_sub_inplace(&mut acc_add_dft[j], i, &vmp_res[j], i);
236                        });
237                    });
238                }
239
240            // Non trivial case: rotation between polynomials
241            // In this case we can't directly multiply with (X^{-ai} - 1) because of the
242            // ring homomorphism R^{N} -> prod R^{N/extension_factor}, so we split the
243            // computation in two steps: acc_add_dft = (acc * sk) * (-1) + (acc * sk) * X^{-ai}
244            } else {
245                // Sets acc_add_dft[0..ai_lo] += (acc[extension_factor - ai_lo..extension_factor] * sk) * X^{-ai+1}
246                if (ai_hi + 1) & (two_n - 1) != 0 {
247                    for (i, j) in (0..ai_lo).zip(extension_factor - ai_lo..extension_factor) {
248                        (0..cols).for_each(|k| {
249                            module.svp_apply_dft_to_dft(&mut vmp_xai, 0, &x_pow_a[ai_hi + 1], 0, &vmp_res[j], k);
250                            module.vec_znx_dft_add_inplace(&mut acc_add_dft[i], k, &vmp_xai, 0);
251                            module.vec_znx_dft_sub_inplace(&mut acc_add_dft[i], k, &vmp_res[i], k);
252                        });
253                    }
254                }
255
256                // Sets acc_add_dft[ai_lo..extension_factor] += (acc[0..extension_factor - ai_lo] * sk) * X^{-ai}
257                if ai_hi != 0 {
258                    // Sets acc_add_dft[ai_lo..extension_factor] += (acc[0..extension_factor - ai_lo] * sk) * X^{-ai}
259                    for (i, j) in (ai_lo..extension_factor).zip(0..extension_factor - ai_lo) {
260                        (0..cols).for_each(|k| {
261                            module.svp_apply_dft_to_dft(&mut vmp_xai, 0, &x_pow_a[ai_hi], 0, &vmp_res[j], k);
262                            module.vec_znx_dft_add_inplace(&mut acc_add_dft[i], k, &vmp_xai, 0);
263                            module.vec_znx_dft_sub_inplace(&mut acc_add_dft[i], k, &vmp_res[i], k);
264                        });
265                    }
266                }
267            }
268        });
269
270        {
271            let (mut acc_add_big, scratch7) = scratch_5.take_vec_znx_big(n_glwe, 1, brk.size());
272
273            (0..extension_factor).for_each(|j| {
274                (0..cols).for_each(|i| {
275                    module.vec_znx_idft_apply(&mut acc_add_big, 0, &acc_add_dft[j], i, scratch7);
276                    module.vec_znx_big_add_small_inplace(&mut acc_add_big, 0, &acc[j], i);
277                    module.vec_znx_big_normalize(base2k, &mut acc[j], i, base2k, &acc_add_big, 0, scratch7);
278                });
279            });
280        }
281    });
282
283    (0..cols).for_each(|i| {
284        module.vec_znx_copy(res.data_mut(), i, &acc[0], i);
285    });
286}
287
288fn execute_block_binary<DataRes, DataIn, DataBrk, B: Backend>(
289    module: &Module<B>,
290    res: &mut GLWE<DataRes>,
291    lwe: &LWE<DataIn>,
292    lut: &LookUpTable,
293    brk: &BlindRotationKeyPrepared<DataBrk, CGGI, B>,
294    scratch: &mut Scratch<B>,
295) where
296    DataRes: DataMut,
297    DataIn: DataRef,
298    DataBrk: DataRef,
299    Module<B>: VecZnxBigAllocBytes
300        + VecZnxDftAllocBytes
301        + SvpPPolAllocBytes
302        + VmpApplyDftToDftTmpBytes
303        + VecZnxBigNormalizeTmpBytes
304        + VecZnxIdftApplyTmpBytes
305        + VecZnxIdftApply<B>
306        + VecZnxDftAdd<B>
307        + VecZnxDftAddInplace<B>
308        + VecZnxDftApply<B>
309        + VecZnxDftZero<B>
310        + SvpApplyDftToDft<B>
311        + VecZnxDftSubInplace<B>
312        + VecZnxBigAddSmallInplace<B>
313        + VecZnxRotate
314        + VecZnxAddInplace
315        + VecZnxSubInplace
316        + VecZnxNormalize<B>
317        + VecZnxNormalizeInplace<B>
318        + VecZnxCopy
319        + VecZnxMulXpMinusOneInplace<B>
320        + VmpApplyDftToDft<B>
321        + VecZnxBigNormalize<B>,
322    Scratch<B>: TakeVecZnxDftSlice<B> + TakeVecZnxDft<B> + TakeVecZnxBig<B> + TakeVecZnxSlice + ScratchAvailable + TakeVecZnx,
323{
324    let n_glwe: usize = brk.n_glwe().into();
325    let mut lwe_2n: Vec<i64> = vec![0i64; (lwe.n() + 1).into()]; // TODO: from scratch space
326    let mut out_mut: GLWE<&mut [u8]> = res.to_mut();
327    let lwe_ref: LWE<&[u8]> = lwe.to_ref();
328    let two_n: usize = n_glwe << 1;
329    let base2k: usize = brk.base2k().into();
330    let dnum: usize = brk.dnum().into();
331
332    let cols: usize = (out_mut.rank() + 1).into();
333
334    mod_switch_2n(
335        2 * lut.domain_size(),
336        &mut lwe_2n,
337        &lwe_ref,
338        lut.rotation_direction(),
339    );
340
341    let a: &[i64] = &lwe_2n[1..];
342    let b: i64 = lwe_2n[0];
343
344    out_mut.data_mut().zero();
345
346    // Initialize out to X^{b} * LUT(X)
347    module.vec_znx_rotate(b, out_mut.data_mut(), 0, &lut.data[0], 0);
348
349    let block_size: usize = brk.block_size();
350
351    // ACC + [sum DFT(X^ai -1) * (DFT(ACC) x BRKi)]
352
353    let (mut acc_dft, scratch_1) = scratch.take_vec_znx_dft(n_glwe, cols, dnum);
354    let (mut vmp_res, scratch_2) = scratch_1.take_vec_znx_dft(n_glwe, cols, brk.size());
355    let (mut acc_add_dft, scratch_3) = scratch_2.take_vec_znx_dft(n_glwe, cols, brk.size());
356    let (mut vmp_xai, scratch_4) = scratch_3.take_vec_znx_dft(n_glwe, 1, brk.size());
357
358    let x_pow_a: &Vec<SvpPPol<Vec<u8>, B>>;
359    if let Some(b) = &brk.x_pow_a {
360        x_pow_a = b
361    } else {
362        panic!("invalid key: x_pow_a has not been initialized")
363    }
364
365    izip!(
366        a.chunks_exact(block_size),
367        brk.data.chunks_exact(block_size)
368    )
369    .for_each(|(ai, ski)| {
370        (0..cols).for_each(|j| {
371            module.vec_znx_dft_apply(1, 0, &mut acc_dft, j, out_mut.data_mut(), j);
372        });
373
374        module.vec_znx_dft_zero(&mut acc_add_dft);
375
376        izip!(ai.iter(), ski.iter()).for_each(|(aii, skii)| {
377            let ai_pos: usize = ((aii + two_n as i64) & (two_n - 1) as i64) as usize;
378
379            // vmp_res = DFT(acc) * BRK[i]
380            module.vmp_apply_dft_to_dft(&mut vmp_res, &acc_dft, skii.data(), scratch_4);
381
382            // DFT(X^ai -1) * (DFT(acc) * BRK[i])
383            (0..cols).for_each(|i| {
384                module.svp_apply_dft_to_dft(&mut vmp_xai, 0, &x_pow_a[ai_pos], 0, &vmp_res, i);
385                module.vec_znx_dft_add_inplace(&mut acc_add_dft, i, &vmp_xai, 0);
386                module.vec_znx_dft_sub_inplace(&mut acc_add_dft, i, &vmp_res, i);
387            });
388        });
389
390        {
391            let (mut acc_add_big, scratch_5) = scratch_4.take_vec_znx_big(n_glwe, 1, brk.size());
392
393            (0..cols).for_each(|i| {
394                module.vec_znx_idft_apply(&mut acc_add_big, 0, &acc_add_dft, i, scratch_5);
395                module.vec_znx_big_add_small_inplace(&mut acc_add_big, 0, out_mut.data_mut(), i);
396                module.vec_znx_big_normalize(
397                    base2k,
398                    out_mut.data_mut(),
399                    i,
400                    base2k,
401                    &acc_add_big,
402                    0,
403                    scratch_5,
404                );
405            });
406        }
407    });
408}
409
410fn execute_standard<DataRes, DataIn, DataBrk, B: Backend>(
411    module: &Module<B>,
412    res: &mut GLWE<DataRes>,
413    lwe: &LWE<DataIn>,
414    lut: &LookUpTable,
415    brk: &BlindRotationKeyPrepared<DataBrk, CGGI, B>,
416    scratch: &mut Scratch<B>,
417) where
418    DataRes: DataMut,
419    DataIn: DataRef,
420    DataBrk: DataRef,
421    Module<B>: VecZnxBigAllocBytes
422        + VecZnxDftAllocBytes
423        + SvpPPolAllocBytes
424        + VmpApplyDftToDftTmpBytes
425        + VecZnxBigNormalizeTmpBytes
426        + VecZnxIdftApplyTmpBytes
427        + VecZnxIdftApply<B>
428        + VecZnxDftAdd<B>
429        + VecZnxDftAddInplace<B>
430        + VecZnxDftApply<B>
431        + VecZnxDftZero<B>
432        + SvpApplyDftToDft<B>
433        + VecZnxDftSubInplace<B>
434        + VecZnxBigAddSmallInplace<B>
435        + VecZnxRotate
436        + VecZnxAddInplace
437        + VecZnxSubInplace
438        + VecZnxNormalize<B>
439        + VecZnxNormalizeInplace<B>
440        + VecZnxCopy
441        + VecZnxMulXpMinusOneInplace<B>
442        + VmpApplyDftToDft<B>
443        + VmpApplyDftToDftAdd<B>
444        + VecZnxIdftApplyConsume<B>
445        + VecZnxBigNormalize<B>
446        + VecZnxNormalizeTmpBytes,
447    Scratch<B>: TakeVecZnxDftSlice<B> + TakeVecZnxDft<B> + TakeVecZnxBig<B> + TakeVecZnxSlice + ScratchAvailable + TakeVecZnx,
448{
449    #[cfg(debug_assertions)]
450    {
451        assert_eq!(
452            res.n(),
453            brk.n(),
454            "res.n(): {} != brk.n(): {}",
455            res.n(),
456            brk.n()
457        );
458        assert_eq!(
459            lut.domain_size(),
460            brk.n_glwe().as_usize(),
461            "lut.n(): {} != brk.n(): {}",
462            lut.domain_size(),
463            brk.n_glwe().as_usize()
464        );
465        assert_eq!(
466            res.rank(),
467            brk.rank(),
468            "res.rank(): {} != brk.rank(): {}",
469            res.rank(),
470            brk.rank()
471        );
472        assert_eq!(
473            lwe.n(),
474            brk.n_lwe(),
475            "lwe.n(): {} != brk.data.len(): {}",
476            lwe.n(),
477            brk.n_lwe()
478        );
479    }
480
481    let mut lwe_2n: Vec<i64> = vec![0i64; (lwe.n() + 1).into()]; // TODO: from scratch space
482    let mut out_mut: GLWE<&mut [u8]> = res.to_mut();
483    let lwe_ref: LWE<&[u8]> = lwe.to_ref();
484
485    mod_switch_2n(
486        2 * lut.domain_size(),
487        &mut lwe_2n,
488        &lwe_ref,
489        lut.rotation_direction(),
490    );
491
492    let a: &[i64] = &lwe_2n[1..];
493    let b: i64 = lwe_2n[0];
494
495    out_mut.data_mut().zero();
496
497    // Initialize out to X^{b} * LUT(X)
498    module.vec_znx_rotate(b, out_mut.data_mut(), 0, &lut.data[0], 0);
499
500    // ACC + [sum DFT(X^ai -1) * (DFT(ACC) x BRKi)]
501    let (mut acc_tmp, scratch_1) = scratch.take_glwe_ct(&out_mut);
502
503    // TODO: see if faster by skipping normalization in external product and keeping acc in big coeffs
504    // TODO: first iteration can be optimized to be a gglwe product
505    izip!(a.iter(), brk.data.iter()).for_each(|(ai, ski)| {
506        // acc_tmp = sk[i] * acc
507        acc_tmp.external_product(module, &out_mut, ski, scratch_1);
508
509        // acc_tmp = (sk[i] * acc) * (X^{ai} - 1)
510        acc_tmp.mul_xp_minus_one_inplace(module, *ai, scratch_1);
511
512        // acc = acc + (sk[i] * acc) * (X^{ai} - 1)
513        out_mut.add_inplace(module, &acc_tmp);
514    });
515
516    // We can normalize only at the end because we add normalized values in [-2^{base2k-1}, 2^{base2k-1}]
517    // on top of each others, thus ~ 2^{63-base2k} additions are supported before overflow.
518    out_mut.normalize_inplace(module, scratch_1);
519}
520
521pub fn mod_switch_2n(n: usize, res: &mut [i64], lwe: &LWE<&[u8]>, rot_dir: LookUpTableRotationDirection) {
522    let base2k: usize = lwe.base2k().into();
523
524    let log2n: usize = usize::BITS as usize - (n - 1).leading_zeros() as usize + 1;
525
526    res.copy_from_slice(lwe.data().at(0, 0));
527
528    match rot_dir {
529        LookUpTableRotationDirection::Left => {
530            res.iter_mut().for_each(|x| *x = -*x);
531        }
532        LookUpTableRotationDirection::Right => {}
533    }
534
535    if base2k > log2n {
536        let diff: usize = base2k - (log2n - 1); // additional -1 because we map to [-N/2, N/2) instead of [0, N)
537        res.iter_mut().for_each(|x| {
538            *x = div_round_by_pow2(x, diff);
539        })
540    } else {
541        let rem: usize = base2k - (log2n % base2k);
542        let size: usize = log2n.div_ceil(base2k);
543        (1..size).for_each(|i| {
544            if i == size - 1 && rem != base2k {
545                let k_rem: usize = base2k - rem;
546                izip!(lwe.data().at(0, i).iter(), res.iter_mut()).for_each(|(x, y)| {
547                    *y = (*y << k_rem) + (x >> rem);
548                });
549            } else {
550                izip!(lwe.data().at(0, i).iter(), res.iter_mut()).for_each(|(x, y)| {
551                    *y = (*y << base2k) + x;
552                });
553            }
554        })
555    }
556}
557
558#[inline(always)]
559fn div_round_by_pow2(x: &i64, k: usize) -> i64 {
560    (x + (1 << (k - 1))) >> k
561}