I tried to implement the sieve of Erastosthenes as fast as possible.
I didn't thought it would make a big difference, but having a bit-wise
vector<bool> instead of a vector<char> makes the code twice as fast
for my AMD 7950X PC and an upper bound of one billion. The additional
overhead for the bit-wise access is less than the additional overhead
of fetching more memory.
I also tried Eulers sieve but as I expected filtering the numbers by
modulo operations with the recent prime numbers takes much more time
than simply resetting a bit; the resulting code becomes mulltiple
times slower.

#include <iostream>
#include <vector>
#include <charconv>
#include <string_view>
#include <algorithm>
#include <fstream>
#include <cctype>
#include <cstring>

#if defined(_MSC_VER)
#pragma warning(disable: 4996)
#endif

using namespace std;

int main( int argc, char **argv )
{ constexpr bool USE_BOOL = true;
constexpr size_t BUF_SIZE = 0x100000;
try
{
size_t max;
if( argc < 2 )
return EXIT_FAILURE;
char const *sizeStr = argv[1];
bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
sizeStr += 2 * hex;
if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
return EXIT_FAILURE;
using bool_t = conditional_t<USE_BOOL, bool, char>;
vector<bool_t> sieve( max, true );
for( size_t p = 2; p < max; )
{
for( size_t m = 2 * p; m < max; m += p )
sieve[m] = false;
while( ++p < max && !sieve[p] );
}
if( argc >= 3 && !*argv[2] )
return EXIT_SUCCESS;
ofstream ofs;
ofs.exceptions( ofstream::failbit | ofstream::badbit );
ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
union ndi_t { }; // non-default-initialzied
vector<ndi_t> buf;
buf.resize( BUF_SIZE + 32 );
char
*bufBegin = (char *)buf.data(),
*bufEnd = bufBegin;
for( size_t i = 2; i != max; ++i )
if( sieve[i] )
{
to_chars_result tcr = to_chars( bufEnd, bufEnd + 32, i );
(bufEnd = tcr.ptr + 1)[-1] = '\n';
if( bufEnd - bufBegin >= BUF_SIZE )
ofs << string_view( bufBegin, bufEnd ),
bufEnd = bufBegin;
}
ofs << string_view( bufBegin, bufEnd );
}
catch( bad_alloc const & )
{
cout << "out of memory" << endl;
}
catch( ios_base::failure const & )
{
cout << "I/O error" << endl;
}
}

I also optimized the file output because simply doing a ofs << i << endl
is multiple times slower than buffering a megabyte of data until pushing
it to the file. The numbers printed to the buffer are directly printed
inside the buffer to make the code more efficient. So if the upper bound
is one billion the time taken to print the numbers to the file is within
the measurement error.

On Wednesday, August 16, 2023 at 9:38:10 PM UTC+8, Bonita Montero wrote:
> I tried to implement the sieve of Erastosthenes as fast as possible.
> I didn't thought it would make a big difference, but having a bit-wise
> vector<bool> instead of a vector<char> makes the code twice as fast
> for my AMD 7950X PC and an upper bound of one billion. The additional
> overhead for the bit-wise access is less than the additional overhead
> of fetching more memory.
> I also tried Eulers sieve but as I expected filtering the numbers by
> modulo operations with the recent prime numbers takes much more time
> than simply resetting a bit; the resulting code becomes mulltiple
> times slower.
>
> #include <iostream>
> #include <vector>
> #include <charconv>
> #include <string_view>
> #include <algorithm>
> #include <fstream>
> #include <cctype>
> #include <cstring>
>
> #if defined(_MSC_VER)
> #pragma warning(disable: 4996)
> #endif
>
> using namespace std;
>
> int main( int argc, char **argv )
> {
> constexpr bool USE_BOOL = true;
> constexpr size_t BUF_SIZE = 0x100000;
> try
> {
> size_t max;
> if( argc < 2 )
> return EXIT_FAILURE;
> char const *sizeStr = argv[1];
> bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
> sizeStr += 2 * hex;
> if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
> sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
> return EXIT_FAILURE;
> using bool_t = conditional_t<USE_BOOL, bool, char>;
> vector<bool_t> sieve( max, true );
> for( size_t p = 2; p < max; )
> {
> for( size_t m = 2 * p; m < max; m += p )
> sieve[m] = false;
> while( ++p < max && !sieve[p] );
> }
> if( argc >= 3 && !*argv[2] )
> return EXIT_SUCCESS;
> ofstream ofs;
> ofs.exceptions( ofstream::failbit | ofstream::badbit );
> ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
> union ndi_t { }; // non-default-initialzied
> vector<ndi_t> buf;
> buf.resize( BUF_SIZE + 32 );
> char
> *bufBegin = (char *)buf.data(),
> *bufEnd = bufBegin;
> for( size_t i = 2; i != max; ++i )
> if( sieve[i] )
> {
> to_chars_result tcr = to_chars( bufEnd, bufEnd + 32, i );
> (bufEnd = tcr.ptr + 1)[-1] = '\n';
> if( bufEnd - bufBegin >= BUF_SIZE )
> ofs << string_view( bufBegin, bufEnd ),
> bufEnd = bufBegin;
> }
> ofs << string_view( bufBegin, bufEnd );
> }
> catch( bad_alloc const & )
> {
> cout << "out of memory" << endl;
> }
> catch( ios_base::failure const & )
> {
> cout << "I/O error" << endl;
> }
> }
>
> I also optimized the file output because simply doing a ofs << i << endl
> is multiple times slower than buffering a megabyte of data until pushing
> it to the file. The numbers printed to the buffer are directly printed
> inside the buffer to make the code more efficient. So if the upper bound
> is one billion the time taken to print the numbers to the file is within
> the measurement error.

Confirmed, but you code don't compile.

--------------------------------- t.cpp
#include <Wy.stdio.h>
#include <Wy.time.h>
#include <CSCall/VLInt.h>
#include <vector>

using namespace Wy;

// [Syn] Print prime numbers less or equal to n
//
// Algo: sieve
//
void print_prime1(size_t n)
{ Errno r;
VLInt<unsigned int> barr; // bit=0 =prime

for(size_t i=2; i<=n; ++i) {
if(barr.bit(i)) {
continue;
}
// Cond: i is a prime number
for(size_t j=2*i; j<=n; j+=i) {
if((r=barr.set_bit(j,true))!=Ok) {
WY_THROW(r);
}
}
}
};

void print_prime2(size_t n)
{ Errno r;
Array<char> barr(n+1,0,0); // bit=0 =prime

for(size_t i=2; i<=n; ++i) {
if(barr[i]) {
continue;
}
// Cond: i is a prime number
for(size_t j=2*i; j<=n; j+=i) {
barr[j]=1;
}
}
};

void print_prime3(size_t n)
{ Errno r;
std::vector<char> barr(n+1,0); // bit=0 =prime

for(size_t i=2; i<=n; ++i) {
if(barr.at(i)) {
continue;
}
// Cond: i is a prime number
for(size_t j=2*i; j<=n; j+=i) {
barr.at(j)=1;
}
}
};

void print_prime4(size_t n)
{ Errno r;
std::vector<bool> barr(n+1,false); // bit=0 =prime

for(size_t i=2; i<=n; ++i) {
if(barr.at(i)) {
continue;
}
// Cond: i is a prime number
for(size_t j=2*i; j<=n; j+=i) {
barr.at(j)=true;
}
}
};

int main()
try {
constexpr size_t Cnt=100000000L;
Timespec tm0;

tm0=now(CLOCK_THREAD_CPUTIME_ID);
print_prime1(Cnt);
cout << "VLInt<unsigned>: " << (now(CLOCK_THREAD_CPUTIME_ID)-tm0) << WY_ENDL;

tm0=now(CLOCK_THREAD_CPUTIME_ID);
print_prime2(Cnt);
cout << "Array<char>: " << (now(CLOCK_THREAD_CPUTIME_ID)-tm0) << WY_ENDL;

tm0=now(CLOCK_THREAD_CPUTIME_ID);
print_prime3(Cnt);
cout << "vector<char>: " << (now(CLOCK_THREAD_CPUTIME_ID)-tm0) << WY_ENDL;

tm0=now(CLOCK_THREAD_CPUTIME_ID);
print_prime4(Cnt);
cout << "vector<bool>: " << (now(CLOCK_THREAD_CPUTIME_ID)-tm0) << WY_ENDL;

cout << "Ok" WY_ENDL;
return 0;
} catch(const Errno& e) {
cerr << wrd(e) << WY_ENDL;
return -1;
} catch(...) {
cerr << "unknown stack unwind" << WY_ENDL;
throw;
};
----------------------------------

