Talk:Legendre prime counting function

From Rosetta Code

several questions

Several questions: In principle this function has infinite recursion, so it should have a cut-off for n smaller than some number. 2? Then it is not clear how primes are counted I presume 2 is the first prime, 3 the second prime etc. Perhaps this could be made clear in the intro. The case of '1' is also a bit strange (i see 0 and 1 in the results, probably because some ambiguity in how the question is now posed/explained).

I checked and indeed pi(1) should be 0 (fixed that). Yes, the first prime would be 2.--Wherrera (talk) 17:27, 5 August 2021 (UTC)

There are a few ways to exit the recursion:

One is to use the base case that pi(n) = 0 for n < 2.

Another method would be to keep a table of small primes (say primes <= 127) and count from the array when n is small.

the third method is to count from the sieved primes that you need to generate anyway, using a binary search, as the sieved primes could be large enough for binary search to be preferred.

In the examples, the CoffeeScript uses the identity pi(n) = 0 when n < 2. The Erlang version performs a binary search on the sieved primes (method 3.) --Davgot (talk) 19:07, 5 August 2021 (UTC)

Thanks for clarifying, and thanks for the update to the description, this makes it a lot more clear. Updated my Mathematica solution to have pi(1) = 0.

V and C entries giving the wrong answer for a trillion numbers

The V entry gives a count of 35,460,428,370 primes as does the C entry (which I've translated from V).

However, according to this independent source, the correct answer is 37,607,912,018.

Another resource to help check results for large ranges is Kim Walisch's primecount. --GordonBGood (talk) 04:16, 11 November 2022 (UTC)

Curiously, the Go and Wren entries which I've also translated from V give the correct answer for a trillion numbers and the Go entry carries on giving the correct answers up to 10^15 numbers (Wren up to 10^13 - I didn't bother after that). However, both V and C produce a segmentation fault when I try 10^13.

I encountered segmentation faults for larger ranges on my V contribution and didn't look into it too much as I attributed it to V's "beta" status; interesting on the C/C++ stack segmentation problem as I first encountered this idea in C++ at this link and translated it to other languages starting with Nim, reorganizing, renaming variables, and adding comments, to make it more understandable (the competitive programmers from which this is taken seem to think that the algorithm is Meissel-Lehmer, likely due to the computational complexity, but as explained in the Nim text, it is not and its only a coincidence that it has about the same computational complexity).
I tested something like this C++ version and didn't get segmentation faults to 1e16 although there is an accuracy/precision problem at this range due to using floating point divisions as explained in the text; depending on the conversion to integer algorithm, this could cause segmentation faults so I converted the "half" and "divide" lambda's as follows:
  const auto divide = [](int64_t n, int64_t d) -> int64_t { return double(n) / d; };
  const auto half = [](int64_t n) -> int64_t { return (n - 1) >> 1; };
Since both V and C/C++ can use 64-bit integers as index variables, this and similar fixes wherever indexing is used fixes the segmentation problem and the accuracy problem is as you have analysed below. --GordonBGood (talk) 04:16, 11 November 2022 (UTC)

It eventually dawned on me that the reason for this is that the 'int' type in V and also C (at least in GCC) is 4 bytes, whereas in Go (and I think Nim) it's 8 bytes on a 64-bit system. Wren effectively has 53 bit integers.

Consequently, this line in the V example is producing an overflow:

mut ans := larges[0] + i64(((rilmt + 2 * (nbps - 1)) * (rilmt - 1) / 2))

If we change it to:

mut ans := larges[0] + i64(rilmt + 2 * (nbps - 1)) * i64(rilmt - 1) / 2

then the correct number of primes is returned for a trillion numbers.

Exactly, one of the problems with V is knowing when the compiler will promote smaller integers into larger ones; in this case it needed the prompt inside the brackets to cause the whole adjustment calculation to be done in 64-bit integers, and if one forces the conversion of any of the interior variables to 64-bits, then the external conversion notation isn't actually necessary as the internal result will already be in 64-bits. --GordonBGood (talk) 04:16, 11 November 2022 (UTC)

There's still a segmentation fault when I try to go above a trillion which I haven't been able to track down though, as the Go example is working, probably the easiest way to fix this is to change 'int' to 'i64' throughout. --PureFox (talk) 22:21, 10 November 2022 (UTC)

I think you'll find that the segmentation fault is the lack of precision in converting to a 32-bit index; in both V and C/C++ one can index by 64-bit integers rather than by 32-bit ones and converting the that will fix the segmentation problem. I've fixed the V version and will leave it up to you to do your C and/or C++ versions... --GordonBGood (talk) 04:16, 11 November 2022 (UTC)
BTW, thanks for "twigging" me on this issue leading to the resolution of both of the problems of the incorrect results and the segmentation faults for higher ranges. --GordonBGood (talk) 06:37, 11 November 2022 (UTC)
I've aligned my C example with the changes you've made to the V example and it's now working fine for a trillion numbers and beyond. As the Go and Wren examples were already working, I've just changed a couple of variable names rather than messing with the types.
It's all making perfect sense now :) --PureFox (talk) 09:51, 11 November 2022 (UTC)

A minor point but I've just noticed that the V entry (and hence the translations of it) don't deal with the trivial cases n = 3 to 8 inclusive. To solve it without messing with the algorithm I suggest inserting the following as the second line of the count_primes_to function:

if n < u64(9) { return i64((n + 1) / 2) }
--PureFox (talk) 20:50, 12 November 2022 (UTC)

If we do that, we can simplify the first line to:

if n < u64(2) { return i64(0) }

--PureFox (talk) 21:21, 12 November 2022 (UTC)