poulpy/poulpy-backend/src/cpu_fft64_avx/znx_avx/automorphism.rs
poulpy-backend/src/cpu_fft64_avx/znx_avx/automorphism.rsRaw
1use core::arch::x86_64::*;
2
3#[inline]
4fn inv_mod_pow2(p: usize, bits: u32) -> usize {
5 // Compute p^{-1} mod 2^bits (p must be odd) through Hensel lifting.
6 debug_assert!(p % 2 == 1);
7 let mut x: usize = 1usize; // inverse mod 2
8 let mut i: u32 = 1;
9 while i < bits {
10 // x <- x * (2 - p*x) mod 2^(2^i) (wrapping arithmetic)
11 x = x.wrapping_mul(2usize.wrapping_sub(p.wrapping_mul(x)));
12 i <<= 1;
13 }
14 x & ((1usize << bits) - 1)
15}
16
17/// # Safety
18/// Caller must ensure the CPU supports AVX2 (e.g., via `is_x86_feature_detected!("avx2")`);
19/// all inputs must have the same length and must not alias.
20#[cfg(target_arch = "x86_64")]
21#[target_feature(enable = "avx2", enable = "fma")]
22pub fn znx_automorphism_avx(p: i64, res: &mut [i64], a: &[i64]) {
23 debug_assert_eq!(res.len(), a.len());
24 let n: usize = res.len();
25 if n == 0 {
26 return;
27 }
28 debug_assert!(n.is_power_of_two(), "n must be power of two");
29 debug_assert!(p & 1 == 1, "p must be odd (invertible mod 2n)");
30
31 if n < 4 {
32 use poulpy_hal::reference::znx::znx_automorphism_ref;
33
34 znx_automorphism_ref(p, res, a);
35 return;
36 }
37
38 unsafe {
39 let two_n: usize = n << 1;
40 let span: usize = n >> 2;
41 let bits: u32 = (two_n as u64).trailing_zeros();
42 let mask_2n: usize = two_n - 1;
43 let mask_1n: usize = n - 1;
44
45 // p mod 2n (positive)
46 let p_2n: usize = (((p & mask_2n as i64) + two_n as i64) as usize) & mask_2n;
47
48 // p^-1 mod 2n
49 let inv: usize = inv_mod_pow2(p_2n, bits);
50
51 // Broadcast constants
52 let n_minus1_vec: __m256i = _mm256_set1_epi64x((n as i64) - 1);
53 let mask_2n_vec: __m256i = _mm256_set1_epi64x(mask_2n as i64);
54 let mask_1n_vec: __m256i = _mm256_set1_epi64x(mask_1n as i64);
55
56 // Lane offsets [0, inv, 2*inv, 3*inv] (mod 2n)
57 let lane_offsets: __m256i = _mm256_set_epi64x(
58 ((inv * 3) & mask_2n) as i64,
59 ((inv * 2) & mask_2n) as i64,
60 inv as i64,
61 0i64,
62 );
63
64 // t_base = (j * inv) mod 2n.
65 let mut t_base: usize = 0;
66 let step: usize = (inv << 2) & mask_2n;
67
68 let mut rr: *mut __m256i = res.as_mut_ptr() as *mut __m256i;
69 let aa: *const i64 = a.as_ptr();
70
71 for _ in 0..span {
72 // t_vec = (t_base + [0, inv, 2*inv, 3*inv]) & (2n-1)
73 let t_base_vec: __m256i = _mm256_set1_epi64x(t_base as i64);
74 let t_vec: __m256i = _mm256_and_si256(_mm256_add_epi64(t_base_vec, lane_offsets), mask_2n_vec);
75
76 // idx = t_vec & (n-1)
77 let idx_vec: __m256i = _mm256_and_si256(t_vec, mask_1n_vec);
78
79 // sign = t >= n ? -1 : 0 (mask of all-ones where negate)
80 let sign_mask: __m256i = _mm256_cmpgt_epi64(t_vec, n_minus1_vec);
81
82 // gather a[idx] (scale = 8 bytes per i64)
83 let vals: __m256i = _mm256_i64gather_epi64(aa, idx_vec, 8);
84
85 // Conditional negate: (vals ^ sign_mask) - sign_mask
86 let vals_x: __m256i = _mm256_xor_si256(vals, sign_mask);
87 let out: __m256i = _mm256_sub_epi64(vals_x, sign_mask);
88
89 // store to res[j..j+4]
90 _mm256_storeu_si256(rr, out);
91
92 // advance
93 rr = rr.add(1);
94 t_base = (t_base + step) & mask_2n;
95 }
96 }
97}
98
99/// # Safety
100/// Caller must ensure the CPU supports AVX2 (e.g., via `is_x86_feature_detected!("avx2")`);
101/// all inputs must have the same length and must not alias.
102#[cfg(all(test, any(target_arch = "x86_64", target_arch = "x86")))]
103mod tests {
104 use poulpy_hal::reference::znx::znx_automorphism_ref;
105
106 use super::*;
107
108 #[target_feature(enable = "avx2", enable = "fma")]
109 fn test_znx_automorphism_internal() {
110 let a: [i64; 16] = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16];
111
112 let p: i64 = -5;
113
114 let mut r0: Vec<i64> = vec![0i64; a.len()];
115 let mut r1: Vec<i64> = vec![0i64; a.len()];
116
117 znx_automorphism_ref(p, &mut r0, &a);
118 znx_automorphism_avx(p, &mut r1, &a);
119
120 assert_eq!(r0, r1);
121 }
122
123 #[test]
124 fn test_znx_automorphism_avx() {
125 if !std::is_x86_feature_detected!("avx2") {
126 eprintln!("skipping: CPU lacks avx2");
127 return;
128 };
129 unsafe {
130 test_znx_automorphism_internal();
131 }
132 }
133}