[]$ g++ t.cpp -lwy
[]$ ./a.out
VLInt<unsigned>: 6.162083886
Array<char>: 3.555815058
vector<char>: 5.307075961
vector<bool>: 30.257717265

I am surprised std::vector's at(size_t) is unexpected slow.

Note: Wy::Array's operator[] checks index range, should be equ. to std::vector's at(size_t)

Am 16.08.2023 um 19:57 schrieb wij:
> On Wednesday, August 16, 2023 at 9:38:10 PM UTC+8, Bonita Montero wrote:
>> I tried to implement the sieve of Erastosthenes as fast as possible.
>> I didn't thought it would make a big difference, but having a bit-wise
>> vector<bool> instead of a vector<char> makes the code twice as fast
>> for my AMD 7950X PC and an upper bound of one billion. The additional
>> overhead for the bit-wise access is less than the additional overhead
>> of fetching more memory.
>> I also tried Eulers sieve but as I expected filtering the numbers by
>> modulo operations with the recent prime numbers takes much more time
>> than simply resetting a bit; the resulting code becomes mulltiple
>> times slower.
>>
>> #include <iostream>
>> #include <vector>
>> #include <charconv>
>> #include <string_view>
>> #include <algorithm>
>> #include <fstream>
>> #include <cctype>
>> #include <cstring>
>>
>> #if defined(_MSC_VER)
>> #pragma warning(disable: 4996)
>> #endif
>>
>> using namespace std;
>>
>> int main( int argc, char **argv )
>> {
>> constexpr bool USE_BOOL = true;
>> constexpr size_t BUF_SIZE = 0x100000;
>> try
>> {
>> size_t max;
>> if( argc < 2 )
>> return EXIT_FAILURE;
>> char const *sizeStr = argv[1];
>> bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
>> sizeStr += 2 * hex;
>> if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
>> sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
>> return EXIT_FAILURE;
>> using bool_t = conditional_t<USE_BOOL, bool, char>;
>> vector<bool_t> sieve( max, true );
>> for( size_t p = 2; p < max; )
>> {
>> for( size_t m = 2 * p; m < max; m += p )
>> sieve[m] = false;
>> while( ++p < max && !sieve[p] );
>> }
>> if( argc >= 3 && !*argv[2] )
>> return EXIT_SUCCESS;
>> ofstream ofs;
>> ofs.exceptions( ofstream::failbit | ofstream::badbit );
>> ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
>> union ndi_t { }; // non-default-initialzied
>> vector<ndi_t> buf;
>> buf.resize( BUF_SIZE + 32 );
>> char
>> *bufBegin = (char *)buf.data(),
>> *bufEnd = bufBegin;
>> for( size_t i = 2; i != max; ++i )
>> if( sieve[i] )
>> {
>> to_chars_result tcr = to_chars( bufEnd, bufEnd + 32, i );
>> (bufEnd = tcr.ptr + 1)[-1] = '\n';
>> if( bufEnd - bufBegin >= BUF_SIZE )
>> ofs << string_view( bufBegin, bufEnd ),
>> bufEnd = bufBegin;
>> }
>> ofs << string_view( bufBegin, bufEnd );
>> }
>> catch( bad_alloc const & )
>> {
>> cout << "out of memory" << endl;
>> }
>> catch( ios_base::failure const & )
>> {
>> cout << "I/O error" << endl;
>> }
>> }
>>
>> I also optimized the file output because simply doing a ofs << i << endl
>> is multiple times slower than buffering a megabyte of data until pushing
>> it to the file. The numbers printed to the buffer are directly printed
>> inside the buffer to make the code more efficient. So if the upper bound
>> is one billion the time taken to print the numbers to the file is within
>> the measurement error.
>
> Confirmed, but you code don't compile.

-std=c++20

>
> --------------------------------- t.cpp
> #include <Wy.stdio.h>
> #include <Wy.time.h>
> #include <CSCall/VLInt.h>
> #include <vector>
>
> using namespace Wy;
>
> // [Syn] Print prime numbers less or equal to n
> //
> // Algo: sieve
> //
> void print_prime1(size_t n)
> {
> Errno r;
> VLInt<unsigned int> barr; // bit=0 =prime
>
> for(size_t i=2; i<=n; ++i) {
> if(barr.bit(i)) {
> continue;
> }
> // Cond: i is a prime number
> for(size_t j=2*i; j<=n; j+=i) {
> if((r=barr.set_bit(j,true))!=Ok) {
> WY_THROW(r);
> }
> }
> }
> };
>
> void print_prime2(size_t n)
> {
> Errno r;
> Array<char> barr(n+1,0,0); // bit=0 =prime
>
> for(size_t i=2; i<=n; ++i) {
> if(barr[i]) {
> continue;
> }
> // Cond: i is a prime number
> for(size_t j=2*i; j<=n; j+=i) {
> barr[j]=1;
> }
> }
> };
>
> void print_prime3(size_t n)
> {
> Errno r;
> std::vector<char> barr(n+1,0); // bit=0 =prime
>
> for(size_t i=2; i<=n; ++i) {
> if(barr.at(i)) {
> continue;
> }
> // Cond: i is a prime number
> for(size_t j=2*i; j<=n; j+=i) {
> barr.at(j)=1;
> }
> }
> };
>
> void print_prime4(size_t n)
> {
> Errno r;
> std::vector<bool> barr(n+1,false); // bit=0 =prime
>
> for(size_t i=2; i<=n; ++i) {
> if(barr.at(i)) {
> continue;
> }
> // Cond: i is a prime number
> for(size_t j=2*i; j<=n; j+=i) {
> barr.at(j)=true;
> }
> }
> };

Show me the times to print primes below 10000000.

16.08.2023 20:57 wij kirjutas:
>
> []$ g++ t.cpp -lwy
>
> I am surprised std::vector's at(size_t) is unexpected slow.

You are missing -O2 or -O3. There is no point to talk about the speed of
unoptimized code.

On Thursday, August 17, 2023 at 2:45:12 AM UTC+8, Paavo Helde wrote:
> 16.08.2023 20:57 wij kirjutas:
> >
> > []$ g++ t.cpp -lwy
> >
> > I am surprised std::vector's at(size_t) is unexpected slow.
> You are missing -O2 or -O3. There is no point to talk about the speed of
> unoptimized code.

[]$ g++ t.cpp -lwy -O2
[]$ ./a.out
VLInt<unsigned>: 1.259642293
Array<char>: 1.525368272
vector<char>: 1.484200452
vector<bool>: 1.533252238

----------------
-O2 makes sense. But the problem now is why VLInt is faster? (VLInt allocates memory on demand)

