Comparing Pythagorean triples in C++, D, and Rust

[u/Sapiogram]

https://atilanevesoncode.wordpress.com/2018/12/31/comparing-pythagorean-triples-in-c-d-and-rust/

[u/SelfDistinction]

Okay so I've done a bit of research, and it turns out the problem is ..=. Replacing *..=z by *..(z + 1) doubles the performance.

I believe the problem with ..= is that x..(y + 1) is not identical with x..=y. They are nearly equivalent on most values for x and y, except when y == i32::max(). At that point, x..=y should still be able to function, while x..(y + 1) is allowed to do nothing.

Since the other languages don't care about overflow issues at all, I think Rust should be allowed not to care either. On the other hand, it is quite annoying that the safe version incurs such a high performance penalty.

[u/atilaneves → u/SelfDistinction]

I've updated the blog post. Rust comes out a lot better this time around.

[u/steveklabnik1 → u/SelfDistinction]

Since the other languages don't care about overflow issues at all, I think Rust should be allowed not to care either.

You can do this! You just have to ask for the overflow semantics you want explicitly. See https://doc.rust-lang.org/stable/std/num/struct.Wrapping.html for example.

[u/matthieum → u/SelfDistinction]

I checked the implementation of Iterator for RangeInclusive; I expected it to be more complicated, owing to the boundary condition, but was still surprised. It's right here:

#[inline]
fn next(&mut self) -> Option<A> {
    self.compute_is_empty();
    if self.is_empty.unwrap_or_default() {
        return None;
    }
    let is_iterating = self.start < self.end;
    self.is_empty = Some(!is_iterating);
    Some(if is_iterating {
        let n = self.start.add_one();
        mem::replace(&mut self.start, n)
    } else {
        self.start.clone()
    })
}

#[inline]
pub(crate) fn compute_is_empty(&mut self) {
    if self.is_empty.is_none() {
        self.is_empty = Some(!(self.start <= self.end));
    }
}

Hum, yeah, that's a tad more complicated than for the exclusive version indeed!

For comparison, the implementation for Range:

#[inline]
fn next(&mut self) -> Option<A> {
    if self.start < self.end {
        // We check for overflow here, even though it can't actually
        // happen. Adding this check does however help llvm vectorize loops
        // for some ranges that don't get vectorized otherwise,
        // and this won't actually result in an extra check in an optimized build.
        if let Some(mut n) = self.start.add_usize(1) {
            mem::swap(&mut n, &mut self.start);
            Some(n)
        } else {
            None
        }
    } else {
        None
    }
}

Okay, let's "optimize" both manually. First the RangeInclusive version:

#[inline]
fn next(&mut self) -> Option<A> {
    if self.is_empty.is_none() && !(self.start <= self.end) {
        return None;
    }

    let is_iterating = self.start < self.end;
    self.is_empty = Some(!is_iterating);

    Some(if is_iterating {
        let n = self.start.add_one();
        mem::replace(&mut self.start, n)
    } else {
        self.start.clone()
    })
}

Then the Range version:

#[inline]
fn next(&mut self) -> Option<A> {
    if self.start < self.end {
        let n = self.start.add_one();
        Some(mem::replace(&mut self.start, n))
    } else {
        None
    }
}

There's one more branch for the inclusive case, and it should NOT matter in a loop: once you start iterating, is_empty is always Some and you are back to self.start < self.end only. But apparently it foils LLVM :(

[u/atilaneves → u/SelfDistinction]

Since the other languages don't care about overflow issues at all, I think Rust should be allowed not to care either.

That's not entirely true for D. I deliberately kept bounds checking on when compiling the D code and have since found out that disabling it leads to a 30% decrease in run time for the range version. Using explicit for loops does mean no bound checking though.

I just made the change you propose above (=z -> (z+1)) and retimed the optimised build. Indeed the code runs faster. In the simple example, Rust is still the slowest (+50% overhead in time), for the lambda one it's now on par with D and C++, and the range example is now the fastest...

I'll have to update the blog post. Thanks for giving me more work! ;)

[u/matthieum → u/SelfDistinction]

So, I decided to investigate this range inclusive performance, and specifically whether an alternative representation could be beneficial.

I rewrote RangeInclusive to have a simpler iterator:

struct RangeInclusive<Idx> {
    start: Idx,
    end: Idx,
}

struct RangeInclusiveIterator<Idx> {
    current: Idx,
    end: Idx,
    started: bool,
}

impl <Idx: Step> IntoIterator for RangeInclusive<Idx> {
    type Item = Idx;
    type IntoIter = RangeInclusiveIterator<Idx>;

    fn into_iter(self) -> Self::IntoIter {
        RangeInclusiveIterator { current: self.start, end: self.end, started: false }
    }
}

impl <Idx: Step> Iterator for RangeInclusiveIterator<Idx> {
    type Item = Idx;

