@@ -44,7 +44,7 @@ for i in 0..p0.n() {
         p0.0[i] = i as u64
     }
 
-    r.automorphism::<false>(&p0, 2 * r.n - 1, nth_root, &mut p1);
+    r.a_apply_automorphism_into_b::<false>(&p0, 2 * r.n - 1, nth_root, &mut p1);
 
     println!("{:?}", p1);
 }

@@ -1,6 +1,7 @@
-use crate::modulus::WordOps;
+use crate::modulus::{REDUCEMOD, WordOps};
 use crate::poly::Poly;
 use crate::ring::Ring;
+use crate::modulus::{ONCE, ScalarOperations};
 
 /// Returns a lookup table for the automorphism X^{i} -> X^{i * k mod nth_root}.
 /// Method will panic if n or nth_root are not power-of-two.
@@ -44,7 +45,41 @@ index
 }
 
 impl Ring<u64> {
-    pub fn automorphism<const NTT: bool>(
+
+    // b <- auto(a)
+    pub fn a_apply_automorphism_into_b<const NTT: bool>(
+        &self,
+        a: &Poly<u64>,
+        gal_el: usize,
+        nth_root: usize,
+        b: &mut Poly<u64>,
+    ){
+        self.apply_automorphism_core::<0, ONCE, NTT>(a, gal_el, nth_root, b)
+    }
+
+    // b <- REDUCEMOD(b + auto(a))
+    pub fn a_apply_automorphism_add_b_into_b<const REDUCE:REDUCEMOD, const NTT: bool>(
+        &self,
+        a: &Poly<u64>,
+        gal_el: usize,
+        nth_root: usize,
+        b: &mut Poly<u64>,
+    ){
+        self.apply_automorphism_core::<1, REDUCE, NTT>(a, gal_el, nth_root, b)
+    }
+
+    // b <- REDUCEMOD(b - auto(a))
+    pub fn a_apply_automorphism_sub_b_into_b<const REDUCE:REDUCEMOD, const NTT:bool>(
+    &self,
+    a: &Poly<u64>,
+    gal_el: usize,
+    nth_root: usize,
+    b: &mut Poly<u64>,
+    ){
+        self.apply_automorphism_core::<2, REDUCE, NTT>(a, gal_el, nth_root, b)
+    }
+
+    fn apply_automorphism_core<const MOD:u8, const REDUCE:REDUCEMOD, const NTT: bool>(
         &self,
         a: &Poly<u64>,
         gal_el: usize,
@@ -81,7 +116,12 @@ a_vec.iter().enumerate().for_each(|(i, ai)| {
                 let i_rev: usize = 2 * i.reverse_bits_msb(log_nth_root_half) + 1;
                 let gal_el_i: usize = (((gal_el * i_rev) & mask) - 1) >> 1;
                 let idx: usize = gal_el_i.reverse_bits_msb(log_nth_root_half);
-                b_vec[idx] = *ai;
+                match MOD{
+                    0 =>{b_vec[idx] = *ai}
+                    1=>{self.modulus.sa_add_sb_into_sb::<REDUCE>(ai, &mut b_vec[idx])}
+                    2=>{self.modulus.sa_sub_sb_into_sa::<1, REDUCE>(ai, &mut b_vec[idx])}
+                    _=>{panic!("invalid const MOD should be 0, 1, or 2 but is {}", MOD)}
+                }
             });
         } else {
             let n: usize = a.n();
@@ -92,12 +132,49 @@ a_vec.iter().enumerate().for_each(|(i, ai)| {
                 let gal_el_i: usize = i * gal_el;
                 let sign: u64 = ((gal_el_i >> log_n) & 1) as u64;
                 let i_out: usize = gal_el_i & mask;
-                b_vec[i_out] = ai * (sign ^ 1) | (q - ai) * sign
+                let v: u64 = ai * (sign ^ 1) | (q - ai) * sign;
+                match MOD{
+                    0 =>{b_vec[i_out] = v}
+                    1=>{self.modulus.sa_add_sb_into_sb::<REDUCE>(&v, &mut b_vec[i_out])}
+                    2=>{self.modulus.sa_sub_sb_into_sa::<1, REDUCE>(&v, &mut b_vec[i_out])}
+                    _=>{panic!("invalid const MOD should be 0, 1, or 2 but is {}", MOD)}
+                }
             });
         }
     }
 
-    pub fn automorphism_from_index<const NTT: bool>(
+    // b <- auto(a)
+    pub fn a_apply_automorphism_from_index_into_b<const NTT: bool>(
+        &self,
+        a: &Poly<u64>,
+        idx: &[usize],
+        b: &mut Poly<u64>,
+    ){
+        self.automorphism_from_index_core::<0, ONCE, NTT>(a, idx, b)
+    }
+
+    // b <- REDUCEMOD(b + auto(a))
+    pub fn a_apply_automorphism_from_index_add_b_into_b<const REDUCE:REDUCEMOD, const NTT: bool>(
+        &self,
+        a: &Poly<u64>,
+        idx: &[usize],
+        b: &mut Poly<u64>,
+    ){
+        self.automorphism_from_index_core::<1, REDUCE, NTT>(a, idx, b)
+    }
+
+    // b <- REDUCEMOD(b - auto(a))
+    pub fn a_apply_automorphism_from_index_sub_b_into_b<const REDUCE:REDUCEMOD, const NTT:bool>(
+    &self,
+    a: &Poly<u64>,
+    idx: &[usize],
+    b: &mut Poly<u64>,
+    ){
+        self.automorphism_from_index_core::<2, REDUCE, NTT>(a, idx, b)
+    }
+
+    // b <- auto(a) if OVERWRITE else b <- REDUCEMOD(b + auto(a))
+    fn automorphism_from_index_core<const MOD:u8, const REDUCE:REDUCEMOD, const NTT: bool>(
         &self,
         a: &Poly<u64>,
         idx: &[usize],
@@ -115,7 +192,12 @@ let a_vec: &Vec<u64> = &a.0;
 
         if NTT {
             a_vec.iter().enumerate().for_each(|(i, ai)| {
-                b_vec[idx[i]] = *ai;
+                match MOD{
+                    0 =>{ b_vec[idx[i]] = *ai}
+                    1 =>{self.modulus.sa_add_sb_into_sb::<REDUCE>(ai, &mut b_vec[idx[i]])}
+                    2=>{self.modulus.sa_sub_sb_into_sa::<1, REDUCE>(ai, &mut b_vec[idx[i]])}
+                    _=>{panic!("invalid const MOD should be 0, 1, or 2 but is {}", MOD)}
+                }
             });
         } else {
             let n: usize = a.n();
@@ -123,7 +205,13 @@ let mask: usize = n - 1;
             let q: u64 = self.modulus.q();
             a_vec.iter().enumerate().for_each(|(i, ai)| {
                 let sign: u64 = (idx[i] >> usize::BITS - 1) as u64;
-                b_vec[idx[i] & mask] = ai * (sign ^ 1) | (q - ai) * sign;
+                let v: u64 = ai * (sign ^ 1) | (q - ai) * sign;
+                match MOD{
+                    0 =>{b_vec[idx[i] & mask] = v}
+                    1 =>{self.modulus.sa_add_sb_into_sb::<REDUCE>(&v, &mut b_vec[idx[i] & mask])}
+                    2=>{self.modulus.sa_sub_sb_into_sa::<1, REDUCE>(&v, &mut b_vec[idx[i] & mask])}
+                    _=>{panic!("invalid const MOD should be 0, 1, or 2 but is {}", MOD)}
+                }
             });
         }
     }

@@ -1,6 +1,8 @@
 pub mod automorphism;
 pub mod packing;
+pub mod trace;
 pub mod rescaling_rns;
 pub mod ring;
 pub mod ring_rns;
 pub mod sampling;
+pub mod utils;

@@ -3,41 +3,8 @@ use crate::modulus::ONCE;
 use crate::poly::Poly;
 use crate::ring::Ring;
 use std::cmp::min;
-use std::collections::HashSet;
 
 impl Ring<u64> {
-    // Generates a vector storing {X^{2^0}, X^{2^1}, .., X^{2^log_n}}.
-    pub fn gen_x_pow_2<const NTT: bool, const INV: bool>(&self, log_n: usize) -> Vec<Poly<u64>> {
-        let mut x_pow: Vec<Poly<u64>> = Vec::<Poly<u64>>::with_capacity(log_n);
-
-        (0..log_n).for_each(|i| {
-            let mut idx: usize = 1 << i;
-
-            if INV {
-                idx = self.n() - idx;
-            }
-
-            x_pow.push(self.new_poly());
-
-            if i == 0 {
-                x_pow[i].0[idx] = self.modulus.montgomery.one();
-                self.ntt_inplace::<false>(&mut x_pow[i]);
-            } else {
-                let (left, right) = x_pow.split_at_mut(i);
-                self.a_mul_b_montgomery_into_c::<ONCE>(&left[i - 1], &left[i - 1], &mut right[0]);
-            }
-        });
-
-        if INV {
-            self.a_neg_into_a::<1, ONCE>(&mut x_pow[0]);
-        }
-
-        if !NTT {
-            x_pow.iter_mut().for_each(|x| self.intt_inplace::<false>(x));
-        }
-
-        x_pow
-    }
 
     pub fn pack<const ZEROGARBAGE: bool, const NTT: bool>(
         &self,
@@ -60,8 +27,6 @@ }
         });
 
         let gap: usize = max_gap(&indices);
-
-        let set: HashSet<_> = indices.into_iter().collect();
 
         if !ZEROGARBAGE {
             if gap > 0 {
@@ -79,7 +44,7 @@ .modulus
             .barrett
             .prepare(self.modulus.inv(1 << (log_end - log_start)));
 
-        set.iter().for_each(|i| {
+        indices.iter().for_each(|i| {
             if let Some(poly) = polys[*i].as_mut() {
                 if !NTT {
                     self.ntt_inplace::<true>(poly);
@@ -111,15 +76,14 @@ if let Some(poly_lo) = polys_lo[j].as_mut() {
                     let gal_el: usize = self.galois_element((1 << i) >> 1, i == 0, log_nth_root);
 
                     if !polys_hi[j].is_none() {
-                        self.automorphism::<true>(&tmpa, gal_el, 2 << self.log_n(), &mut tmpb);
-                        self.a_add_b_into_b::<ONCE>(&tmpb, poly_lo);
+                        self.a_apply_automorphism_add_b_into_b::<ONCE, true>(&tmpa, gal_el, 2 << self.log_n(), poly_lo);
                     } else {
-                        self.automorphism::<true>(poly_lo, gal_el, nth_root, &mut tmpa);
+                        self.a_apply_automorphism_into_b::<true>(poly_lo, gal_el, nth_root, &mut tmpa);
                         self.a_add_b_into_b::<ONCE>(&tmpa, poly_lo);
                     }
                 } else if let Some(poly_hi) = polys_hi[j].as_mut() {
                     let gal_el: usize = self.galois_element((1 << i) >> 1, i == 0, log_nth_root);
-                    self.automorphism::<true>(poly_hi, gal_el, nth_root, &mut tmpa);
+                    self.a_apply_automorphism_into_b::<true>(poly_hi, gal_el, nth_root, &mut tmpa);
                     self.a_sub_b_into_a::<1, ONCE>(&tmpa, poly_hi);
                     std::mem::swap(&mut polys_lo[j], &mut polys_hi[j]);
                 }
@@ -136,12 +100,14 @@ }
     }
 }
 
-// Returns the largest gap.
+
+// Returns the largest gap between two values in an ordered array of distinct values.
+// Panics if the array is not ordered or values are not distincts.
 fn max_gap(vec: &[usize]) -> usize {
     let mut gap: usize = usize::MAX;
     for i in 1..vec.len() {
         let (l, r) = (vec[i - 1], vec[i]);
-        assert!(r > l, "invalid input vec: not sorted");
+        assert!(r > l, "invalid input vec: not sorted or collision between indices");
         gap = min(gap, r - l);
         if gap == 1 {
             break;

@@ -0,0 +1,32 @@
+use crate::ring::Ring;
+use crate::poly::Poly;
+use crate::modulus::barrett::Barrett;
+use crate::modulus::ONCE;
+
+impl Ring<u64>{
+    pub fn trace_inplace<const NTT:bool>(&self, log_start: usize, log_end: usize, a: &mut Poly<u64>){
+        assert!(log_end <= self.log_n(), "invalid argument log_end: log_end={} > self.log_n()={}", log_end, self.log_n());
+        assert!(log_end > log_start, "invalid argument log_start: log_start={} > log_end={}", log_start, log_end);
+
+        let log_steps = log_end - log_start;
+
+        if log_steps > 0 {
+            let n_inv: Barrett<u64> = self.modulus.barrett.prepare(self.modulus.inv(1<<log_steps));
+            self.a_mul_b_scalar_barrett_into_a::<ONCE>(&n_inv, a);
+
+            if !NTT{
+                self.ntt_inplace::<false>(a);
+            }
+
+            let mut tmp: Poly<u64> = self.new_poly();
+
+            (log_start..log_end).for_each(|i|{
+
+            });
+        }
+
+        
+
+        
+    }
+}

@@ -0,0 +1,38 @@
+use crate::ring::Ring;
+use crate::poly::Poly;
+use crate::modulus::ONCE;
+
+impl Ring<u64>{
+    // Generates a vector storing {X^{2^0}, X^{2^1}, .., X^{2^log_n}}.
+    pub fn gen_x_pow_2<const NTT: bool, const INV: bool>(&self, log_n: usize) -> Vec<Poly<u64>> {
+        let mut x_pow: Vec<Poly<u64>> = Vec::<Poly<u64>>::with_capacity(log_n);
+
+        (0..log_n).for_each(|i| {
+            let mut idx: usize = 1 << i;
+
+            if INV {
+                idx = self.n() - idx;
+            }
+
+            x_pow.push(self.new_poly());
+
+            if i == 0 {
+                x_pow[i].0[idx] = self.modulus.montgomery.one();
+                self.ntt_inplace::<false>(&mut x_pow[i]);
+            } else {
+                let (left, right) = x_pow.split_at_mut(i);
+                self.a_mul_b_montgomery_into_c::<ONCE>(&left[i - 1], &left[i - 1], &mut right[0]);
+            }
+        });
+
+        if INV {
+            self.a_neg_into_a::<1, ONCE>(&mut x_pow[0]);
+        }
+
+        if !NTT {
+            x_pow.iter_mut().for_each(|x| self.intt_inplace::<false>(x));
+        }
+
+        x_pow
+    }
+}

@@ -39,7 +39,7 @@ if NTT {
         ring.ntt_inplace::<false>(&mut p0);
     }
 
-    ring.automorphism::<NTT>(&p0, 2 * n - 1, nth_root, &mut p1);
+    ring.a_apply_automorphism_into_b::<NTT>(&p0, 2 * n - 1, nth_root, &mut p1);
 
     if NTT {
         ring.intt_inplace::<false>(&mut p1);

@@ -1,63 +0,0 @@
-thread 'rustc' panicked at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\compiler\rustc_query_system\src\query\plumbing.rs:729:9:
-Found unstable fingerprints for thir_abstract_const(math[de55]::dft::ntt::{impl#2}::forward_inplace_core::{closure#2}): Ok(None)
-stack backtrace:
-   0:     0x7ff8de431d43 - std::backtrace_rs::backtrace::dbghelp64::trace
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\..\..\backtrace\src\backtrace\dbghelp64.rs:91
-   1:     0x7ff8de431d43 - std::backtrace_rs::backtrace::trace_unsynchronized
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\..\..\backtrace\src\backtrace\mod.rs:66
-   2:     0x7ff8de431d43 - std::backtrace::Backtrace::create
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\backtrace.rs:331
-   3:     0x7ff8de431c8a - std::backtrace::Backtrace::force_capture
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\backtrace.rs:312
-   4:     0x7ff8dfa55e22 - strncpy
-   5:     0x7ff8de44cb82 - alloc::boxed::impl$30::call
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/alloc\src\boxed.rs:1982
-   6:     0x7ff8de44cb82 - std::panicking::rust_panic_with_hook
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\panicking.rs:809
-   7:     0x7ff8de44c9c9 - std::panicking::begin_panic_handler::closure$0
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\panicking.rs:674
-   8:     0x7ff8de44a1cf - std::sys::backtrace::__rust_end_short_backtrace<std::panicking::begin_panic_handler::closure_env$0,never$>
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\sys\backtrace.rs:170
-   9:     0x7ff8de44c5ce - std::panicking::begin_panic_handler
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\panicking.rs:665
-  10:     0x7ff8e1457e31 - core::panicking::panic_fmt
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/core\src\panicking.rs:76
-  11:     0x7ff8e087eb42 - <rustc_ty_utils[611601a2e2e4762f]::opaque_types::OpaqueTypeCollector as rustc_type_ir[ed289d4b95bf9d77]::visit::TypeVisitor<rustc_middle[680c00e878606981]::ty::context::TyCtxt>>::visit_ty
-  12:     0x7ff8df148500 - rustc_ty_utils[611601a2e2e4762f]::ty::self_ty_of_trait_impl_enabling_order_dep_trait_object_hack
-  13:     0x7ff8df0b6a5f - rustc_ty_utils[611601a2e2e4762f]::ty::self_ty_of_trait_impl_enabling_order_dep_trait_object_hack
-  14:     0x7ff8de046cf9 - rustc_query_impl[5b7e8b0772f52dfc]::query_system
-  15:     0x7ff8de9b3ab6 - rustc_mir_build[c21120f9a6e0100a]::build::mir_build
-  16:     0x7ff8de899beb - rustc_mir_transform[79b20237fda78f87]::mir_built
-  17:     0x7ff8df1f5b8b - rustc_query_impl[5b7e8b0772f52dfc]::plumbing::query_key_hash_verify_all
-  18:     0x7ff8df13db6e - rustc_ty_utils[611601a2e2e4762f]::ty::self_ty_of_trait_impl_enabling_order_dep_trait_object_hack
-  19:     0x7ff8de038988 - rustc_ty_utils[611601a2e2e4762f]::ty::adt_sized_constraint
-  20:     0x7ff8e08f6c2f - <rustc_ty_utils[611601a2e2e4762f]::opaque_types::OpaqueTypeCollector as rustc_type_ir[ed289d4b95bf9d77]::visit::TypeVisitor<rustc_middle[680c00e878606981]::ty::context::TyCtxt>>::visit_ty
-  21:     0x7ff8df0ac364 - rustc_ty_utils[611601a2e2e4762f]::ty::self_ty_of_trait_impl_enabling_order_dep_trait_object_hack
-  22:     0x7ff8df0b3406 - rustc_ty_utils[611601a2e2e4762f]::ty::self_ty_of_trait_impl_enabling_order_dep_trait_object_hack
-  23:     0x7ff8df1fd6c7 - rustc_query_impl[5b7e8b0772f52dfc]::plumbing::query_key_hash_verify_all
-  24:     0x7ff8de4bf29e - rustc_interface[8d61ccc26bf4800c]::passes::resolver_for_lowering_raw
-  25:     0x7ff8db4854be - rustc_interface[8d61ccc26bf4800c]::passes::analysis
-  26:     0x7ff8de03dc5b - rustc_ty_utils[611601a2e2e4762f]::ty::adt_sized_constraint
-  27:     0x7ff8ddf512b7 - rustc_ty_utils[611601a2e2e4762f]::ty::adt_sized_constraint
-  28:     0x7ff8de043603 - rustc_query_impl[5b7e8b0772f52dfc]::query_system
-  29:     0x7ff8db443ce7 - _rust_alloc_error_handler
-  30:     0x7ff8db43f9bf - _rust_alloc_error_handler
-  31:     0x7ff8db4498d7 - _rust_alloc_error_handler
-  32:     0x7ff8de45e7cd - alloc::boxed::impl$28::call_once
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/alloc\src\boxed.rs:1968
-  33:     0x7ff8de45e7cd - alloc::boxed::impl$28::call_once
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/alloc\src\boxed.rs:1968
-  34:     0x7ff8de45e7cd - std::sys::pal::windows::thread::impl$0::new::thread_start
-                               at /rustc/5ec7d6eee7e0f5236ec1559499070eaf836bc608\library/std\src\sys\pal\windows\thread.rs:55
-  35:     0x7ff9be46e8d7 - BaseThreadInitThunk
-  36:     0x7ff9bec7fbcc - RtlUserThreadStart
-
-
-rustc version: 1.84.0-nightly (5ec7d6eee 2024-11-17)
-platform: x86_64-pc-windows-msvc
-
-query stack during panic:
-#0 [thir_abstract_const] building an abstract representation for `dft::ntt::<impl at math\src\dft\ntt.rs:91:1: 91:16>::forward_inplace_core::{closure#2}`
-#1 [mir_built] building MIR for `dft::ntt::<impl at math\src\dft\ntt.rs:91:1: 91:16>::forward_inplace_core::{closure#2}`
-#2 [analysis] running analysis passes on this crate
-end of query stack