On Thursday, August 17, 2023 at 2:22:02 AM UTC+8, Bonita Montero wrote:
> Am 16.08.2023 um 19:57 schrieb wij:
> > On Wednesday, August 16, 2023 at 9:38:10 PM UTC+8, Bonita Montero wrote:
> >> I tried to implement the sieve of Erastosthenes as fast as possible.
> >> I didn't thought it would make a big difference, but having a bit-wise
> >> vector<bool> instead of a vector<char> makes the code twice as fast
> >> for my AMD 7950X PC and an upper bound of one billion. The additional
> >> overhead for the bit-wise access is less than the additional overhead
> >> of fetching more memory.
> >> I also tried Eulers sieve but as I expected filtering the numbers by
> >> modulo operations with the recent prime numbers takes much more time
> >> than simply resetting a bit; the resulting code becomes mulltiple
> >> times slower.
> >>
> >> #include <iostream>
> >> #include <vector>
> >> #include <charconv>
> >> #include <string_view>
> >> #include <algorithm>
> >> #include <fstream>
> >> #include <cctype>
> >> #include <cstring>
> >>
> >> #if defined(_MSC_VER)
> >> #pragma warning(disable: 4996)
> >> #endif
> >>
> >> using namespace std;
> >>
> >> int main( int argc, char **argv )
> >> {
> >> constexpr bool USE_BOOL = true;
> >> constexpr size_t BUF_SIZE = 0x100000;
> >> try
> >> {
> >> size_t max;
> >> if( argc < 2 )
> >> return EXIT_FAILURE;
> >> char const *sizeStr = argv[1];
> >> bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
> >> sizeStr += 2 * hex;
> >> if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
> >> sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
> >> return EXIT_FAILURE;
> >> using bool_t = conditional_t<USE_BOOL, bool, char>;
> >> vector<bool_t> sieve( max, true );
> >> for( size_t p = 2; p < max; )
> >> {
> >> for( size_t m = 2 * p; m < max; m += p )
> >> sieve[m] = false;
> >> while( ++p < max && !sieve[p] );
> >> }
> >> if( argc >= 3 && !*argv[2] )
> >> return EXIT_SUCCESS;
> >> ofstream ofs;
> >> ofs.exceptions( ofstream::failbit | ofstream::badbit );
> >> ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
> >> union ndi_t { }; // non-default-initialzied
> >> vector<ndi_t> buf;
> >> buf.resize( BUF_SIZE + 32 );
> >> char
> >> *bufBegin = (char *)buf.data(),
> >> *bufEnd = bufBegin;
> >> for( size_t i = 2; i != max; ++i )
> >> if( sieve[i] )
> >> {
> >> to_chars_result tcr = to_chars( bufEnd, bufEnd + 32, i );
> >> (bufEnd = tcr.ptr + 1)[-1] = '\n';
> >> if( bufEnd - bufBegin >= BUF_SIZE )
> >> ofs << string_view( bufBegin, bufEnd ),
> >> bufEnd = bufBegin;
> >> }
> >> ofs << string_view( bufBegin, bufEnd );
> >> }
> >> catch( bad_alloc const & )
> >> {
> >> cout << "out of memory" << endl;
> >> }
> >> catch( ios_base::failure const & )
> >> {
> >> cout << "I/O error" << endl;
> >> }
> >> }
> >>
> >> I also optimized the file output because simply doing a ofs << i << endl
> >> is multiple times slower than buffering a megabyte of data until pushing
> >> it to the file. The numbers printed to the buffer are directly printed
> >> inside the buffer to make the code more efficient. So if the upper bound
> >> is one billion the time taken to print the numbers to the file is within
> >> the measurement error.
> >
> > Confirmed, but you code don't compile.
> -std=c++20
> >
> > --------------------------------- t.cpp
> > #include <Wy.stdio.h>
> > #include <Wy.time.h>
> > #include <CSCall/VLInt.h>
> > #include <vector>
> >
> > using namespace Wy;
> >
> > // [Syn] Print prime numbers less or equal to n
> > //
> > // Algo: sieve
> > //
> > void print_prime1(size_t n)
> > {
> > Errno r;
> > VLInt<unsigned int> barr; // bit=0 =prime
> >
> > for(size_t i=2; i<=n; ++i) {
> > if(barr.bit(i)) {
> > continue;
> > }
> > // Cond: i is a prime number
> > for(size_t j=2*i; j<=n; j+=i) {
> > if((r=barr.set_bit(j,true))!=Ok) {
> > WY_THROW(r);
> > }
> > }
> > }
> > };
> >
> > void print_prime2(size_t n)
> > {
> > Errno r;
> > Array<char> barr(n+1,0,0); // bit=0 =prime
> >
> > for(size_t i=2; i<=n; ++i) {
> > if(barr[i]) {
> > continue;
> > }
> > // Cond: i is a prime number
> > for(size_t j=2*i; j<=n; j+=i) {
> > barr[j]=1;
> > }
> > }
> > };
> >
> > void print_prime3(size_t n)
> > {
> > Errno r;
> > std::vector<char> barr(n+1,0); // bit=0 =prime
> >
> > for(size_t i=2; i<=n; ++i) {
> > if(barr.at(i)) {
> > continue;
> > }
> > // Cond: i is a prime number
> > for(size_t j=2*i; j<=n; j+=i) {
> > barr.at(j)=1;
> > }
> > }
> > };
> >
> > void print_prime4(size_t n)
> > {
> > Errno r;
> > std::vector<bool> barr(n+1,false); // bit=0 =prime
> >
> > for(size_t i=2; i<=n; ++i) {
> > if(barr.at(i)) {
> > continue;
> > }
> > // Cond: i is a prime number
> > for(size_t j=2*i; j<=n; j+=i) {
> > barr.at(j)=true;
> > }
> > }
> > };
> Show me the times to print primes below 10000000.

[]$ g++ t12.cpp -lwy -std=c++20
[]$ ./a.out
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 211 2
....(cut)
889 9999901 9999907 9999929 9999931 9999937 9999943 9999971 9999973 9999991
pi(10000000)= 664579
print_prime: 2.740338666
Ok

--------------------- t12.cpp
#include <Wy.stdio.h>
#include <Wy.time.h>
#include <CSCall/VLInt.h>

using namespace Wy;

// [Syn] Print prime numbers less or equal to n
//
// Algo: sieve
//
void print_prime(size_t n)
{ Errno r;
VLInt<unsigned int> barr; // bit=0 =prime

for(size_t i=2; i<=n; ++i) {
if(barr.bit(i)) {
continue;
}
// Cond: i is a prime number
for(size_t j=2*i; j<=n; j+=i) {
if((r=barr.set_bit(j,true))!=Ok) {
WY_THROW(r);
}
}
}

// print prime number (from 2)
size_t pcnt=0;
for(size_t i=2; i<=n; ++i) {
if(barr.bit(i)==false) {
cout << wrd(i) << ' ';
++pcnt;
}
}
cout << WY_ENDL;
cout << "pi(" << n << ")= " << pcnt << WY_ENDL;
};

int main()
try {
Timespec tm0;

tm0=now(CLOCK_THREAD_CPUTIME_ID);
print_prime(10000000);
cout << "print_prime: " << (now(CLOCK_THREAD_CPUTIME_ID)-tm0) << WY_ENDL;
cout << "Ok" WY_ENDL;
return 0;
} catch(const Errno& e) {
cerr << wrd(e) << WY_ENDL;
return -1;
} catch(...) {
cerr << "unknown stack unwind" << WY_ENDL;
throw;
};

On 2023-08-16 3:37 PM, Bonita Montero wrote:
> I tried to implement the sieve of Erastosthenes as fast as possible.
> I didn't thought it would make a big difference, but having a bit-wise
> vector<bool> instead of a vector<char> makes the code twice as fast
> for my AMD 7950X PC and an upper bound of one billion. The additional
> overhead for the bit-wise access is less than the additional overhead
> of fetching more memory.

An apparent mystery, because I get only a marginal difference, down at
10% or less.

On my computer the code below is slightly faster than yours for a
10'000'000 sieve, respectively (typical results, three runs, only the
vector allocation + sieve code measured, not file i/o)

0.0477802 seconds.
0.0466974 seconds.
0.0472202 seconds.

.... for my code, versus

0.0502888 seconds.
0.0482821 seconds.
0.0505414 seconds.

.... for yours with the `std::string_view` instantiations fixed,
producing the same file contents.

#include <fmt/core.h>
#include "cppm-Timer.hpp"

#include <vector>
using std::ofstream,
std::vector;

#include <stdlib.h> // EXIT_FAILURE
#include <stdio.h> // fopen, fclose

auto main() -> int
{
constexpr int n = 10'000'000;

cppm::Timer timer;
auto is_composite = vector<bool>( n + 1 );
for( int i = 2; i < n; ) {
while( is_composite[i] ) { ++i; }
const int prime = i;
for( int j = 2*prime; j < n; j += prime ) { is_composite[j]
= true; }
++i;
}
timer.stop();
fmt::print( "{} seconds.\n", timer.seconds() );

auto f = fopen( "primes.alf.txt", "w" );
if( not f ) { return EXIT_FAILURE; }
for( int i = 2; i < n; ++i ) {
if( not is_composite[i] ) {
fmt::print( f, "{:d}\n", i );
}
}
fclose( f );
}

The timer code:

#pragma once

#include <chrono>
#include <optional>
#include <type_traits>

