Float

Bite 297 showed signum reading the sign bit. copysign is the other half: it writes one — take this magnitude, give it that number’s sign, one call.

The branch

You have a fixed magnitude and a direction living in another variable — a step size, a force, a correction — and out comes the same conditional every time:

1
2
3

let delta: f64 = -3.2;

let step = if delta < 0.0 { -0.25 } else { 0.25 };

One call

1
2
3
4
5

let step = 0.25_f64.copysign(delta); // -0.25

assert_eq!(3.5_f64.copysign(-1.0), -3.5);
assert_eq!((-3.5_f64).copysign(1.0), 3.5);
assert_eq!(7.0_f64.copysign(2.0), 7.0);

The receiver contributes the magnitude, the argument contributes the sign. Yesterday’s “step toward the target” loop works for floats too — magnitude from you, direction from the gap:

1
2
3
4
5
6
7

let mut pos = 1.0_f64;
let target = 2.0;

while (target - pos).abs() >= 0.25 {
    pos += 0.25_f64.copysign(target - pos);
}
assert_eq!(pos, 2.0); // 0.25 divides 1.0 exactly

It really is the sign bit

Like f64::signum, copysign works on the raw sign bit — which is exactly what makes it more honest than the if:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

// the branch calls -0.0 "positive"...
let x = -0.0_f64;
let s = if x < 0.0 { -1.0 } else { 1.0 };
assert_eq!(s, 1.0);

// ...copysign actually looks at the bit
assert!(1.0_f64.copysign(-0.0).is_sign_negative());

// and a NaN's magnitude stays NaN, sign still applied
assert!(f64::NAN.copysign(-1.0).is_sign_negative());

Anywhere you’re bolting a direction onto a magnitude — drag opposing velocity, rounding away from zero, mirroring a coordinate — copysign is the branch-free, -0.0-correct way to say it.

a * b + c rounds twice and may compile to two separate operations. a.mul_add(b, c) computes the whole thing at full precision, rounds once, and on a CPU with FMA folds into a single instruction.

Two roundings vs. one

Floating-point math rounds after every operation. Writing a * b + c first rounds the product a * b to the nearest f64, then rounds again after adding c. That intermediate rounding throws away bits before they ever reach the sum:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

let a = 1.0_f64;
let b = 2.0_f64;
let c = 3.0_f64;

// Two roundings: (a*b) rounded, then (+ c) rounded
let separate = a * b + c;

// One rounding: a*b + c evaluated at full precision, rounded once
let fused = a.mul_add(b, c);

assert_eq!(separate, 5.0);
assert_eq!(fused, 5.0);

For tidy values they agree. The difference shows up when the product needs more bits than an f64 can hold and the addition would have recovered them — exactly the catastrophic-cancellation cases that wreck numeric code.

Where it pays off: Horner’s method

Polynomial evaluation is built out of multiply-then-add, so it’s the textbook home for mul_add. To evaluate 2x² + 3x + 4 with Horner’s method you nest the operations, and each step is one fused step:

1
2
3
4
5
6

fn poly(x: f64) -> f64 {
    // ((2*x) + 3)*x + 4
    2.0_f64.mul_add(x, 3.0).mul_add(x, 4.0)
}

assert_eq!(poly(2.0), 18.0); // 2*4 + 3*2 + 4

Each mul_add keeps full precision through the chain instead of rounding at every * and +. The same pattern carries dot products, which are a sum of products:

1
2
3
4
5
6
7

fn dot(a: &[f64], b: &[f64]) -> f64 {
    a.iter()
        .zip(b)
        .fold(0.0, |acc, (&x, &y)| x.mul_add(y, acc))
}

assert_eq!(dot(&[1.0, 2.0, 3.0], &[4.0, 5.0, 6.0]), 32.0);

The honest caveat

mul_add is faster only when the target CPU has a hardware fused-multiply-add instruction (modern x86-64 with FMA, AArch64, most others you’ll deploy to). Where it doesn’t, the standard library has to emulate the exact single-rounding semantics in software, and that emulation is slower than a plain a * b + c. So this is a hot-path tool: reach for it in tight numeric loops on hardware you control (or behind target-feature/target-cpu flags), and lean on it freely when you need the extra accuracy regardless of speed.

The bottom line

mul_add gives you a single rounding step and, on capable hardware, a single instruction for a * b + c. Use it in polynomial evaluation, dot products, and any multiply-accumulate loop where precision or throughput matters — and remember it can regress on FMA-less targets, so measure when speed is the goal.

Ever wondered what the next representable floating point number after 1.0 is? Since Rust 1.86, f64::next_up and f64::next_down let you step through the number line one float at a time.

The problem

Floating point numbers aren’t evenly spaced — the gap between representable values grows as the magnitude increases. Before next_up / next_down, figuring out the next neighbor required bit-level manipulation of the IEEE 754 representation:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

fn main() {
    // The hard way (before 1.86): manually decode bits
    let x: f64 = 1.0;
    let bits = x.to_bits();
    let next_bits = bits + 1;
    let next = f64::from_bits(next_bits);

    assert!(next > x);
    assert_eq!(next, 1.0000000000000002);
}

Error-prone, unreadable, and doesn’t handle edge cases like negative numbers, zero, or special values.

The clean way

next_up returns the smallest f64 greater than self. next_down returns the largest f64 less than self:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15

fn main() {
    let x: f64 = 1.0;

    let up = x.next_up();
    let down = x.next_down();

    assert!(up > x);
    assert!(down < x);
    assert_eq!(up, 1.0000000000000002);
    assert_eq!(down, 0.9999999999999998);

    // There's no float between x and its neighbors
    assert_eq!(up.next_down(), x);
    assert_eq!(down.next_up(), x);
}

They handle all the edge cases you’d rather not think about — negative numbers, subnormals, and infinity:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

fn main() {
    // Works across zero
    assert_eq!(0.0_f64.next_up(), 5e-324);   // smallest positive f64
    assert_eq!(0.0_f64.next_down(), -5e-324); // largest negative f64

    // Infinity is the boundary
    assert_eq!(f64::MAX.next_up(), f64::INFINITY);
    assert_eq!(f64::INFINITY.next_up(), f64::INFINITY);

    // NaN stays NaN
    assert!(f64::NAN.next_up().is_nan());
}

Practical use: precision-aware comparisons

The gap between adjacent floats is called an ULP (unit in the last place). next_up lets you build tolerance-aware comparisons without guessing at epsilon values:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21

fn almost_equal(a: f64, b: f64, max_ulps: u32) -> bool {
    if a == b { return true; }

    let mut current = a;
    for _ in 0..max_ulps {
        current = if a < b { current.next_up() } else { current.next_down() };
        if current == b { return true; }
    }
    false
}

fn main() {
    let a = 0.1 + 0.2;
    let b = 0.3;

    // They're not equal...
    assert_ne!(a, b);

    // ...but they're within 1 ULP of each other
    assert!(almost_equal(a, b, 1));
}

Also available on f32 with the same API. These methods are const fn, so you can use them in const contexts too.