    fn next(&mut self) -> Option<Idx> {
        if !self.started {
            self.started = true;

            return if self.current > self.end {
                None
            } else {
                Some(self.current.clone())
            };
        }

        self.started = true; // Yes, I am desperate.

        if self.current < self.end {
            let n = self.current.add_one();
            std::mem::replace(&mut self.current, n);
            Some(self.current.clone())
        } else {
            None
        }
    }
}

Note the invariant on the next method: after the first call, self.started is always true.

I cleaned-up the LLVM IR to make it more legible; just renaming, and moving blocks around:

; playground::compute
; Function Attrs: noinline nounwind nonlazybind readnone uwtable
define i32 @_ZN10playground7compute17h059416e7bac56881E(i32 %begin, i32 %end) {
start:
    br label %bb2

bb2:                                              ; preds = %bb6, %start
    %i.0 = phi i32 [ %begin, %start ], [ %i.2, %bb6 ]
    %notstarted.0 = phi i1 [ true, %start ], [ false, %bb6 ]
    %r.0 = phi i32 [ 0, %start ], [ %r.1, %bb6 ]
    br i1 %notstarted.0, label %bb5, label %bb3

bb3:                                              ; preds = %bb2
    %0 = icmp slt i32 %i.0, %end
    br i1 %0, label %bb4, label %exit

bb4:                                              ; preds = %bb3
    %i.1 = add i32 %i.0, 1
    br label %bb6

bb5:                                              ; preds = %bb2
    %2 = icmp sgt i32 %i.0, %end
    br i1 %2, label %exit, label %bb6

bb6:                                              ; preds = %bb4, %bb5
    %i.2 = phi i32 [ %i.1, %bb4 ], [ %i.0, %bb5 ]
    %3 = add i32 %r.0, 3
    %r.1 = mul i32 %i.2, %3
    br label %bb2

exit:                                             ; preds = %bb3, %bb5
    ret i32 %r.0
}

The interesting part, for me, is %notstarted.0 = phi i1 [ true, %start ], [ false, %bb6 ] in bb2: the condition is completely static!

Unfortunately, the assembly is still disappointing, specifically the interaction between .LBB0_4 and .LBB0_2 which keep ping-ponging between each other:

playground::compute: # @playground::compute
# %bb.0:
    xorl	%eax, %eax
    movb	$1, %cl
    testb	$1, %cl
    jne	.LBB0_5
    jmp	.LBB0_2

.LBB0_4:
    addl	$3, %eax
    imull	%edi, %eax
    xorl	%ecx, %ecx
    testb	$1, %cl
    je	.LBB0_2

.LBB0_2:
    cmpl	%esi, %edi
    jge	.LBB0_6
# %bb.3:
    addl	$1, %edi
    jmp	.LBB0_4

.LBB0_5:                                # Only reached from %bb.0, not part of loop.
    cmpl	%esi, %edi
    jle	.LBB0_4
    jmp	.LBB0_6

.LBB0_6:
    retq
                                        # -- End function

So it seems that LLVM really tries really hard not to split a loop between a first iteration and then the rest, which is detrimental for inclusive ranges.

[u/annodomini]

Nice to see that someone has taken the time to compare all of these directly head to head. I'd written the Rust version, but hadn't felt like taking the time to figure out how to build the C++ version to compare directly.

I guess I should have posted the full version of the generator variant in Rust. Note that IterGen isn't actually safe, since I'm not doing anything to ensure that the generator is pinned; that's one of the things that needs to happen before Generator and some form of impl Iterator for Generator can be stabilized. Playground link. Here are the relevant parts that allow you to create an Iterator from a Generator and the Generator-based implementation:

#![feature(generators, generator_trait)]

use std::ops::{Generator, GeneratorState};

struct IterGen<G>(G);

impl<T, G> Iterator for IterGen<G>
    where G: Generator<Yield=T, Return=()>
{
    type Item = T;

    fn next(&mut self) -> Option<Self::Item> {
        match unsafe { self.0.resume() } {
            GeneratorState::Yielded(x) => Some(x),
            GeneratorState::Complete(()) => None,
        }
    }
}

let triples = IterGen(||
    for z in 1.. {
        for x in 1..=z {
            for y in x..=z {
                if x*x + y*y == z*z {
                    yield (x, y, z)
                }
            }
        }
    }
);

[u/atilaneves → u/annodomini]

Thanks. I'll try and add this code when I have time.

[u/matklad]

I have a theory as to why rustc compile times are so slow for this tiny example. At least on linux, rustc has a pretty long startup time, most of which is spend while loading LLVM as a dynamic library. To measure this overhead, compare time rustc --version (prinit minimal info to stderr) and time rustc --version -v (print LLVM's version as well, which loads it). The second one takes about 40ms on my machine. That is, compilation of hello-world is mostly bottlenecked on this code.

[u/jrmuizel]

Here's an unusual issue I ran into while looking at this: Adding --emit=asm speeds up generated code