namespace cppm {
namespace chrono = std::chrono;
using std::optional, // <optional>
std::conditional_t; // <type_traits>

using Duration = chrono::duration<double>;

class Timer
{
using Hrc = chrono::high_resolution_clock; using Sc =
chrono::steady_clock;

public:
using Clock = conditional_t<Hrc::is_steady, Hrc, Sc>;
using Time_point = chrono::time_point<Clock>;

private:
Time_point m_start_time;
optional<Time_point> m_opt_end_time;

public:
Timer(): m_start_time( Clock::now() ), m_opt_end_time() {}

void stop() { if( not m_opt_end_time ) { m_opt_end_time =
Clock::now(); } }

auto duration() const -> Duration { return
m_opt_end_time.value() - m_start_time; }
auto seconds() const -> double { return
duration().count(); }
};
} // namespace cppm

I have some other more advanced prime sieve code that I made for an
example, but I think it's much slower, it's just a way that can go on
and on forever, or more precisely till the end of the integer type's
range or it runs out of memory for a map of the deduces primes.

- Alf

On 2023-08-16 11:17 PM, Alf P. Steinbach wrote:
> ...

Unfortunately the `using std::ofstream`, left in the code due to
imperfect editing, compiled -- with Visual C++. Gah! Not used.

- Alf (wishing for old Modula-2, just with better syntax etc.)

On 2023-08-16 11:17 PM, Alf P. Steinbach wrote:
> On 2023-08-16 3:37 PM, Bonita Montero wrote:
>> I tried to implement the sieve of Erastosthenes as fast as possible.
>> I didn't thought it would make a big difference, but having a bit-wise
>> vector<bool> instead of a vector<char> makes the code twice as fast
>> for my AMD 7950X PC and an upper bound of one billion. The additional
>> overhead for the bit-wise access is less than the additional overhead
>> of fetching more memory.
>
>
> An apparent mystery, because I get only a marginal difference, down at
> 10% or less.
>
> On my computer the code below is slightly faster than yours for a
> 10'000'000 sieve, respectively (typical results, three runs, only the
> vector allocation + sieve code measured, not file i/o)
>
>     0.0477802 seconds.
>     0.0466974 seconds.
>     0.0472202 seconds.
>
> ... for my code, versus
>
>     0.0502888 seconds.
>     0.0482821 seconds.
>     0.0505414 seconds.
>
> ... for yours with the `std::string_view` instantiations fixed,
> producing the same file contents.
>
>
>     #include <fmt/core.h>
>     #include "cppm-Timer.hpp"
>
>     #include <vector>
>     using   std::ofstream,
>             std::vector;
>
>     #include <stdlib.h>         // EXIT_FAILURE
>     #include <stdio.h>          // fopen, fclose
>
>     auto main() -> int
>     {
>         constexpr int n = 10'000'000;
>
>         cppm::Timer timer;
>         auto is_composite = vector<bool>( n + 1 );
>         for( int i = 2; i < n; ) {
>             while( is_composite[i] ) { ++i; }
>             const int prime = i;
>             for( int j = 2*prime; j < n; j += prime ) { is_composite[j]
> = true; }
>             ++i;
>         }
>         timer.stop();
>         fmt::print( "{} seconds.\n", timer.seconds() );
>
>         auto f = fopen( "primes.alf.txt", "w" );
>         if( not f ) { return EXIT_FAILURE; }
>         for( int i = 2; i < n; ++i ) {
>             if( not is_composite[i] ) {
>                 fmt::print( f, "{:d}\n", i );
>             }
>         }
>         fclose( f );
>     }
>
>
> The timer code:
>
>
>     #pragma once
>
>     #include <chrono>
>     #include <optional>
>     #include <type_traits>
>
>     namespace cppm {
>         namespace chrono = std::chrono;
>         using   std::optional,          // <optional>
>                 std::conditional_t;     // <type_traits>
>
>         using Duration      = chrono::duration<double>;
>
>         class Timer
>         {
>             using Hrc = chrono::high_resolution_clock;  using Sc =
> chrono::steady_clock;
>
>         public:
>             using Clock         = conditional_t<Hrc::is_steady, Hrc, Sc>;
>             using Time_point    = chrono::time_point<Clock>;
>
>         private:
>             Time_point              m_start_time;
>             optional<Time_point>    m_opt_end_time;
>
>         public:
>             Timer(): m_start_time( Clock::now() ), m_opt_end_time() {}
>
>             void stop() { if( not m_opt_end_time ) { m_opt_end_time =
> Clock::now(); } }
>
>             auto duration() const   -> Duration { return
> m_opt_end_time.value() - m_start_time; }
>             auto seconds() const    -> double   { return
> duration().count(); }
>         };
>     }  // namespace cppm
>
>
> I have some other more advanced prime sieve code that I made for an
> example, but I think it's much slower, it's just a way that can go on
> and on forever, or more precisely till the end of the integer type's
> range or it runs out of memory for a map of the deduces primes.

Well, I've become stupid in my old days. Reduced running time by say 25%
or so by very simple optimization. I did not think of that until /after/
I'd posted.

New running time, three tries:

0.0365686 seconds.
0.0362102 seconds.
0.0359283 seconds.

New better shiny sparkling code:

#include <fmt/core.h>
#include "cppm-Timer.hpp"

#include <vector>
using std::vector;

#include <stdlib.h> // EXIT_FAILURE
#include <stdio.h> // fopen, fclose

auto main() -> int
{
constexpr int n = 10'000'000;
constexpr int sqrt_n = 3'162;

cppm::Timer timer;
auto is_composite = vector<bool>( n + 1 );
for( int i = 2; i < n; ) {
while( is_composite[i] ) { ++i; }
const int prime = i;
if( prime < sqrt_n ) {
for( int j = prime*prime; j < n; j += prime ) {
is_composite[j] = true;
}
}
++i;
}
timer.stop();
fmt::print( "{} seconds.\n", timer.seconds() );

auto f = fopen( "primes.alf.txt", "w" );
if( not f ) { return EXIT_FAILURE; }
for( int i = 2; i < n; ++i ) {
if( not is_composite[i] ) {
fmt::print( f, "{:d}\n", i );
}
}
fclose( f );
}

Clue: for efficiency, avoid totally needless work.

- Alf

