Sieve of Eratosthenes with Accelerate

Sieve.swift

// Calculate prime numbers in a given range
// using Sieve of Eratosthenes

import Accelerate

var primes: [Int] = []

let range = 0...999
var numbers = range.map(Float.init)

for n in range.dropLast() {
    guard numbers[n] > 1 else { continue }
    defer { primes.append(n) }
  
    // Clear vector elements with given stride
    // (i.e. set value to 0)
    vDSP_vclr(&numbers, 
              numericCast(n), 
              numericCast(numbers.count / n))
}

print(primes)

But it doesn't work if range is from 0...100, or 0...1000. It thinks 99 and 999 are prime.

@willie Right you are! I noticed that as well, but haven't had a chance to dig into the apparent off-by-one error. Maybe it's as simple as taking the ceiling of the quotient?

I'm not quite sure what you mean when you talk about "taking the ceiling of the quotient" in this context

My Swift game is very weak, but my suspicion is that since vDSP_vclr is zero-based, and there are n elements (and not n-1 elements like in C, where you have an array of 0 to n-1), that in some cases the tail end is not getting cleared. Adding +1 to numbers.count ( https://gist.github.com/mattt/d725a54f9b5aa6a3eaa425f4e1f18f78#file-sieve-swift-L19 ) before the division seems to solve the bug locally, but I'm concerned that it might be a buffer overrun...

This is weird.

Put:

	let cnt = (numbers.count) / n
	print(numbers.count, n, cnt)
	numbers[0] = 1.0
	print(numbers)
	vDSP_vclr(&numbers,
			  numericCast(n),
			  numericCast(cnt))
	print(numbers)

results in this for n=3

11 3 3
[1.0, 1.0, 0.0, 3.0, 0.0, 5.0, 0.0, 7.0, 0.0, 9.0, 10.0]
[0.0, 1.0, 0.0, 0.0, 0.0, 5.0, 0.0, 7.0, 0.0, 9.0, 10.0]

It looks like vDSP_vclr clears numbers[0] regardless of stride and elements?

My understanding of the stride parameter is that a value of 2 , means "clear every 2nd (even) element", 3 means "clear every 3rd element)", etc. So yes, it clears numbers[0] each time, but that's fine because numbers < 2 aren't considered.

To clarify my point from before, I wonder if always rounding up numbers.count / n (like with ceil) would solve the problem. Adding 1 as you described would have the same effect (at least for any divisions that are rounded down.

Here is what I figured out:

1. numbers.count is the last number in the range + 1, so the number of elements calculation should be: (numbers.count - 1) / n
2. there is no need to dropLast
3. to account for the zeroth element behavior, it's easy enough to start vDSP_vclr at &numbers[n].

I'll post a modified version in a few.

Try this: https://gist.github.com/willie/c589832bc009350578acdbfda4dbcd2c