From e4be228a5b50acd898e34cace369a533d844c635 Mon Sep 17 00:00:00 2001 From: Jens Tirsvad Nielsen Date: Tue, 20 Oct 2026 00:07:17 +0800 Subject: [PATCH 1/3] Add AVX2 prime check with cross-platform Makefile Trial division by 6k+/-1 using double-precision AVX2 lanes so the quotient is exact for 32-bit n. Makefile builds on Linux and Windows (MinGW gcc). --- Makefile | 36 ++++++++++++++++++++++++++++++++ prime_avx2.c | 59 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 95 insertions(+) create mode 100644 Makefile create mode 100644 prime_avx2.c diff --git a/Makefile b/Makefile new file mode 100644 index 0000000..d8a7449 --- /dev/null +++ b/Makefile @@ -0,0 +1,36 @@ +# Builds prime_avx2 on Linux (gcc) and Windows (MinGW-w64 / MSYS2 gcc). +# make build +# make run build and run +# make clean remove build output + +CC ?= gcc +CFLAGS ?= -O3 -mavx2 -Wall -Wextra +LDLIBS += -lm + +ifeq ($(OS),Windows_NT) + EXE := .exe + RM := del /Q /F + # make on Windows may use cmd.exe, which needs backslashes + RUN := .\prime_avx2$(EXE) + RMTARGET := prime_avx2$(EXE) +else + EXE := + RM := rm -f + RUN := ./prime_avx2 + RMTARGET := prime_avx2 +endif + +TARGET := prime_avx2$(EXE) + +.PHONY: all run clean + +all: $(TARGET) + +$(TARGET): prime_avx2.c + $(CC) $(CFLAGS) $< -o $@ $(LDLIBS) + +run: $(TARGET) + $(RUN) + +clean: + -$(RM) $(RMTARGET) diff --git a/prime_avx2.c b/prime_avx2.c new file mode 100644 index 0000000..4daec44 --- /dev/null +++ b/prime_avx2.c @@ -0,0 +1,59 @@ +#include +#include +#include +#include +#include +#include // AVX2 intrinsics + +// Returns true if n is prime, using AVX2 trial division by 6k +/- 1. +// Supports n up to UINT32_MAX. +// +// Four divisors are tested per iteration in double precision. A double holds +// every 32-bit integer exactly and IEEE division is correctly rounded, so +// floor(n / d) is exact here (a 24-bit float quotient is not). +bool is_prime_avx2(uint64_t n) { + if (n <= 1) return false; + if (n <= 3) return true; + if (n % 2 == 0 || n % 3 == 0) return false; + if (n > UINT32_MAX) return false; // out of supported range + + const double limit = floor(sqrt((double)n)); + + __m256d v_n = _mm256_set1_pd((double)n); + __m256d v_limit = _mm256_set1_pd(limit); + + // Divisors 5, 7, 11, 13 (6k +/- 1); each iteration advances by 12. + __m128i v_div_i = _mm_setr_epi32(5, 7, 11, 13); + const __m128i v_step = _mm_set1_epi32(12); + + for (double base = 5.0; base <= limit; base += 12.0) { + __m256d v_d = _mm256_cvtepi32_pd(v_div_i); + __m256d v_q = _mm256_floor_pd(_mm256_div_pd(v_n, v_d)); + __m256d v_rem = _mm256_sub_pd(v_n, _mm256_mul_pd(v_q, v_d)); + + __m256d is_zero = _mm256_cmp_pd(v_rem, _mm256_setzero_pd(), _CMP_EQ_OQ); + __m256d in_range = _mm256_cmp_pd(v_d, v_limit, _CMP_LE_OQ); + + // A lane only counts if it divides n and the divisor is <= sqrt(n) + if (_mm256_movemask_pd(_mm256_and_pd(is_zero, in_range)) != 0) + return false; + + v_div_i = _mm_add_epi32(v_div_i, v_step); + } + + return true; +} + +int main(void) { + // Known large prime: 1,000,000,007 + uint64_t num = 1000000007; + + printf("Checking if %" PRIu64 " is prime using AVX2...\n", num); + + if (is_prime_avx2(num)) + printf("%" PRIu64 " is a PRIME number.\n", num); + else + printf("%" PRIu64 " is a COMPOSITE number.\n", num); + + return 0; +} -- 2.54.0 From 42774d3d3a89dce1cca6c6a6e417b90b841186df Mon Sep 17 00:00:00 2001 From: Jens Tirsvad Nielsen Date: Tue, 20 Oct 2026 01:07:36 +0800 Subject: [PATCH 2/3] Build into dist/, list primes to 1,000,000 with timing, add wiki --- .gitignore | 1 + Makefile | 25 +++--- docs/wiki/Build-and-Run.md | 24 ++++++ docs/wiki/Description.md | 34 ++++++++ docs/wiki/Home.md | 10 +++ docs/wiki/Results.md | 14 ++++ docs/wiki/_Sidebar.md | 4 + prime_avx2.c => src/prime_avx2.c | 131 +++++++++++++++++-------------- 8 files changed, 173 insertions(+), 70 deletions(-) create mode 100644 .gitignore create mode 100644 docs/wiki/Build-and-Run.md create mode 100644 docs/wiki/Description.md create mode 100644 docs/wiki/Home.md create mode 100644 docs/wiki/Results.md create mode 100644 docs/wiki/_Sidebar.md rename prime_avx2.c => src/prime_avx2.c (65%) diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..849ddff --- /dev/null +++ b/.gitignore @@ -0,0 +1 @@ +dist/ diff --git a/Makefile b/Makefile index d8a7449..2072a51 100644 --- a/Makefile +++ b/Makefile @@ -1,36 +1,39 @@ -# Builds prime_avx2 on Linux (gcc) and Windows (MinGW-w64 / MSYS2 gcc). +# Builds prime_avx2 into dist/ on Linux (gcc) and Windows (MinGW-w64 / MSYS2 gcc). # make build # make run build and run -# make clean remove build output +# make clean remove dist/ CC ?= gcc CFLAGS ?= -O3 -mavx2 -Wall -Wextra LDLIBS += -lm +DIST := dist + ifeq ($(OS),Windows_NT) EXE := .exe - RM := del /Q /F # make on Windows may use cmd.exe, which needs backslashes - RUN := .\prime_avx2$(EXE) - RMTARGET := prime_avx2$(EXE) + RUN := $(DIST)\prime_avx2$(EXE) + RMDIR := if exist $(DIST) rmdir /S /Q $(DIST) else EXE := - RM := rm -f - RUN := ./prime_avx2 - RMTARGET := prime_avx2 + RUN := ./$(DIST)/prime_avx2 + RMDIR := rm -rf $(DIST) endif -TARGET := prime_avx2$(EXE) +TARGET := $(DIST)/prime_avx2$(EXE) .PHONY: all run clean all: $(TARGET) -$(TARGET): prime_avx2.c +$(TARGET): src/prime_avx2.c | $(DIST) $(CC) $(CFLAGS) $< -o $@ $(LDLIBS) +$(DIST): + mkdir $(DIST) + run: $(TARGET) $(RUN) clean: - -$(RM) $(RMTARGET) + -$(RMDIR) diff --git a/docs/wiki/Build-and-Run.md b/docs/wiki/Build-and-Run.md new file mode 100644 index 0000000..6535a6e --- /dev/null +++ b/docs/wiki/Build-and-Run.md @@ -0,0 +1,24 @@ +# Build and run + +**English** + +The `Makefile` is in the repository root and puts the executable in `dist/`. +It needs `gcc` (Linux, or MinGW-w64 / MSYS2 on Windows) and a CPU with AVX2. + +```bash +make +``` + +```bash +make run +``` + +```bash +make clean +``` + +To save the primes and still see the time, redirect only standard output: + +```bash +./dist/prime_avx2 > primes.txt +``` diff --git a/docs/wiki/Description.md b/docs/wiki/Description.md new file mode 100644 index 0000000..ebf60a4 --- /dev/null +++ b/docs/wiki/Description.md @@ -0,0 +1,34 @@ +# Description + +**English** + +## Purpose + +The program shows how fast a CPU is. How fast it looks depends on how the +program is written, so the number is a measure of the CPU and the code +together, not of the CPU alone. + +## What it does + +It tests every number from 0 to 1,000,000 with `is_prime_avx2`, which uses +AVX2 trial division by 6k ± 1 (four divisors per step). It prints each prime +on its own line and writes a summary line to standard error: + +```text +Found 78498 primes from 0 to 1000000 in 13.806 ms +``` + +78,498 is the known count of primes up to 1,000,000, so the summary line also +works as a correctness check. + +## Why the time depends on the code + +The same CPU gives very different times for different programs: + +- the algorithm (trial division, a sieve, a segmented sieve), +- the instruction set used (scalar code, AVX2, AVX-512), +- compiler and flags (`-O3 -mavx2`, MSVC `/O2 /arch:AVX2`), +- what is measured (only the search, or also printing to a console). + +Here only the search is timed. Output is collected first and printed +afterwards, so console speed does not change the number. diff --git a/docs/wiki/Home.md b/docs/wiki/Home.md new file mode 100644 index 0000000..e3e4085 --- /dev/null +++ b/docs/wiki/Home.md @@ -0,0 +1,10 @@ +# prime-numbers-c + +**English** + +A small C program that lists every prime number from 0 to 1,000,000 and +reports how long it took. + +- [Description](Description): what the program is for +- [Results](Results): measured times on real hardware +- [Build and run](Build-and-Run): `make`, `make run` diff --git a/docs/wiki/Results.md b/docs/wiki/Results.md new file mode 100644 index 0000000..c07f493 --- /dev/null +++ b/docs/wiki/Results.md @@ -0,0 +1,14 @@ +# Results + +**English** + +Time to find all primes from 0 to 1,000,000. Only the search is timed. + +| CPU | Result | +| --- | --- | +| 13th Gen Intel(R) Core(TM) i5-13600K | Found 78498 primes from 0 to 1000000 in 20.977 ms | + +A result is only comparable with another when the code, compiler and flags +are the same. To add a result, run the program on your machine (see +[Build and run](Build-and-Run)) and add a row with your CPU name and the +summary line. diff --git a/docs/wiki/_Sidebar.md b/docs/wiki/_Sidebar.md new file mode 100644 index 0000000..b002eee --- /dev/null +++ b/docs/wiki/_Sidebar.md @@ -0,0 +1,4 @@ +- [Home](Home) +- [Description](Description) +- [Results](Results) +- [Build and run](Build-and-Run) diff --git a/prime_avx2.c b/src/prime_avx2.c similarity index 65% rename from prime_avx2.c rename to src/prime_avx2.c index 4daec44..b6216fd 100644 --- a/prime_avx2.c +++ b/src/prime_avx2.c @@ -1,59 +1,72 @@ -#include -#include -#include -#include -#include -#include // AVX2 intrinsics - -// Returns true if n is prime, using AVX2 trial division by 6k +/- 1. -// Supports n up to UINT32_MAX. -// -// Four divisors are tested per iteration in double precision. A double holds -// every 32-bit integer exactly and IEEE division is correctly rounded, so -// floor(n / d) is exact here (a 24-bit float quotient is not). -bool is_prime_avx2(uint64_t n) { - if (n <= 1) return false; - if (n <= 3) return true; - if (n % 2 == 0 || n % 3 == 0) return false; - if (n > UINT32_MAX) return false; // out of supported range - - const double limit = floor(sqrt((double)n)); - - __m256d v_n = _mm256_set1_pd((double)n); - __m256d v_limit = _mm256_set1_pd(limit); - - // Divisors 5, 7, 11, 13 (6k +/- 1); each iteration advances by 12. - __m128i v_div_i = _mm_setr_epi32(5, 7, 11, 13); - const __m128i v_step = _mm_set1_epi32(12); - - for (double base = 5.0; base <= limit; base += 12.0) { - __m256d v_d = _mm256_cvtepi32_pd(v_div_i); - __m256d v_q = _mm256_floor_pd(_mm256_div_pd(v_n, v_d)); - __m256d v_rem = _mm256_sub_pd(v_n, _mm256_mul_pd(v_q, v_d)); - - __m256d is_zero = _mm256_cmp_pd(v_rem, _mm256_setzero_pd(), _CMP_EQ_OQ); - __m256d in_range = _mm256_cmp_pd(v_d, v_limit, _CMP_LE_OQ); - - // A lane only counts if it divides n and the divisor is <= sqrt(n) - if (_mm256_movemask_pd(_mm256_and_pd(is_zero, in_range)) != 0) - return false; - - v_div_i = _mm_add_epi32(v_div_i, v_step); - } - - return true; -} - -int main(void) { - // Known large prime: 1,000,000,007 - uint64_t num = 1000000007; - - printf("Checking if %" PRIu64 " is prime using AVX2...\n", num); - - if (is_prime_avx2(num)) - printf("%" PRIu64 " is a PRIME number.\n", num); - else - printf("%" PRIu64 " is a COMPOSITE number.\n", num); - - return 0; -} +#include +#include +#include +#include +#include +#include +#include // AVX2 intrinsics + +// Returns true if n is prime, using AVX2 trial division by 6k +/- 1. +// Supports n up to UINT32_MAX. +// +// Four divisors are tested per iteration in double precision. A double holds +// every 32-bit integer exactly and IEEE division is correctly rounded, so +// floor(n / d) is exact here (a 24-bit float quotient is not). +bool is_prime_avx2(uint64_t n) { + if (n <= 1) return false; + if (n <= 3) return true; + if (n % 2 == 0 || n % 3 == 0) return false; + if (n > UINT32_MAX) return false; // out of supported range + + const double limit = floor(sqrt((double)n)); + + __m256d v_n = _mm256_set1_pd((double)n); + __m256d v_limit = _mm256_set1_pd(limit); + + // Divisors 5, 7, 11, 13 (6k +/- 1); each iteration advances by 12. + __m128i v_div_i = _mm_setr_epi32(5, 7, 11, 13); + const __m128i v_step = _mm_set1_epi32(12); + + for (double base = 5.0; base <= limit; base += 12.0) { + __m256d v_d = _mm256_cvtepi32_pd(v_div_i); + __m256d v_q = _mm256_floor_pd(_mm256_div_pd(v_n, v_d)); + __m256d v_rem = _mm256_sub_pd(v_n, _mm256_mul_pd(v_q, v_d)); + + __m256d is_zero = _mm256_cmp_pd(v_rem, _mm256_setzero_pd(), _CMP_EQ_OQ); + __m256d in_range = _mm256_cmp_pd(v_d, v_limit, _CMP_LE_OQ); + + // A lane only counts if it divides n and the divisor is <= sqrt(n) + if (_mm256_movemask_pd(_mm256_and_pd(is_zero, in_range)) != 0) + return false; + + v_div_i = _mm_add_epi32(v_div_i, v_step); + } + + return true; +} + +#define LIMIT 1000000u + +// Seconds elapsed between two timespec_get(TIME_UTC) samples. +static double elapsed(const struct timespec *a, const struct timespec *b) { + return (double)(b->tv_sec - a->tv_sec) + (double)(b->tv_nsec - a->tv_nsec) / 1e9; +} + +int main(void) { + static uint32_t primes[LIMIT / 2]; // pi(1e6) = 78498, so this is ample + size_t count = 0; + struct timespec t0, t1; + + // Time the computation only; printing is done afterwards. + timespec_get(&t0, TIME_UTC); + for (uint32_t n = 0; n <= LIMIT; n++) + if (is_prime_avx2(n)) primes[count++] = n; + timespec_get(&t1, TIME_UTC); + + for (size_t i = 0; i < count; i++) + printf("%" PRIu32 "\n", primes[i]); + + fprintf(stderr, "Found %zu primes from 0 to %u in %.3f ms\n", + count, LIMIT, elapsed(&t0, &t1) * 1000.0); + return 0; +} -- 2.54.0 From 8c5359a56741004cff5ce0ba99727fcd76cfdd1f Mon Sep 17 00:00:00 2001 From: Jens Tirsvad Nielsen Date: Tue, 20 Oct 2026 01:07:51 +0800 Subject: [PATCH 3/3] Use LF line endings in src/prime_avx2.c and Makefile --- src/prime_avx2.c | 144 +++++++++++++++++++++++------------------------ 1 file changed, 72 insertions(+), 72 deletions(-) diff --git a/src/prime_avx2.c b/src/prime_avx2.c index b6216fd..a54d76e 100644 --- a/src/prime_avx2.c +++ b/src/prime_avx2.c @@ -1,72 +1,72 @@ -#include -#include -#include -#include -#include -#include -#include // AVX2 intrinsics - -// Returns true if n is prime, using AVX2 trial division by 6k +/- 1. -// Supports n up to UINT32_MAX. -// -// Four divisors are tested per iteration in double precision. A double holds -// every 32-bit integer exactly and IEEE division is correctly rounded, so -// floor(n / d) is exact here (a 24-bit float quotient is not). -bool is_prime_avx2(uint64_t n) { - if (n <= 1) return false; - if (n <= 3) return true; - if (n % 2 == 0 || n % 3 == 0) return false; - if (n > UINT32_MAX) return false; // out of supported range - - const double limit = floor(sqrt((double)n)); - - __m256d v_n = _mm256_set1_pd((double)n); - __m256d v_limit = _mm256_set1_pd(limit); - - // Divisors 5, 7, 11, 13 (6k +/- 1); each iteration advances by 12. - __m128i v_div_i = _mm_setr_epi32(5, 7, 11, 13); - const __m128i v_step = _mm_set1_epi32(12); - - for (double base = 5.0; base <= limit; base += 12.0) { - __m256d v_d = _mm256_cvtepi32_pd(v_div_i); - __m256d v_q = _mm256_floor_pd(_mm256_div_pd(v_n, v_d)); - __m256d v_rem = _mm256_sub_pd(v_n, _mm256_mul_pd(v_q, v_d)); - - __m256d is_zero = _mm256_cmp_pd(v_rem, _mm256_setzero_pd(), _CMP_EQ_OQ); - __m256d in_range = _mm256_cmp_pd(v_d, v_limit, _CMP_LE_OQ); - - // A lane only counts if it divides n and the divisor is <= sqrt(n) - if (_mm256_movemask_pd(_mm256_and_pd(is_zero, in_range)) != 0) - return false; - - v_div_i = _mm_add_epi32(v_div_i, v_step); - } - - return true; -} - -#define LIMIT 1000000u - -// Seconds elapsed between two timespec_get(TIME_UTC) samples. -static double elapsed(const struct timespec *a, const struct timespec *b) { - return (double)(b->tv_sec - a->tv_sec) + (double)(b->tv_nsec - a->tv_nsec) / 1e9; -} - -int main(void) { - static uint32_t primes[LIMIT / 2]; // pi(1e6) = 78498, so this is ample - size_t count = 0; - struct timespec t0, t1; - - // Time the computation only; printing is done afterwards. - timespec_get(&t0, TIME_UTC); - for (uint32_t n = 0; n <= LIMIT; n++) - if (is_prime_avx2(n)) primes[count++] = n; - timespec_get(&t1, TIME_UTC); - - for (size_t i = 0; i < count; i++) - printf("%" PRIu32 "\n", primes[i]); - - fprintf(stderr, "Found %zu primes from 0 to %u in %.3f ms\n", - count, LIMIT, elapsed(&t0, &t1) * 1000.0); - return 0; -} +#include +#include +#include +#include +#include +#include +#include // AVX2 intrinsics + +// Returns true if n is prime, using AVX2 trial division by 6k +/- 1. +// Supports n up to UINT32_MAX. +// +// Four divisors are tested per iteration in double precision. A double holds +// every 32-bit integer exactly and IEEE division is correctly rounded, so +// floor(n / d) is exact here (a 24-bit float quotient is not). +bool is_prime_avx2(uint64_t n) { + if (n <= 1) return false; + if (n <= 3) return true; + if (n % 2 == 0 || n % 3 == 0) return false; + if (n > UINT32_MAX) return false; // out of supported range + + const double limit = floor(sqrt((double)n)); + + __m256d v_n = _mm256_set1_pd((double)n); + __m256d v_limit = _mm256_set1_pd(limit); + + // Divisors 5, 7, 11, 13 (6k +/- 1); each iteration advances by 12. + __m128i v_div_i = _mm_setr_epi32(5, 7, 11, 13); + const __m128i v_step = _mm_set1_epi32(12); + + for (double base = 5.0; base <= limit; base += 12.0) { + __m256d v_d = _mm256_cvtepi32_pd(v_div_i); + __m256d v_q = _mm256_floor_pd(_mm256_div_pd(v_n, v_d)); + __m256d v_rem = _mm256_sub_pd(v_n, _mm256_mul_pd(v_q, v_d)); + + __m256d is_zero = _mm256_cmp_pd(v_rem, _mm256_setzero_pd(), _CMP_EQ_OQ); + __m256d in_range = _mm256_cmp_pd(v_d, v_limit, _CMP_LE_OQ); + + // A lane only counts if it divides n and the divisor is <= sqrt(n) + if (_mm256_movemask_pd(_mm256_and_pd(is_zero, in_range)) != 0) + return false; + + v_div_i = _mm_add_epi32(v_div_i, v_step); + } + + return true; +} + +#define LIMIT 1000000u + +// Seconds elapsed between two timespec_get(TIME_UTC) samples. +static double elapsed(const struct timespec *a, const struct timespec *b) { + return (double)(b->tv_sec - a->tv_sec) + (double)(b->tv_nsec - a->tv_nsec) / 1e9; +} + +int main(void) { + static uint32_t primes[LIMIT / 2]; // pi(1e6) = 78498, so this is ample + size_t count = 0; + struct timespec t0, t1; + + // Time the computation only; printing is done afterwards. + timespec_get(&t0, TIME_UTC); + for (uint32_t n = 0; n <= LIMIT; n++) + if (is_prime_avx2(n)) primes[count++] = n; + timespec_get(&t1, TIME_UTC); + + for (size_t i = 0; i < count; i++) + printf("%" PRIu32 "\n", primes[i]); + + fprintf(stderr, "Found %zu primes from 0 to %u in %.3f ms\n", + count, LIMIT, elapsed(&t0, &t1) * 1000.0); + return 0; +} -- 2.54.0