On 2023-08-16 11:39 PM, Alf P. Steinbach wrote:
> On 2023-08-16 11:17 PM, Alf P. Steinbach wrote:
>> On 2023-08-16 3:37 PM, Bonita Montero wrote:
>>> I tried to implement the sieve of Erastosthenes as fast as possible.
>>> I didn't thought it would make a big difference, but having a bit-wise
>>> vector<bool> instead of a vector<char> makes the code twice as fast
>>> for my AMD 7950X PC and an upper bound of one billion. The additional
>>> overhead for the bit-wise access is less than the additional overhead
>>> of fetching more memory.
>>
>>
>> An apparent mystery, because I get only a marginal difference, down at
>> 10% or less.
>>
>> On my computer the code below is slightly faster than yours for a
>> 10'000'000 sieve, respectively (typical results, three runs, only the
>> vector allocation + sieve code measured, not file i/o)
>>
>>      0.0477802 seconds.
>>      0.0466974 seconds.
>>      0.0472202 seconds.
>>
>> ... for my code, versus
>>
>>      0.0502888 seconds.
>>      0.0482821 seconds.
>>      0.0505414 seconds.
>>
>> ... for yours with the `std::string_view` instantiations fixed,
>> producing the same file contents.
>>
>>
>>      #include <fmt/core.h>
>>      #include "cppm-Timer.hpp"
>>
>>      #include <vector>
>>      using   std::ofstream,
>>              std::vector;
>>
>>      #include <stdlib.h>         // EXIT_FAILURE
>>      #include <stdio.h>          // fopen, fclose
>>
>>      auto main() -> int
>>      {
>>          constexpr int n = 10'000'000;
>>
>>          cppm::Timer timer;
>>          auto is_composite = vector<bool>( n + 1 );
>>          for( int i = 2; i < n; ) {
>>              while( is_composite[i] ) { ++i; }
>>              const int prime = i;
>>              for( int j = 2*prime; j < n; j += prime ) {
>> is_composite[j] = true; }
>>              ++i;
>>          }
>>          timer.stop();
>>          fmt::print( "{} seconds.\n", timer.seconds() );
>>
>>          auto f = fopen( "primes.alf.txt", "w" );
>>          if( not f ) { return EXIT_FAILURE; }
>>          for( int i = 2; i < n; ++i ) {
>>              if( not is_composite[i] ) {
>>                  fmt::print( f, "{:d}\n", i );
>>              }
>>          }
>>          fclose( f );
>>      }
>>
>>
>> The timer code:
>>
>>
>>      #pragma once
>>
>>      #include <chrono>
>>      #include <optional>
>>      #include <type_traits>
>>
>>      namespace cppm {
>>          namespace chrono = std::chrono;
>>          using   std::optional,          // <optional>
>>                  std::conditional_t;     // <type_traits>
>>
>>          using Duration      = chrono::duration<double>;
>>
>>          class Timer
>>          {
>>              using Hrc = chrono::high_resolution_clock;  using Sc =
>> chrono::steady_clock;
>>
>>          public:
>>              using Clock         = conditional_t<Hrc::is_steady, Hrc,
>> Sc>;
>>              using Time_point    = chrono::time_point<Clock>;
>>
>>          private:
>>              Time_point              m_start_time;
>>              optional<Time_point>    m_opt_end_time;
>>
>>          public:
>>              Timer(): m_start_time( Clock::now() ), m_opt_end_time() {}
>>
>>              void stop() { if( not m_opt_end_time ) { m_opt_end_time =
>> Clock::now(); } }
>>
>>              auto duration() const   -> Duration { return
>> m_opt_end_time.value() - m_start_time; }
>>              auto seconds() const    -> double   { return
>> duration().count(); }
>>          };
>>      }  // namespace cppm
>>
>>
>> I have some other more advanced prime sieve code that I made for an
>> example, but I think it's much slower, it's just a way that can go on
>> and on forever, or more precisely till the end of the integer type's
>> range or it runs out of memory for a map of the deduces primes.
>
> Well, I've become stupid in my old days. Reduced running time by say 25%
> or so by very simple optimization. I did not think of that until /after/
> I'd posted.
>
> New running time, three tries:
>
>     0.0365686 seconds.
>     0.0362102 seconds.
>     0.0359283 seconds.
>
> New better shiny sparkling code:
>
>     #include <fmt/core.h>
>     #include "cppm-Timer.hpp"
>
>     #include <vector>
>     using   std::vector;
>
>     #include <stdlib.h>         // EXIT_FAILURE
>     #include <stdio.h>          // fopen, fclose
>
>     auto main() -> int
>     {
>         constexpr int n = 10'000'000;
>         constexpr int sqrt_n = 3'162;
>
>         cppm::Timer timer;
>         auto is_composite = vector<bool>( n + 1 );
>         for( int i = 2; i < n; ) {
>             while( is_composite[i] ) { ++i; }
>             const int prime = i;
>             if( prime < sqrt_n ) {
>                 for( int j = prime*prime; j < n; j += prime ) {
>                     is_composite[j] = true;
>                 }
>             }
>             ++i;
>         }
>         timer.stop();
>         fmt::print( "{} seconds.\n", timer.seconds() );
>
>         auto f = fopen( "primes.alf.txt", "w" );
>         if( not f ) { return EXIT_FAILURE; }
>         for( int i = 2; i < n; ++i ) {
>             if( not is_composite[i] ) {
>                 fmt::print( f, "{:d}\n", i );
>             }
>         }
>         fclose( f );
>     }
>
> Clue: for efficiency, avoid totally needless work.

Oh, I'm REALLY old.

Even less slow, now clocking in at

0.0260816 seconds.
0.0267381 seconds.
0.0269437 seconds.

Code:

#include <fmt/core.h>
#include "cppm-Timer.hpp"

#include <vector>
using std::vector;

#include <stdlib.h> // EXIT_FAILURE
#include <stdio.h> // fopen, fclose

auto main() -> int
{
constexpr int n = 10'000'000;
constexpr int sqrt_n = 3'162;

cppm::Timer timer;
auto is_composite = vector<bool>( n + 1 );
for( int i = 4; i < n; i += 2 ) { is_composite[i] = true; }
for( int i = 3; i < n; ) {
while( is_composite[i] ) { ++i; }
const int prime = i;
if( prime <= sqrt_n ) {
for( int j = prime*prime; j < n; j += 2*prime ) {
is_composite[j] = true;
}
}
++i;
}
timer.stop();
fmt::print( "{} seconds.\n", timer.seconds() );

auto f = fopen( "primes.alf.txt", "w" );
if( not f ) { return EXIT_FAILURE; }
for( int i = 2; i < n; ++i ) {
if( not is_composite[i] ) {
fmt::print( f, "{:d}\n", i );
}
}
fclose( f );
}

But every little optimization adds complexity and verbosity.

- Alf

On 8/16/2023 2:59 PM, Alf P. Steinbach wrote:
[...]

Fwiw, keep in mind what your system is currently doing when you run a
benchmark. Personally, I shutdown as many unnecessary processes as I can
before I run any benchmark. Bare bones.

On 8/16/2023 2:17 PM, Alf P. Steinbach wrote:
> On 2023-08-16 3:37 PM, Bonita Montero wrote:
>> I tried to implement the sieve of Erastosthenes as fast as possible.
>> I didn't thought it would make a big difference, but having a bit-wise
>> vector<bool> instead of a vector<char> makes the code twice as fast
>> for my AMD 7950X PC and an upper bound of one billion. The additional
>> overhead for the bit-wise access is less than the additional overhead
>> of fetching more memory
[redacted]

Here's mine. I get rid of the vector<bool> by storing the prime factors
in their own vector, and only checking for divisibility by primes.

It can be sped up by checking for divisibility only by values less than
or equal to sqrt(candidate value).

I only check up to 1 million, but you can obviously change max_prime.
Also change to long long for values greater than 2 billion.

--- cut here --

#include <iostream>
#include <ostream>
#include <vector>
#include <cmath>
#include <cstdlib>

const long max_prime = 1000000;
const long max_factor =
static_cast<long>(std::sqrt(static_cast<double>(max_prime))) + 1;

int main(int, char *[])
{ std::vector<long> factors;

std::cout << "2" << std::endl;
for (long n = 3; n < max_prime ; n += 2 )
{
bool composite;
composite = false;
for (std::vector<long>::iterator it = factors.begin();
!composite && it != factors.end();
++it)
{
if ((n % *it) == 0)
composite = true;
}
if (!composite)
{
std::cout << n << std::endl;

if (n < max_factor)
factors.push_back(n);
}
}

return EXIT_SUCCESS;
}

red floyd <no.spam.here@its.invalid> writes:

> On 8/16/2023 2:17 PM, Alf P. Steinbach wrote:
>
>> On 2023-08-16 3:37 PM, Bonita Montero wrote:
>>
>>> I tried to implement the sieve of Erastosthenes as fast as
>>> possible. I didn't thought it would make a big difference, but
>>> having a bit-wise vector<bool> instead of a vector<char> makes the
>>> code twice as fast for my AMD 7950X PC and an upper bound of one
>>> billion. The additional overhead for the bit-wise access is less
>>> than the additional overhead of fetching more memory
>
> [redacted]
>
> Here's mine. I get rid of the vector<bool> by storing the prime
> factors in their own vector, and only checking for divisibility by
> primes. [...]

Here are some suggestions for everyone.

Except for 2, 3, and 5, numbers can be considered mod 30, and
there are (conveniently) only 8 candidates: 1, 7, 11, 13, 17,
19, 23, and 29. For numbers 7 and up, every prime must be
one of those candidates (mod 30).

By storing the results for blocks of 30 in a byte, we can
store results for all 32 bit primes in 143165577 bytes.

There are some extra tricks taking advantage of what happens
with the bit numbers. For example 7 * 37 has two 7's, which
is 49, which is 19 mod 30, so the bit for that product will
be in the 19 position. Being clever will avoid any shifting
or modding.

I suggest computing all 32 bit primes for the benchmark
tests. Shouldn't take too long.

Am 17.08.2023 um 03:06 schrieb red floyd:

> Here's mine.  I get rid of the vector<bool> by storing the prime factors
> in their own vector, and only checking for divisibility by primes.

If I modify my code to use the same, called Euler's sieve, my code
becomes 19 times slower on my 7950X with an upper bound of 1E6.

Am 17.08.2023 um 08:42 schrieb Bonita Montero:
> Am 17.08.2023 um 03:06 schrieb red floyd:
>
>> Here's mine.  I get rid of the vector<bool> by storing the prime factors
>> in their own vector, and only checking for divisibility by primes.
>
> If I modify my code to use the same, called Euler's sieve, my code
> becomes 19 times slower on my 7950X with an upper bound of 1E6.
>

