Add a generic power function

The patch adds a `pow` function for types implementing `One`, `Mul` and
`Clone` trait.

The patch also renames f32 and f64 pow into powf in order to still have
a way to easily have float powers. It uses llvms intrinsics.

The pow implementation for all num types uses the exponentiation by
square.

Fixes bug #11499
This commit is contained in:
Flavio Percoco 2014-01-12 21:02:59 +01:00
parent 5fdc81262a
commit ed7e576d9c
10 changed files with 82 additions and 54 deletions

View File

@ -290,7 +290,7 @@ be distributed on the available cores.
fn partial_sum(start: uint) -> f64 { fn partial_sum(start: uint) -> f64 {
let mut local_sum = 0f64; let mut local_sum = 0f64;
for num in range(start*100000, (start+1)*100000) { for num in range(start*100000, (start+1)*100000) {
local_sum += (num as f64 + 1.0).pow(&-2.0); local_sum += (num as f64 + 1.0).powf(&-2.0);
} }
local_sum local_sum
} }
@ -326,7 +326,7 @@ a single large vector of floats. Each task needs the full vector to perform its
use extra::arc::Arc; use extra::arc::Arc;
fn pnorm(nums: &~[f64], p: uint) -> f64 { fn pnorm(nums: &~[f64], p: uint) -> f64 {
nums.iter().fold(0.0, |a,b| a+(*b).pow(&(p as f64)) ).pow(&(1.0 / (p as f64))) nums.iter().fold(0.0, |a,b| a+(*b).powf(&(p as f64)) ).powf(&(1.0 / (p as f64)))
} }
fn main() { fn main() {

View File

@ -653,19 +653,19 @@ fn test<T: Float>(given: T, (numer, denom): (&str, &str)) {
// f32 // f32
test(3.14159265359f32, ("13176795", "4194304")); test(3.14159265359f32, ("13176795", "4194304"));
test(2f32.pow(&100.), ("1267650600228229401496703205376", "1")); test(2f32.powf(&100.), ("1267650600228229401496703205376", "1"));
test(-2f32.pow(&100.), ("-1267650600228229401496703205376", "1")); test(-2f32.powf(&100.), ("-1267650600228229401496703205376", "1"));
test(1.0 / 2f32.pow(&100.), ("1", "1267650600228229401496703205376")); test(1.0 / 2f32.powf(&100.), ("1", "1267650600228229401496703205376"));
test(684729.48391f32, ("1369459", "2")); test(684729.48391f32, ("1369459", "2"));
test(-8573.5918555f32, ("-4389679", "512")); test(-8573.5918555f32, ("-4389679", "512"));
// f64 // f64
test(3.14159265359f64, ("3537118876014453", "1125899906842624")); test(3.14159265359f64, ("3537118876014453", "1125899906842624"));
test(2f64.pow(&100.), ("1267650600228229401496703205376", "1")); test(2f64.powf(&100.), ("1267650600228229401496703205376", "1"));
test(-2f64.pow(&100.), ("-1267650600228229401496703205376", "1")); test(-2f64.powf(&100.), ("-1267650600228229401496703205376", "1"));
test(684729.48391f64, ("367611342500051", "536870912")); test(684729.48391f64, ("367611342500051", "536870912"));
test(-8573.5918555, ("-4713381968463931", "549755813888")); test(-8573.5918555, ("-4713381968463931", "549755813888"));
test(1.0 / 2f64.pow(&100.), ("1", "1267650600228229401496703205376")); test(1.0 / 2f64.powf(&100.), ("1", "1267650600228229401496703205376"));
} }
#[test] #[test]

View File

@ -349,8 +349,8 @@ pub fn write_boxplot(w: &mut io::Writer, s: &Summary, width_hint: uint) {
let (q1,q2,q3) = s.quartiles; let (q1,q2,q3) = s.quartiles;
// the .abs() handles the case where numbers are negative // the .abs() handles the case where numbers are negative
let lomag = (10.0_f64).pow(&(s.min.abs().log10().floor())); let lomag = (10.0_f64).powf(&(s.min.abs().log10().floor()));
let himag = (10.0_f64).pow(&(s.max.abs().log10().floor())); let himag = (10.0_f64).powf(&(s.max.abs().log10().floor()));
// need to consider when the limit is zero // need to consider when the limit is zero
let lo = if lomag == 0.0 { let lo = if lomag == 0.0 {

View File

@ -409,7 +409,7 @@ fn ln_10() -> f32 { 2.30258509299404568401799145468436421 }
fn recip(&self) -> f32 { 1.0 / *self } fn recip(&self) -> f32 { 1.0 / *self }
#[inline] #[inline]
fn pow(&self, n: &f32) -> f32 { pow(*self, *n) } fn powf(&self, n: &f32) -> f32 { pow(*self, *n) }
#[inline] #[inline]
fn sqrt(&self) -> f32 { sqrt(*self) } fn sqrt(&self) -> f32 { sqrt(*self) }
@ -1265,7 +1265,7 @@ fn test_frexp_nowin() {
fn test_integer_decode() { fn test_integer_decode() {
assert_eq!(3.14159265359f32.integer_decode(), (13176795u64, -22i16, 1i8)); assert_eq!(3.14159265359f32.integer_decode(), (13176795u64, -22i16, 1i8));
assert_eq!((-8573.5918555f32).integer_decode(), (8779358u64, -10i16, -1i8)); assert_eq!((-8573.5918555f32).integer_decode(), (8779358u64, -10i16, -1i8));
assert_eq!(2f32.pow(&100.0).integer_decode(), (8388608u64, 77i16, 1i8)); assert_eq!(2f32.powf(&100.0).integer_decode(), (8388608u64, 77i16, 1i8));
assert_eq!(0f32.integer_decode(), (0u64, -150i16, 1i8)); assert_eq!(0f32.integer_decode(), (0u64, -150i16, 1i8));
assert_eq!((-0f32).integer_decode(), (0u64, -150i16, -1i8)); assert_eq!((-0f32).integer_decode(), (0u64, -150i16, -1i8));
assert_eq!(INFINITY.integer_decode(), (8388608u64, 105i16, 1i8)); assert_eq!(INFINITY.integer_decode(), (8388608u64, 105i16, 1i8));

View File

@ -411,7 +411,7 @@ fn ln_10() -> f64 { 2.30258509299404568401799145468436421 }
fn recip(&self) -> f64 { 1.0 / *self } fn recip(&self) -> f64 { 1.0 / *self }
#[inline] #[inline]
fn pow(&self, n: &f64) -> f64 { pow(*self, *n) } fn powf(&self, n: &f64) -> f64 { pow(*self, *n) }
#[inline] #[inline]
fn sqrt(&self) -> f64 { sqrt(*self) } fn sqrt(&self) -> f64 { sqrt(*self) }
@ -1269,7 +1269,7 @@ fn test_frexp_nowin() {
fn test_integer_decode() { fn test_integer_decode() {
assert_eq!(3.14159265359f64.integer_decode(), (7074237752028906u64, -51i16, 1i8)); assert_eq!(3.14159265359f64.integer_decode(), (7074237752028906u64, -51i16, 1i8));
assert_eq!((-8573.5918555f64).integer_decode(), (4713381968463931u64, -39i16, -1i8)); assert_eq!((-8573.5918555f64).integer_decode(), (4713381968463931u64, -39i16, -1i8));
assert_eq!(2f64.pow(&100.0).integer_decode(), (4503599627370496u64, 48i16, 1i8)); assert_eq!(2f64.powf(&100.0).integer_decode(), (4503599627370496u64, 48i16, 1i8));
assert_eq!(0f64.integer_decode(), (0u64, -1075i16, 1i8)); assert_eq!(0f64.integer_decode(), (0u64, -1075i16, 1i8));
assert_eq!((-0f64).integer_decode(), (0u64, -1075i16, -1i8)); assert_eq!((-0f64).integer_decode(), (0u64, -1075i16, -1i8));
assert_eq!(INFINITY.integer_decode(), (4503599627370496u64, 972i16, 1i8)); assert_eq!(INFINITY.integer_decode(), (4503599627370496u64, 972i16, 1i8));

View File

@ -121,38 +121,6 @@ fn checked_mul(&self, v: &int) -> Option<int> {
} }
} }
/// Returns `base` raised to the power of `exponent`
pub fn pow(base: int, exponent: uint) -> int {
if exponent == 0u {
//Not mathemtically true if ~[base == 0]
return 1;
}
if base == 0 { return 0; }
let mut my_pow = exponent;
let mut acc = 1;
let mut multiplier = base;
while(my_pow > 0u) {
if my_pow % 2u == 1u {
acc *= multiplier;
}
my_pow /= 2u;
multiplier *= multiplier;
}
return acc;
}
#[test]
fn test_pow() {
assert!((pow(0, 0u) == 1));
assert!((pow(0, 1u) == 0));
assert!((pow(0, 2u) == 0));
assert!((pow(-1, 0u) == 1));
assert!((pow(1, 0u) == 1));
assert!((pow(-3, 2u) == 9));
assert!((pow(-3, 3u) == -27));
assert!((pow(4, 9u) == 262144));
}
#[test] #[test]
fn test_overflows() { fn test_overflows() {
assert!((::int::max_value > 0)); assert!((::int::max_value > 0));

View File

@ -185,9 +185,9 @@ pub trait Real: Signed
fn recip(&self) -> Self; fn recip(&self) -> Self;
// Algebraic functions // Algebraic functions
/// Raise a number to a power. /// Raise a number to a power.
fn pow(&self, n: &Self) -> Self; fn powf(&self, n: &Self) -> Self;
/// Take the square root of a number. /// Take the square root of a number.
fn sqrt(&self) -> Self; fn sqrt(&self) -> Self;
/// Take the reciprocal (inverse) square root of a number, `1/sqrt(x)`. /// Take the reciprocal (inverse) square root of a number, `1/sqrt(x)`.
@ -263,6 +263,50 @@ pub trait Real: Signed
fn to_radians(&self) -> Self; fn to_radians(&self) -> Self;
} }
/// Raises a value to the power of exp, using
/// exponentiation by squaring.
///
/// # Example
///
/// ```rust
/// use std::num;
///
/// let sixteen = num::pow(2, 4u);
/// assert_eq!(sixteen, 16);
/// ```
#[inline]
pub fn pow<T: Clone+One+Mul<T, T>>(num: T, exp: uint) -> T {
let one: uint = One::one();
let num_one: T = One::one();
if exp.is_zero() { return num_one; }
if exp == one { return num.clone(); }
let mut i: uint = exp;
let mut v: T;
let mut r: T = num_one;
// This if is to avoid cloning self.
if (i & one) == one {
r = r * num;
i = i - one;
}
i = i >> one;
v = num * num;
while !i.is_zero() {
if (i & one) == one {
r = r * v;
i = i - one;
}
i = i >> one;
v = v * v;
}
r
}
/// Raise a number to a power. /// Raise a number to a power.
/// ///
/// # Example /// # Example
@ -270,10 +314,10 @@ pub trait Real: Signed
/// ```rust /// ```rust
/// use std::num; /// use std::num;
/// ///
/// let sixteen: f64 = num::pow(2.0, 4.0); /// let sixteen: f64 = num::powf(2.0, 4.0);
/// assert_eq!(sixteen, 16.0); /// assert_eq!(sixteen, 16.0);
/// ``` /// ```
#[inline(always)] pub fn pow<T: Real>(value: T, n: T) -> T { value.pow(&n) } #[inline(always)] pub fn powf<T: Real>(value: T, n: T) -> T { value.powf(&n) }
/// Take the square root of a number. /// Take the square root of a number.
#[inline(always)] pub fn sqrt<T: Real>(value: T) -> T { value.sqrt() } #[inline(always)] pub fn sqrt<T: Real>(value: T) -> T { value.sqrt() }
/// Take the reciprocal (inverse) square root of a number, `1/sqrt(x)`. /// Take the reciprocal (inverse) square root of a number, `1/sqrt(x)`.
@ -1074,6 +1118,7 @@ pub fn test_num<T:Num + NumCast>(ten: T, two: T) {
mod tests { mod tests {
use prelude::*; use prelude::*;
use super::*; use super::*;
use num;
use i8; use i8;
use i16; use i16;
use i32; use i32;
@ -1634,4 +1679,19 @@ fn test_from_primitive() {
assert_eq!(from_f32(5f32), Some(Value { x: 5 })); assert_eq!(from_f32(5f32), Some(Value { x: 5 }));
assert_eq!(from_f64(5f64), Some(Value { x: 5 })); assert_eq!(from_f64(5f64), Some(Value { x: 5 }));
} }
#[test]
fn test_pow() {
fn assert_pow<T: Eq+Clone+One+Mul<T, T>>(num: T, exp: uint) -> () {
assert_eq!(num::pow(num.clone(), exp),
range(1u, exp).fold(num.clone(), |acc, _| acc * num));
}
assert_eq!(num::pow(3, 0), 1);
assert_eq!(num::pow(5, 1), 5);
assert_pow(-4, 2);
assert_pow(8, 3);
assert_pow(8, 5);
assert_pow(2u64, 50);
}
} }

View File

@ -144,7 +144,7 @@ impl IndependentSample<f64> for GammaSmallShape {
fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 { fn ind_sample<R: Rng>(&self, rng: &mut R) -> f64 {
let Open01(u) = rng.gen::<Open01<f64>>(); let Open01(u) = rng.gen::<Open01<f64>>();
self.large_shape.ind_sample(rng) * num::pow(u, self.inv_shape) self.large_shape.ind_sample(rng) * num::powf(u, self.inv_shape)
} }
} }
impl IndependentSample<f64> for GammaLargeShape { impl IndependentSample<f64> for GammaLargeShape {

View File

@ -31,10 +31,10 @@
use clone::Clone; use clone::Clone;
use kinds::Send; use kinds::Send;
use num::{Real, Round};
use option::{Option, Some, None}; use option::{Option, Some, None};
use sync::arc::UnsafeArc; use sync::arc::UnsafeArc;
use sync::atomics::{AtomicUint,Relaxed,Release,Acquire}; use sync::atomics::{AtomicUint,Relaxed,Release,Acquire};
use uint;
use vec; use vec;
struct Node<T> { struct Node<T> {
@ -64,7 +64,7 @@ fn with_capacity(capacity: uint) -> State<T> {
2u 2u
} else { } else {
// use next power of 2 as capacity // use next power of 2 as capacity
2f64.pow(&((capacity as f64).log2().ceil())) as uint uint::next_power_of_two(capacity)
} }
} else { } else {
capacity capacity

View File

@ -62,7 +62,7 @@ fn main() {
let long_lived_tree = bottom_up_tree(&long_lived_arena, 0, max_depth); let long_lived_tree = bottom_up_tree(&long_lived_arena, 0, max_depth);
let mut messages = range_step(min_depth, max_depth + 1, 2).map(|depth| { let mut messages = range_step(min_depth, max_depth + 1, 2).map(|depth| {
use std::int::pow; use std::num::pow;
let iterations = pow(2, (max_depth - depth + min_depth) as uint); let iterations = pow(2, (max_depth - depth + min_depth) as uint);
do Future::spawn { do Future::spawn {
let mut chk = 0; let mut chk = 0;