Julia for loops slower than while?

Question

I'm playing around with the Julia language, and noticed that the small program I wrote was quite slow. Suspecting that it was somehow related to the for loops, I rewrote it to use while, and got about 15x faster.

I'm sure there's something I'm doing wrong with the ranges etc., but I can't figure out what.

function primes_for()
    num_primes = 0
    for a = 2:3000000
        prime = true
  
        sa = floor(sqrt(a))
        
        for c in 2:sa
             if a % c == 0
                 prime = false
                 break
             end
         end

        if prime
            num_primes += 1
        end
    end
    println("Number of primes is $num_primes")
end

function primes()
    num_primes = 0
    a = 2
    while a < 3000000
        prime = true

        c = 2
        while c * c <= a
            if a % c == 0
                prime = false
                break
            end
            c += 1
        end
  
        if prime
            num_primes += 1
        end

        a += 1
    end
    println("Number of primes is $num_primes")
end

@time primes_for()
@time primes()

Answer 1

As explained in the comments by @Vincent Yu and @Kelly Bundy, this is because sa = floor(sqrt(a)) creates a float. Then c becomes a float, and a % c is slow.

You can replace floor(sqrt(a)) with floor(Int, sqrt(a)), or preferably, I think, with isqrt(a), which returns

Integer square root: the largest integer m such that m*m <= n.

This avoids the (unlikely) event that floor(Int, sqrt(a)) may round down too far, which could happen if sqrt(x^2) = x - ε due to floating point errors.

Edit: Here's a benchmark to demonstrate (note the use of isqrt):

function primes_for2()
    num_primes = 0
    for a = 2:3000000
        prime = true
  
        # sa = floor(sqrt(a))
        sa = isqrt(a)
        
        for c in 2:sa
             if a % c == 0
                 prime = false
                 break
             end
         end

        if prime
            num_primes += 1
        end
    end
    println("Number of primes is $num_primes")
end

1.7.0> @time primes_for()
Number of primes is 216816
  6.705099 seconds (15 allocations: 480 bytes)

1.7.0> @time primes_for2()
Number of primes is 216816
  0.691304 seconds (15 allocations: 480 bytes)

1.7.0> @time primes()
Number of primes is 216816
  0.671784 seconds (15 allocations: 480 bytes)

I can note that each call to isqrt on my computer takes approximately 8ns, and that 3000000 times 8ns is 0.024 seconds. A call to regular sqrt is approximately 1ns.

Answer 2

It's not the for/while that makes the speed difference, it's the sqrt. It doesn't help that sqrt returns float, which promotes all the rest of the code around the sqrt output from integers.

Note that @time is not measuring the while and for loops, but also the code outside those loops.

If you are benchmarking code, the rest of your code needs to be the same, and removing the sqrt is one of the prime optimizations in this algorithm. It's also possible to remove the c * c in the test, but this is trickier.