This is the code:

#include <iostream>
#include <vector>
#include <charconv>
#include <string_view>
#include <algorithm>
#include <fstream>
#include <cctype>
#include <cstring>

using namespace std;

int main( int argc, char **argv )
{ constexpr bool USE_BOOL = true;
constexpr size_t BUF_SIZE = 0x100000;
try
{
size_t max;
if( argc < 2 )
return EXIT_FAILURE;
char const *sizeStr = argv[1];
bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
sizeStr += 2 * hex;
if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
return EXIT_FAILURE;
vector<size_t> primes;
for( size_t p = 2; p < max; ++p )
{
bool has = false;
for( size_t soFar : primes )
if( !(p % soFar) )
{
has = true;
break;
}
if( !has )
primes.emplace_back( p );
}
if( argc >= 3 && !*argv[2] )
return EXIT_SUCCESS;
ofstream ofs;
ofs.exceptions( ofstream::failbit | ofstream::badbit );
ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
union ndi_t { }; // non-default-initialzied
vector<ndi_t> buf;
buf.resize( BUF_SIZE + 32 );
char
*bufBegin = (char *)buf.data(),
*bufEnd = bufBegin,
*const bufFlush = bufBegin + BUF_SIZE;
for( size_t p : primes )
{
to_chars_result tcr = to_chars( bufEnd, bufEnd + 32, p );
bufEnd = tcr.ptr;
*bufEnd++ = '\n';
if( bufEnd >= bufFlush )
ofs << string_view( bufBegin, bufEnd ),
bufEnd = bufBegin;
}
ofs << string_view( bufBegin, bufEnd );
}
catch( bad_alloc const & )
{
cout << "out of memory" << endl;
}
catch( ios_base::failure const & )
{
cout << "I/O error" << endl;
}
}

I tried that with an upper bound of 1E8 and after 2,5 hours of computing
I stopped the process. With my primary code the process ends after 0.4s.

On 2023-08-17 11:38 AM, Bonita Montero wrote:
> Am 17.08.2023 um 08:42 schrieb Bonita Montero:
>> Am 17.08.2023 um 03:06 schrieb red floyd:
>>
>>> Here's mine.  I get rid of the vector<bool> by storing the prime factors
>>> in their own vector, and only checking for divisibility by primes.
>>
>> If I modify my code to use the same, called Euler's sieve, my code
>> becomes 19 times slower on my 7950X with an upper bound of 1E6.
>>
>
> This is the code:

[snip]

>         vector<size_t> primes;
>         for( size_t p = 2; p < max; ++p )
>         {
>             bool has = false;
>             for( size_t soFar : primes )
>                 if( !(p % soFar) )
>                 {
>                     has = true;
>                     break;
>                 }
>             if( !has )
>                 primes.emplace_back( p );
>         }

[snip]

> I tried that with an upper bound of 1E8 and after 2,5 hours of computing
> I stopped the process. With my primary code the process ends after 0.4s.

One less unreasonable way to check against already computed primes is to
generate and update the smallest multiple of each already known prime,
that is larger than the current number to be checked.

It reduces the complexity for each number check to roughly logarithmic
in the number of primes so far, instead of linear. And in total roughly
O(n log n) instead of O(n^2). Quadratic is just impractical, as you found.

Code for the collection-of-smallest-multiples approach:

#include <fmt/core.h>
#include "cppm-Timer.hpp"

#include <queue>
#include <vector>
using std::priority_queue,
std::vector;

#include <stdio.h> // fopen, fclose

template< class Item, class Compare_func >
using Priority_q = priority_queue<Item, std::vector<Item>,
Compare_func>;

struct Multiple
{
int base; int n;
auto value() const -> int{ return n*base; }
};

struct Descending_order
{
auto operator()( const Multiple& a, const Multiple& b ) const
-> bool
{ return a.value() > b.value(); } // Note: direction of
comparison.
};

auto main() -> int
{
constexpr int n = 10'000'000;

cppm::Timer timer;
auto primes = vector<int>();
primes.reserve( n/2 );
const auto produce = [&]( const int prime ) { primes.push_back(
prime ); };

produce( 2 );
auto larger_composites = Priority_q<Multiple, Descending_order>();
larger_composites.push( Multiple{2, 2} );
for( int current = 3; current < n; ++current ) {
bool is_composite = false;
for( ;; ) {
const Multiple next_composite = larger_composites.top();
if( current < next_composite.value() ) {
break;
}
larger_composites.pop();
is_composite = true;
larger_composites.push( Multiple{next_composite.base, 1
+ next_composite.n} );
}
if( not is_composite ) {
produce( current );
larger_composites.push( Multiple{current, 2} );
}
}
timer.stop();
fmt::print( "{} seconds.\n", timer.seconds() );

auto f = fopen( "primes.alf-gen.txt", "w" );
if( not f ) { return EXIT_FAILURE; }
for( const int prime: primes ) { fmt::print( f, "{:d}\n", prime
); }
fclose( f );
}

It runs roughly 100 times slower than the unslowest of the earlier
programs I posted in this thread, using 3.5736102 seconds to produce the
primes up to 10'000'000, but I believe still much faster than the
quadratic time code you showed (I would guess maybe 400 seconds = 6'40
minutes for the 1 billion case?), and has the possible advantage that
one doesn't need to know in advance a cap on the prime size.

It can just run and run, gobbling up memory as it goes.

- Alf

Am 17.08.2023 um 12:50 schrieb Alf P. Steinbach:

> One less unreasonable way to check against already computed primes is to
> generate and update the smallest multiple of each already known prime,
> that is larger than the current number to be checked.

Beginning with the smaller primes aldready seen does make sense
since you can filter the values earlier. F.e. there are more
values which are a multiple of two than a multiple of 11.

I did it the same way with some simpler code:

#include <iostream>
#include <queue>
#include <vector>
#include <charconv>
#include <fstream>

using namespace std;

int main( int argc, char **argv )
{ constexpr size_t BUF_SIZE = 0x100000;
size_t max;
if( argc < 2 )
return EXIT_FAILURE;
char const *sizeStr = argv[1];
bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
sizeStr += 2 * hex;
if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
return EXIT_FAILURE;
struct next_multiple { size_t p, pNext; };
auto nm_less = []( next_multiple const &left, next_multiple const
&right ) { return left.pNext > right.pNext; };
priority_queue<next_multiple, vector<next_multiple>, decltype(nm_less)>
nextMultiples( nm_less );
nextMultiples.emplace( 2, 4 );
auto *hitOrNot = &nextMultiples.top();
for( size_t p = 3; p <= max; ++p )
{
if( p < hitOrNot->pNext ) [[unlikely]]
{
nextMultiples.emplace( p, p + p );
hitOrNot = &nextMultiples.top();
continue;
}
do
{
auto nextIncr = *hitOrNot;
nextMultiples.pop();
nextIncr.pNext += nextIncr.p;
nextMultiples.emplace( nextIncr );
hitOrNot = &nextMultiples.top();
} while( p == hitOrNot->pNext );
}
vector<size_t> primes;
primes.reserve( nextMultiples.size() );
while( !nextMultiples.empty() )
{
auto const next = nextMultiples.top();
primes.emplace_back( next.p );
nextMultiples.pop();
}
sort( primes.begin(), primes.end() );
if( argc >= 3 && !*argv[2] )
return EXIT_FAILURE;
ofstream ofs;
ofs.exceptions( ofstream::failbit | ofstream::badbit );
ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
union ndi_t { }; // non-default-initialzied
vector<ndi_t> buf( BUF_SIZE + 32 );
char
*bufBegin = (char*)buf.data(),
*bufEnd = bufBegin,
*bufFlush = bufBegin + BUF_SIZE;
for( size_t p : primes )
{
to_chars_result tcr = to_chars( bufEnd, bufEnd + 32, p );
bufEnd = tcr.ptr;
*bufEnd++ = '\n';
if( bufEnd >= bufFlush )
ofs << string_view( bufBegin, bufEnd ),
bufEnd = bufBegin;
}
ofs << string_view( bufBegin, bufEnd );
}

If I print all primes <= 1E8 that's still 82 times slower than my
initial solution.

