Open up the Barrett number construction / debugging a bit.

This commit is contained in:
2018-11-06 21:54:32 -08:00
parent de5ff87f9e
commit 1a2f3aaa7a

View File

@@ -1,78 +1,98 @@
macro_rules! barrett_impl { macro_rules! barrett_impl {
($bar: ident, $name: ident, $name64: ident, $dbl: ident, $dbl64: ident) => { ($bar: ident, $name: ident, $name64: ident, $dbl: ident, $dbl64: ident) => {
pub struct $bar { #[derive(PartialEq)]
pub(crate) k: usize, pub struct $bar {
pub(crate) m: $name64, pub(crate) k: usize,
pub(crate) mu: $name64 pub(crate) m: $name64,
} pub(crate) mu: $name64
}
impl $bar { impl $bar {
/// Generate a new Barrett number from the given input. This operation /// Generate a new Barrett number from the given input. This operation
/// is relatively slow, so you should only do it if you plan to use /// is relatively slow, so you should only do it if you plan to use
/// reduce() multiple times with the ame number. /// reduce() multiple times with the ame number.
pub fn new(m: $name) -> $bar { pub fn new(m: $name) -> $bar {
// Step #1: Figure out k // Step #1: Figure out k
let mut k = 0; let mut k = 0;
for i in 0..m.value.len() { for i in 0..m.value.len() {
if m.value[i] != 0 { if m.value[i] != 0 {
k = i; k = i;
} }
} }
k += 1; k += 1;
// Step #2: Compute b // Step #2: Compute b
let mut b = $dbl64::zero(); let mut b = $dbl64::zero();
b.value[2*k] = 1; b.value[2*k] = 1;
// Step #3: Divide b by m. // Step #3: Divide b by m.
let bigm = $dbl64::from(&m); let bigm = $dbl64::from(&m);
let quot = b / &bigm; let quot = b / &bigm;
let resm = $name64::from(&m); let resm = $name64::from(&m);
let mu = $name64::from(&quot); let mu = $name64::from(&quot);
// Done! // Done!
$bar { k: k, m: resm, mu: mu } $bar { k: k, m: resm, mu: mu }
} }
// Reduce the input by this value; in other words, perform a mod /// Generate a new Barrett number from its component parts. You
// operation. /// probably don't want to use it; it's mostly here for debugging
pub fn reduce(&self, x: &$dbl) -> $name { /// purposes, so that tests won't have to compute this all the
// 1. q1←⌊x/bk1⌋, q2←q1 · μ, q3←⌊q2/bk+1⌋. /// time.
let q1: $name64 = $name64::from(x >> ((self.k - 1) * 64)); pub fn from_components(k: usize, m: $name64, mu: $name64) -> $bar {
let q2 = q1 * &self.mu; $bar { k: k, m: m, mu: mu }
let q3: $name64 = $name64::from(q2 >> ((self.k + 1) * 64)); }
// 2. r1←x mod bk+1, r2←q3 · m mod bk+1, r←r1 r2.
let mut r: $dbl64 = $dbl64::from(x);
r.mask(self.k + 1);
let mut r2: $dbl64 = $dbl64::from(q3 * &self.m);
r2.mask(self.k + 1);
let went_negative = &r < &r2;
r -= &r2;
// 3. If r<0 then r←r+bk+1.
if went_negative {
let mut bk1 = $dbl64::zero();
bk1.value[self.k+1] = 1;
// this may overflow, and we should probably be OK with it.
r += &bk1;
}
// 4. While r≥m do: r←rm.
let m2 = $dbl64::from(&self.m);
while &r > &m2 {
r -= &m2;
}
// Done!
$name::from(&r)
}
}
impl $name { // Reduce the input by this value; in other words, perform a mod
/// Generate a Barrett number from this number. // operation.
pub fn generate_barrett(&self) -> $bar { #[inline(always)]
$bar::new(self.clone()) pub fn reduce(&self, x: &$dbl) -> $name {
} // 1. q1←⌊x/bk1⌋, q2←q1 · μ, q3←⌊q2/bk+1⌋.
let q1: $name64 = $name64::from(x >> ((self.k - 1) * 64));
let q2 = q1 * &self.mu;
let q3: $name64 = $name64::from(q2 >> ((self.k + 1) * 64));
// 2. r1←x mod bk+1, r2←q3 · m mod bk+1, r←r1 r2.
let mut r: $dbl64 = $dbl64::from(x);
r.mask(self.k + 1);
let mut r2: $dbl64 = $dbl64::from(q3 * &self.m);
r2.mask(self.k + 1);
let went_negative = &r < &r2;
r -= &r2;
// 3. If r<0 then r←r+bk+1.
if went_negative {
let mut bk1 = $dbl64::zero();
bk1.value[self.k+1] = 1;
// this may overflow, and we should probably be OK with it.
r += &bk1;
}
// 4. While r≥m do: r←rm.
let m2 = $dbl64::from(&self.m);
while &r > &m2 {
r -= &m2;
}
// Done!
$name::from(&r)
}
}
#[cfg(test)] impl $name {
pub(crate) fn new_barrett(k: usize, m: $name64, mu: $name64) -> $bar { /// Generate a Barrett number from this number.
$bar{ k: k, m: m, mu: mu } pub fn generate_barrett(&self) -> $bar {
} $bar::new(self.clone())
} }
#[cfg(test)]
pub(crate) fn new_barrett(k: usize, m: $name64, mu: $name64) -> $bar {
$bar{ k: k, m: m, mu: mu }
}
}
impl fmt::Debug for $bar {
fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
f.debug_struct(stringify!($name))
.field("k", &self.k)
.field("m", &self.m)
.field("mu", &self.mu)
.finish()
}
}
}; };
} }