On Thursday, August 17, 2023 at 5:17:54 AM UTC+8, Alf P. Steinbach wrote:
> On 2023-08-16 3:37 PM, Bonita Montero wrote:
> > I tried to implement the sieve of Erastosthenes as fast as possible.
> > I didn't thought it would make a big difference, but having a bit-wise
> > vector<bool> instead of a vector<char> makes the code twice as fast
> > for my AMD 7950X PC and an upper bound of one billion. The additional
> > overhead for the bit-wise access is less than the additional overhead
> > of fetching more memory.
> An apparent mystery, because I get only a marginal difference, down at
> 10% or less.
>
> On my computer the code below is slightly faster than yours for a
> 10'000'000 sieve, respectively (typical results, three runs, only the
> vector allocation + sieve code measured, not file i/o)
>
> 0.0477802 seconds.
> 0.0466974 seconds.
> 0.0472202 seconds.
>
> ... for my code, versus
>
> 0.0502888 seconds.
> 0.0482821 seconds.
> 0.0505414 seconds.
>
> ... for yours with the `std::string_view` instantiations fixed,
> producing the same file contents.
>
>
> #include <fmt/core.h>
> #include "cppm-Timer.hpp"
>
> #include <vector>
> using std::ofstream,
> std::vector;
>
> #include <stdlib.h> // EXIT_FAILURE
> #include <stdio.h> // fopen, fclose
>
> auto main() -> int
> {
> constexpr int n = 10'000'000;
>
> cppm::Timer timer;
> auto is_composite = vector<bool>( n + 1 );
> for( int i = 2; i < n; ) {
> while( is_composite[i] ) { ++i; }
> const int prime = i;
> for( int j = 2*prime; j < n; j += prime ) { is_composite[j]
> = true; }
> ++i;
> }
> timer.stop();
> fmt::print( "{} seconds.\n", timer.seconds() );
>
> auto f = fopen( "primes.alf.txt", "w" );
> if( not f ) { return EXIT_FAILURE; }
> for( int i = 2; i < n; ++i ) {
> if( not is_composite[i] ) {
> fmt::print( f, "{:d}\n", i );
> }
> }
> fclose( f );
> }
>
>
> The timer code:
>
>
> #pragma once
>
> #include <chrono>
> #include <optional>
> #include <type_traits>
>
> namespace cppm {
> namespace chrono = std::chrono;
> using std::optional, // <optional>
> std::conditional_t; // <type_traits>
>
> using Duration = chrono::duration<double>;
>
> class Timer
> {
> using Hrc = chrono::high_resolution_clock; using Sc =
> chrono::steady_clock;
>
> public:
> using Clock = conditional_t<Hrc::is_steady, Hrc, Sc>;
> using Time_point = chrono::time_point<Clock>;
>
> private:
> Time_point m_start_time;
> optional<Time_point> m_opt_end_time;
>
> public:
> Timer(): m_start_time( Clock::now() ), m_opt_end_time() {}
>
> void stop() { if( not m_opt_end_time ) { m_opt_end_time =
> Clock::now(); } }
>
> auto duration() const -> Duration { return
> m_opt_end_time.value() - m_start_time; }
> auto seconds() const -> double { return
> duration().count(); }
> };
> } // namespace cppm
>
>
> I have some other more advanced prime sieve code that I made for an
> example, but I think it's much slower, it's just a way that can go on
> and on forever, or more precisely till the end of the integer type's
> range or it runs out of memory for a map of the deduces primes.
>
>
> - Alf

Code Comment:
1. The use of 'not' is impressive, bringing dead to life.
2. The creation of class Timer may be questionable. Or it is just specifically
to simply and emphasize the main codes.
3. I don't know what "fmt::print( f, "{:d}\n", i );" really mean, but in general
I think cstring format in its current usage (e.g. "foo %s= %d\n") does not
provide real benefit for the executing time spent interpreting it.

On 2023-08-17 10:43 PM, wij wrote:
> On Thursday, August 17, 2023 at 5:17:54 AM UTC+8, Alf P. Steinbach wrote:
[snip]
>
> Code Comment:
> 1. The use of 'not' is impressive, bringing dead to life.

I guess that means you didn't like it, but I do.

It may be an acquired taste, like coffee and beer. And I came to C++
from Pascals and Modula-2 and such. But try it: you may come to like it.

> 2. The creation of class Timer may be questionable. Or it is just specifically
> to simply and emphasize the main codes.

Factoring out that functionality in a header makes it easy to reuse.

And it was reused.

But even with no hope of reuse you should not be afraid to define useful
abstractions, as opposed to just using others' abstractions.

> 3. I don't know what "fmt::print( f, "{:d}\n", i );" really mean, but in general
> I think cstring format in its current usage (e.g. "foo %s= %d\n") does not
> provide real benefit for the executing time spent interpreting it.

With the {fmt} library's `fmt::print`, which will become C++23's
`std::print`, the format string is interpreted at compile time.

That's almost like a lint implemented via C++ template system, plus that
it's just about the fastest string formatting around, plus handles
Unicode (guaranteed by the standard) if the output device can do it.

The design of the compile time format string handling does however not
adhere to the principle of least surprise, and that caused the
standard's wording to be a bit sabotage-like for some time. That was
eventually fixed. But means that somewhere I have this C++20 code:

#if __cpp_lib_format >= (2022*100 + 7)
#define myprefix_FMTLIB_API_PROVIDES_FORMAT_STRING_TYPE
constexpr bool api_provides_format_string_type = true;

- Alf

On Friday, August 18, 2023 at 8:00:49 AM UTC+8, Alf P. Steinbach wrote:
> On 2023-08-17 10:43 PM, wij wrote:

> > 3. I don't know what "fmt::print( f, "{:d}\n", i );" really mean, but in general
> > I think cstring format in its current usage (e.g. "foo %s= %d\n") does not
> > provide real benefit for the executing time spent interpreting it.
> With the {fmt} library's `fmt::print`, which will become C++23's
> `std::print`, the format string is interpreted at compile time.
>
> That's almost like a lint implemented via C++ template system, plus that
> it's just about the fastest string formatting around, plus handles
> Unicode (guaranteed by the standard) if the output device can do it.
>
> The design of the compile time format string handling does however not
> adhere to the principle of least surprise, and that caused the
> standard's wording to be a bit sabotage-like for some time. That was
> eventually fixed. But means that somewhere I have this C++20 code:
>
> #if __cpp_lib_format >= (2022*100 + 7)
> #define myprefix_FMTLIB_API_PROVIDES_FORMAT_STRING_TYPE
> constexpr bool api_provides_format_string_type = true;
>
>
> - Alf

Not only time, but many kind of resources are used by std::print.

My idea of the 'format' issue was adding a format member to class String should
be enough, e.g. "String& String::fmt(..)", where the .. means unclear.

So, the usecase (example) printf("IC0123 qty:%5d, price:%.2f\n",qty,prc) can be
translated like the following:

result << String("IC0123 qty:") << wrd(qty).fmt("5") << ", price:"
<< wrd(prc).fmt(".2") << "\n";

('resutl' can be String or any device)

The main point is that such provision does not deserve any complexity adding to
the language, so burden to the user as well. Of course, who knows what the stdc++'s
head is thinking about.

On Friday, August 18, 2023 at 8:01:58 PM UTC+8, wij wrote:
> On Friday, August 18, 2023 at 8:00:49 AM UTC+8, Alf P. Steinbach wrote:
> > On 2023-08-17 10:43 PM, wij wrote:
>
> > > 3. I don't know what "fmt::print( f, "{:d}\n", i );" really mean, but in general
> > > I think cstring format in its current usage (e.g. "foo %s= %d\n") does not
> > > provide real benefit for the executing time spent interpreting it.
> > With the {fmt} library's `fmt::print`, which will become C++23's
> > `std::print`, the format string is interpreted at compile time.
> >
> > That's almost like a lint implemented via C++ template system, plus that
> > it's just about the fastest string formatting around, plus handles
> > Unicode (guaranteed by the standard) if the output device can do it.
> >
> > The design of the compile time format string handling does however not
> > adhere to the principle of least surprise, and that caused the
> > standard's wording to be a bit sabotage-like for some time. That was
> > eventually fixed. But means that somewhere I have this C++20 code:
> >
> > #if __cpp_lib_format >= (2022*100 + 7)
> > #define myprefix_FMTLIB_API_PROVIDES_FORMAT_STRING_TYPE
> > constexpr bool api_provides_format_string_type = true;
> >
> >
> > - Alf
> Not only time, but many kind of resources are used by std::print.
>
> My idea of the 'format' issue was adding a format member to class String should
> be enough, e.g. "String& String::fmt(..)", where the .. means unclear.
>
> So, the usecase (example) printf("IC0123 qty:%5d, price:%.2f\n",qty,prc) can be
> translated like the following:
>
> result << String("IC0123 qty:") << wrd(qty).fmt("5") << ", price:"
> << wrd(prc).fmt(".2") << "\n";
>
> ('resutl' can be String or any device)
>
> The main point is that such provision does not deserve any complexity adding to
> the language, so burden to the user as well. Of course, who knows what the stdc++'s
> head is thinking about.

Just came to me, libwy already had the basic capability of 'format' :

wrd(i).indent('0',5) // left align, width=5, fill-character='0';
wrd(i).indent(5,' ') // right align, width=5, fill-character=' ';
wrd(val,10,2) // float number=val, radix=10, fraction digits=2
wrd(val, 0,2) // float number=val, scientific notation, fraction digits=2

I never remember all these rules of format specifiers:
https://cplusplus.com/reference/cstdio/printf/

As to the example shown "constexpr int n = 10'000'000;", there are many kind of
preferences, e.g. "10,000,000" or "1000 0000",... or 'unicode'. IMO, these are
better be done by programmers, not in 'standard'.

I tried to combine the bitmap-approach with the heap-approach,
having a bitmap that fits into the cache and which is a window
inside the rows of numbers.

#include <iostream>
#include <vector>
#include <charconv>
#include <string_view>
#include <algorithm>
#include <fstream>
#include <cctype>
#include <cstring>
#include <queue>

#if defined(_MSC_VER)
#pragma warning(disable: 26495) // uninitialized member
#endif

using namespace std;

int main( int argc, char **argv )
{ try
{
constexpr size_t
BITMAP_SIZE = 0x80000,
BITS = BITMAP_SIZE * 8;
size_t max;
if( argc < 2 )
return EXIT_FAILURE;
char const *sizeStr = argv[1];
bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
sizeStr += 2 * hex;
if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
return EXIT_FAILURE;
union ndi_t { uint8_t c; ndi_t() {} };
vector<ndi_t> bitmap( BITMAP_SIZE );
struct next_multiple { size_t p, pNext; };
auto nm_less = []( next_multiple const &left, next_multiple const
&right ) { return left.pNext > right.pNext; };
priority_queue<next_multiple, vector<next_multiple>, decltype(nm_less)>
multiples( nm_less );
for( size_t w = 2; w <= max; w += BITS )
{
for_each( bitmap.begin(), bitmap.end(), [&]( ndi_t &n ) { n.c = 0xFF;
} );
size_t wEnd = w + BITS <= max ? w + BITS : max;
auto clearInterval = [&]( next_multiple &nm )
{
size_t b = nm.pNext - w;
for( ; b < wEnd - w; b += nm.p )
bitmap[b / 8].c &= ~(1u << b % 8);
nm.pNext = w + b;
};
if( multiples.size() ) [[likely]]
for( next_multiple const *pNm; (pNm = &multiples.top())->pNext < wEnd; )
{
auto nm = *pNm;
multiples.pop();
clearInterval( nm );
multiples.push( nm );
}
for( size_t b = 0, e = wEnd - w; b != e; ++b )
if( bitmap[b / 8].c >> b % 8 & 1 )
{
next_multiple nm { w + b, 2 * (w + b) };
clearInterval( nm );
multiples.push( nm );
}
}
/*while (multiples.size())
cout << multiples.top().p << endl,
multiples.pop();*/
}
catch( bad_alloc const & )
{
cout << "out of memory" << endl;
}

}

This is faster than the pure heap approach but still not as fast
as the pure bitmap approach.

On Mon, 21 Aug 2023 11:00:15 +0200
Bonita Montero <Bonita.Montero@gmail.com> wrote:
>I tried to combine the bitmap-approach with the heap-approach,
>having a bitmap that fits into the cache and which is a window
>inside the rows of numbers.

If only erastothenes had had C++20 he too could have come up with some
hopeless complicated solution no one understood. Or maybe he would have
done some simple baby code in C...

#include <stdio.h>
#include <stdlib.h>
#include <strings.h>
#include <math.h>

int main(int argc, char **argv)
{ char *primes;
int num;
int sqr;
int i;
int j;

if (argc != 2 || (num = atoi(argv[1])) < 1)
{
printf("Usage: %s <max number>\n",argv[0]);
return 1;
}
primes = (char *)malloc(num);

bzero(primes,num);

sqr = (int)ceil(sqrt(num));

for(i=3;i <= sqr;i+=2)
{
if (!primes[i])
{
for(j=i*2;j <= num;j+=i) primes[j] = 1;
}
}

for(i=1;i <= num;i+=2)
{
if (!primes[i])
{
printf("%d\n",i);
if (i == 1) puts("2");
}
}
return 0;
}

Am 21.08.2023 um 19:19 schrieb Muttley@dastardlyhq.com:
> On Mon, 21 Aug 2023 11:00:15 +0200
> Bonita Montero <Bonita.Montero@gmail.com> wrote:
>> I tried to combine the bitmap-approach with the heap-approach,
>> having a bitmap that fits into the cache and which is a window
>> inside the rows of numbers.
>
> If only erastothenes had had C++20 he too could have come up with some
> hopeless complicated solution no one understood. Or maybe he would have
> done some simple baby code in C...

On my PC printing the all primes below 1E8 takes 0.827s with your
code. With my code it takes about 0.36s. Take this as a reference:

#include <iostream>
#include <vector>
#include <charconv>
#include <string_view>
#include <algorithm>
#include <fstream>
#include <cctype>
#include <cstring>

#if defined(_MSC_VER)
#pragma warning(disable: 26495) // uninitialized member
#endif

using namespace std;

int main( int argc, char **argv )
{ constexpr size_t BUF_SIZE = 0x100000;
try
{
size_t max;
if( argc < 2 )
return EXIT_FAILURE;
char const *sizeStr = argv[1];
bool hex = sizeStr[0] == '0' && (sizeStr[1] == 'x' || sizeStr[1] == 'X');
sizeStr += 2 * hex;
if( from_chars_result fcr = from_chars( sizeStr, sizeStr + strlen(
sizeStr ), max, !hex ? 10 : 16 ); (bool)fcr.ec || *fcr.ptr )
return EXIT_FAILURE;
vector<uint8_t> sieve( (max + 1 + 7) / 8, (uint8_t)-1 );
ofstream ofs;
union ndi_t { char c; ndi_t() {} }; // non-default-initialzied
using ndi_vec = vector<ndi_t>;
using ndi_vec_it = ndi_vec::iterator;
constexpr size_t AHEAD = 32;
ndi_vec printBuf;
ndi_vec_it bufBegin, bufEnd, bufFlush;
if( argc < 3 || *argv[2] )
{
ofs.exceptions( ofstream::failbit | ofstream::badbit );
ofs.open( argc >= 3 ? argv[2] : "primes.txt", ofstream::trunc );
printBuf.resize( BUF_SIZE + AHEAD - 1 );
bufBegin = printBuf.begin();
bufEnd = bufBegin;
bufFlush = bufBegin + BUF_SIZE;
}
auto print = [&]() { ofs << string_view( &bufBegin->c, &bufEnd->c ); };
for( size_t p = 2; p <= max; )
{
for( size_t m = 2 * p; m <= max; m += p )
sieve[m / 8] &= ~(1u << m % 8);
while( ++p <= max && !(1 & sieve[p / 8] >> p % 8) );
if( !ofs.is_open() )
continue;
to_chars_result tcr = to_chars( &bufEnd->c, &to_address( bufEnd +
(AHEAD - 1) )->c, p );
bufEnd = tcr.ptr - &bufBegin->c + bufBegin; // NOP
bufEnd++->c = '\n';
if( bufEnd >= bufFlush )
print(),
bufEnd = bufBegin;

}
print();
return EXIT_SUCCESS;
}
catch( bad_alloc const & )
{
cout << "out of memory" << endl;
}
catch( ios_base::failure const & )
{
cout << "I/O error" << endl;
}
}

The complexity is all because of speed. Using a byte-array is much
slower and using C-streams is also slow. I'm using a bit-array and
I'm using my own output buffering.