From 67babd37ab140e00abdf491111523cd0f915e309 Mon Sep 17 00:00:00 2001
From: Heisenberg Vader <23ucs596@lnmiit.ac.in>
Date: Fri, 10 Oct 2025 11:43:38 +0000
Subject: [PATCH 1/6] Add solution for Problem 810

---
 project_euler/problem_810/__init__.py |   0
 project_euler/problem_810/sol1.py     | 168 ++++++++++++++++++++++++++
 2 files changed, 168 insertions(+)
 create mode 100644 project_euler/problem_810/__init__.py
 create mode 100644 project_euler/problem_810/sol1.py

diff --git a/project_euler/problem_810/__init__.py b/project_euler/problem_810/__init__.py
new file mode 100644
index 000000000000..e69de29bb2d1
diff --git a/project_euler/problem_810/sol1.py b/project_euler/problem_810/sol1.py
new file mode 100644
index 000000000000..c293e89f859a
--- /dev/null
+++ b/project_euler/problem_810/sol1.py
@@ -0,0 +1,168 @@
+"""
+Project Euler Problem 810: https://projecteuler.net/problem=810
+
+We use x ⊕ y for the bitwise XOR of x and y.
+Define the XOR-product of x and y, denoted by x ⊗ y,
+similar to a long multiplication in base 2,
+except the intermediate results are XORed instead of usual integer addition.
+
+For example, 7 ⊗ 3 = 9, or in base 2, 111_2 ⊗ 11_2 = 1001_2:
+        111
+    ⊗   11
+    -------
+        111
+       111
+    -------
+       1001
+An XOR-Prime is an integer n greater than 1 that is not an
+XOR-product of two integers greater than 1.
+The above example shows that 9 is not an XOR-prime.
+Similarly, 5 = 3 ⊗ 3 is not an XOR-prime.
+The first few XOR-primes are 2, 3, 7, 11, 13, ... and the 10th XOR-prime is 41.
+
+Find the 5,000,000.th XOR-prime.
+
+References:
+http://en.wikipedia.org/wiki/M%C3%B6bius_function
+
+"""
+
+from bitarray import bitarray
+
+
+def xor_multiply(a: int, b: int) -> int:
+    """
+    Perform XOR multiplication of two integers, equivalent to polynomial
+    multiplication in GF(2).
+
+    >>> xor_multiply(3, 5)  # (011) ⊗ (101)
+    15
+    """
+    res = 0
+    while b:
+        if b & 1:
+            res ^= a
+        a <<= 1
+        b >>= 1
+    return res
+
+
+def divisors(n: int) -> set[int]:
+    """
+    Return all divisors of n (excluding 0).
+
+    >>> divisors(12)
+    {1, 2, 3, 4, 6, 12}
+    """
+    s = {1}
+    for i in range(2, int(n**0.5) + 1):
+        if n % i == 0:
+            s.add(i)
+            s.add(n // i)
+    s.add(n)
+    return s
+
+
+def mobius_table(n: int) -> list[int]:
+    """
+    Generate Möbius function values from 1 to n.
+
+    >>> mobius_table(10)[:6]
+    [0, 1, -1, -1, 0, -1]
+    """
+    mob = [1] * (n + 1)
+    is_prime = [True] * (n + 1)
+    mob[0] = 0
+
+    for p in range(2, n + 1):
+        if is_prime[p]:
+            mob[p] = -1
+            for j in range(2 * p, n + 1, p):
+                is_prime[j] = False
+                mob[j] *= -1
+            p2 = p * p
+            if p2 <= n:
+                for j in range(p2, n + 1, p2):
+                    mob[j] = 0
+    return mob
+
+
+def count_irreducibles(d: int) -> int:
+    """
+    Count the number of irreducible polynomials of degree d over GF(2)
+    using the Möbius function.
+
+    Args:
+        d (int): Polynomial degree.
+
+    Returns:
+        int: Count of irreducible polynomials.
+
+    Example:
+        >>> count_irreducibles(3)
+        2
+    """
+    mob = mobius_table(d)
+    total = 0
+    for div in divisors(d) | {d}:
+        total += mob[div] * (1 << (d // div))
+    return total // d
+
+
+def find_xor_prime(rank: int) -> int:
+    """
+    Find the Nth XOR prime using a bitarray-based sieve.
+
+    >>> find_xor_prime(10)
+    41
+    """
+    total, degree = 0, 1
+    while True:
+        count = count_irreducibles(degree)
+        if total + count > rank:
+            break
+        total += count
+        degree += 1
+
+    limit = 1 << (degree + 1)
+
+    sieve = bitarray(limit)
+    sieve.setall(True)
+    sieve[0:2] = False
+
+    for even in range(4, limit, 2):
+        sieve[even] = False
+
+    current = 1
+    for i in range(3, limit, 2):
+        if sieve[i]:
+            current += 1
+            if current == rank:
+                return i
+
+            j = i
+            while True:
+                prod = xor_multiply(i, j)
+                if prod >= limit:
+                    break
+                sieve[prod] = False
+                j += 2
+
+    raise ValueError(
+        "Failed to locate the requested XOR-prime within the computed limit"
+    )
+
+
+def solution(limit: int = 5000001) -> int:
+    """
+    Wrapper for Project Euler-style solution function.
+
+    >>> solution(10)
+    41
+    """
+    result = find_xor_prime(limit)
+    return result
+
+
+if __name__ == "__main__":
+    print(f"{solution(5000000) = }")

From 743830afae36dbe590618b413bf7ce96db1eae7d Mon Sep 17 00:00:00 2001
From: Heisenberg Vader <23ucs596@lnmiit.ac.in>
Date: Fri, 10 Oct 2025 16:18:02 +0000
Subject: [PATCH 2/6] Add solution for problem 810

---
 project_euler/problem_810/sol1.py | 18 ++++--------------
 1 file changed, 4 insertions(+), 14 deletions(-)

diff --git a/project_euler/problem_810/sol1.py b/project_euler/problem_810/sol1.py
index c293e89f859a..5201d563b460 100644
--- a/project_euler/problem_810/sol1.py
+++ b/project_euler/problem_810/sol1.py
@@ -27,8 +27,6 @@
 
 """
 
-from bitarray import bitarray
-
 
 def xor_multiply(a: int, b: int) -> int:
     """
@@ -92,15 +90,8 @@ def count_irreducibles(d: int) -> int:
     Count the number of irreducible polynomials of degree d over GF(2)
     using the Möbius function.
 
-    Args:
-        d (int): Polynomial degree.
-
-    Returns:
-        int: Count of irreducible polynomials.
-
-    Example:
-        >>> count_irreducibles(3)
-        2
+    >>> count_irreducibles(3)
+    2
     """
     mob = mobius_table(d)
     total = 0
@@ -126,9 +117,8 @@ def find_xor_prime(rank: int) -> int:
 
     limit = 1 << (degree + 1)
 
-    sieve = bitarray(limit)
-    sieve.setall(True)
-    sieve[0:2] = False
+    sieve = [True] * limit
+    sieve[0] = sieve[1] = False
 
     for even in range(4, limit, 2):
         sieve[even] = False

From aa30a0b57d774c6f90787a3cb1d911fe0ff5082f Mon Sep 17 00:00:00 2001
From: Heisenberg Vader <23ucs596@lnmiit.ac.in>
Date: Fri, 10 Oct 2025 16:32:30 +0000
Subject: [PATCH 3/6] Add solution for Problem 810

---
 project_euler/problem_810/sol1.py | 40 +++++++++++++++----------------
 1 file changed, 20 insertions(+), 20 deletions(-)

diff --git a/project_euler/problem_810/sol1.py b/project_euler/problem_810/sol1.py
index 5201d563b460..c6730cd6569e 100644
--- a/project_euler/problem_810/sol1.py
+++ b/project_euler/problem_810/sol1.py
@@ -28,7 +28,7 @@
 """
 
 
-def xor_multiply(a: int, b: int) -> int:
+def xor_multiply(op_a: int, op_b: int) -> int:
     """
     Perform XOR multiplication of two integers, equivalent to polynomial
     multiplication in GF(2).
@@ -37,33 +37,33 @@ def xor_multiply(a: int, b: int) -> int:
     15
     """
     res = 0
-    while b:
-        if b & 1:
-            res ^= a
-        a <<= 1
-        b >>= 1
+    while op_b:
+        if op_b & 1:
+            res ^= op_a
+        op_a <<= 1
+        op_b >>= 1
     return res
 
 
-def divisors(n: int) -> set[int]:
+def divisors(num: int) -> set[int]:
     """
-    Return all divisors of n (excluding 0).
+    Return all divisors of `num` (excluding 0).
 
     >>> divisors(12)
     {1, 2, 3, 4, 6, 12}
     """
     s = {1}
-    for i in range(2, int(n**0.5) + 1):
-        if n % i == 0:
+    for i in range(2, int(num**0.5) + 1):
+        if num % i == 0:
             s.add(i)
-            s.add(n // i)
-    s.add(n)
+            s.add(num // i)
+    s.add(num)
     return s
 
 
 def mobius_table(n: int) -> list[int]:
     """
-    Generate Möbius function values from 1 to n.
+    Generate a variant of Möbius function values from 1 to n.
 
     >>> mobius_table(10)[:6]
     [0, 1, -1, -1, 0, -1]
@@ -85,19 +85,19 @@ def mobius_table(n: int) -> list[int]:
     return mob
 
 
-def count_irreducibles(d: int) -> int:
+def count_irreducibles(deg: int) -> int:
     """
-    Count the number of irreducible polynomials of degree d over GF(2)
-    using the Möbius function.
+    Count the number of irreducible polynomials of degree `deg` over GF(2)
+    using the variant of Möbius function.
 
     >>> count_irreducibles(3)
     2
     """
-    mob = mobius_table(d)
+    mob = mobius_table(deg)
     total = 0
-    for div in divisors(d) | {d}:
-        total += mob[div] * (1 << (d // div))
-    return total // d
+    for div in divisors(deg) | {deg}:
+        total += mob[div] * (1 << (deg // div))
+    return total // deg
 
 
 def find_xor_prime(rank: int) -> int:

From d88cb320f4aea310277e09bedbd3fc30cd833ba6 Mon Sep 17 00:00:00 2001
From: Heisenberg Vader <23ucs596@lnmiit.ac.in>
Date: Fri, 10 Oct 2025 16:33:58 +0000
Subject: [PATCH 4/6] Add solution for Problem 810

---
 project_euler/problem_810/sol1.py | 16 ++++++++--------
 1 file changed, 8 insertions(+), 8 deletions(-)

diff --git a/project_euler/problem_810/sol1.py b/project_euler/problem_810/sol1.py
index c6730cd6569e..964b84db1fbe 100644
--- a/project_euler/problem_810/sol1.py
+++ b/project_euler/problem_810/sol1.py
@@ -61,26 +61,26 @@ def divisors(num: int) -> set[int]:
     return s
 
 
-def mobius_table(n: int) -> list[int]:
+def mobius_table(num: int) -> list[int]:
     """
-    Generate a variant of Möbius function values from 1 to n.
+    Generate a variant of Möbius function values from 1 to num.
 
     >>> mobius_table(10)[:6]
     [0, 1, -1, -1, 0, -1]
     """
-    mob = [1] * (n + 1)
-    is_prime = [True] * (n + 1)
+    mob = [1] * (num + 1)
+    is_prime = [True] * (num + 1)
     mob[0] = 0
 
-    for p in range(2, n + 1):
+    for p in range(2, num + 1):
         if is_prime[p]:
             mob[p] = -1
-            for j in range(2 * p, n + 1, p):
+            for j in range(2 * p, num + 1, p):
                 is_prime[j] = False
                 mob[j] *= -1
             p2 = p * p
-            if p2 <= n:
-                for j in range(p2, n + 1, p2):
+            if p2 <= num:
+                for j in range(p2, num + 1, p2):
                     mob[j] = 0
     return mob
 

From ad1e62acf365b2825a6f982333412f153281a183 Mon Sep 17 00:00:00 2001
From: Heisenberg Vader <23ucs596@lnmiit.ac.in>
Date: Fri, 10 Oct 2025 17:27:07 +0000
Subject: [PATCH 5/6] Add solution for Problem 810

---
 project_euler/problem_810/sol1.py | 28 +++++++++++-----------------
 1 file changed, 11 insertions(+), 17 deletions(-)

diff --git a/project_euler/problem_810/sol1.py b/project_euler/problem_810/sol1.py
index 964b84db1fbe..e748626bfabb 100644
--- a/project_euler/problem_810/sol1.py
+++ b/project_euler/problem_810/sol1.py
@@ -27,6 +27,8 @@
 
 """
 
+from array import array
+
 
 def xor_multiply(op_a: int, op_b: int) -> int:
     """
@@ -108,23 +110,17 @@ def find_xor_prime(rank: int) -> int:
     41
     """
     total, degree = 0, 1
-    while True:
-        count = count_irreducibles(degree)
-        if total + count > rank:
-            break
-        total += count
+    while total + count_irreducibles(degree) < rank:
+        total += count_irreducibles(degree)
         degree += 1
 
     limit = 1 << (degree + 1)
 
-    sieve = [True] * limit
-    sieve[0] = sieve[1] = False
-
-    for even in range(4, limit, 2):
-        sieve[even] = False
+    sieve = array("B", [1]) * limit
+    sieve[0] = sieve[1] = 0
 
-    current = 1
-    for i in range(3, limit, 2):
+    current = 0
+    for i in range(2, limit):
         if sieve[i]:
             current += 1
             if current == rank:
@@ -135,12 +131,10 @@ def find_xor_prime(rank: int) -> int:
                 prod = xor_multiply(i, j)
                 if prod >= limit:
                     break
-                sieve[prod] = False
-                j += 2
+                sieve[prod] = 0
+                j += 1
 
-    raise ValueError(
-        "Failed to locate the requested XOR-prime within the computed limit"
-    )
+    raise ValueError("Failed to locate the requested XOR-prime")
 
 
 def solution(limit: int = 5000001) -> int:

From 590e47b0d0741e674cd29964f6f7e335ab4cf3ce Mon Sep 17 00:00:00 2001
From: Heisenberg Vader <23ucs596@lnmiit.ac.in>
Date: Mon, 13 Oct 2025 09:13:36 +0000
Subject: [PATCH 6/6] Add solution for problem 810

---
 project_euler/problem_810/sol1.py | 155 +++++++++++++++---------------
 1 file changed, 77 insertions(+), 78 deletions(-)

diff --git a/project_euler/problem_810/sol1.py b/project_euler/problem_810/sol1.py
index e748626bfabb..c47dae8a5e42 100644
--- a/project_euler/problem_810/sol1.py
+++ b/project_euler/problem_810/sol1.py
@@ -20,133 +20,132 @@
 Similarly, 5 = 3 ⊗ 3 is not an XOR-prime.
 The first few XOR-primes are 2, 3, 7, 11, 13, ... and the 10th XOR-prime is 41.
 
-Find the 5,000,000.th XOR-prime.
+Find the 5,000,000 th XOR-prime.
 
 References:
 http://en.wikipedia.org/wiki/M%C3%B6bius_function
 
 """
 
-from array import array
+import math
 
 
-def xor_multiply(op_a: int, op_b: int) -> int:
+def get_divisors(num: int) -> set[int]:
+    """
+    Return all positive divisors of num.
+    >>> get_divisors(12)
+    {1, 2, 3, 4, 6, 12}
     """
-    Perform XOR multiplication of two integers, equivalent to polynomial
-    multiplication in GF(2).
+    divisors = {1}
+    for i in range(2, int(math.sqrt(num)) + 1):
+        if num % i == 0:
+            divisors.add(i)
+            divisors.add(num // i)
+    divisors.add(num)
+    return divisors
+
 
-    >>> xor_multiply(3, 5)  # (011) ⊗ (101)
+def xor_multiply(op_a: int, op_b: int) -> int:
+    """
+    Perform XOR-based multiplication (polynomial multiplication mod 2).
+    >>> xor_multiply(3, 5)
     15
     """
-    res = 0
+    result = 0
     while op_b:
         if op_b & 1:
-            res ^= op_a
+            result ^= op_a
         op_a <<= 1
         op_b >>= 1
-    return res
-
-
-def divisors(num: int) -> set[int]:
-    """
-    Return all divisors of `num` (excluding 0).
-
-    >>> divisors(12)
-    {1, 2, 3, 4, 6, 12}
-    """
-    s = {1}
-    for i in range(2, int(num**0.5) + 1):
-        if num % i == 0:
-            s.add(i)
-            s.add(num // i)
-    s.add(num)
-    return s
+    return result
 
 
-def mobius_table(num: int) -> list[int]:
+def mobius_table(lim: int, k: int = 2) -> list[int]:
     """
-    Generate a variant of Möbius function values from 1 to num.
-
+    Compute a modified Mobius function table up to `lim`.
     >>> mobius_table(10)[:6]
     [0, 1, -1, -1, 0, -1]
     """
-    mob = [1] * (num + 1)
-    is_prime = [True] * (num + 1)
-    mob[0] = 0
+    mob = [0] + [1] * lim
+    is_prime = [True] * (lim + 1)
+    is_prime[0] = is_prime[1] = False
 
-    for p in range(2, num + 1):
+    for p in range(2, lim + 1):
         if is_prime[p]:
-            mob[p] = -1
-            for j in range(2 * p, num + 1, p):
-                is_prime[j] = False
-                mob[j] *= -1
-            p2 = p * p
-            if p2 <= num:
-                for j in range(p2, num + 1, p2):
-                    mob[j] = 0
+            mob[p] *= -1
+            for mul in range(2 * p, lim + 1, p):
+                is_prime[mul] = False
+                mob[mul] *= -1
+
+            p_pow = p**k
+            if p_pow <= lim:
+                for mul in range(p_pow, lim + 1, p_pow):
+                    mob[mul] = 0
     return mob
 
 
-def count_irreducibles(deg: int) -> int:
+def cnt_irred_pol(num: int) -> int:
     """
-    Count the number of irreducible polynomials of degree `deg` over GF(2)
-    using the variant of Möbius function.
-
-    >>> count_irreducibles(3)
+    Return the number of monic irreducible polynomials of degree num over GF(2).
+    >>> cnt_irred_pol(3)
     2
     """
-    mob = mobius_table(deg)
-    total = 0
-    for div in divisors(deg) | {deg}:
-        total += mob[div] * (1 << (deg // div))
-    return total // deg
+    mob = mobius_table(num)
+    total = sum(mob[d] * (2 ** (num // d)) for d in get_divisors(num))
+    return total // num
 
 
-def find_xor_prime(rank: int) -> int:
+def xor_prime_func(tgt_idx: int) -> int:
     """
-    Find the Nth XOR prime using a bitarray-based sieve.
-
-    >>> find_xor_prime(10)
+    Find the N-th XOR-prime (irreducible polynomial) index approximation.
+    >>> xor_prime_func(10)
     41
     """
     total, degree = 0, 1
-    while total + count_irreducibles(degree) < rank:
-        total += count_irreducibles(degree)
+
+    while True:
+        cnt = cnt_irred_pol(degree)
+        if total + cnt > tgt_idx:
+            break
+        total += cnt
         degree += 1
 
-    limit = 1 << (degree + 1)
+    lim = 1 << (degree + 1)
+    is_prime = [True] * lim
+    is_prime[0] = is_prime[1] = False
 
-    sieve = array("B", [1]) * limit
-    sieve[0] = sieve[1] = 0
+    for even in range(4, lim, 2):
+        is_prime[even] = False
 
-    current = 0
-    for i in range(2, limit):
-        if sieve[i]:
-            current += 1
-            if current == rank:
-                return i
+    cnt = 1
+    for num in range(3, lim, 2):
+        if not is_prime[num]:
+            continue
 
-            j = i
-            while True:
-                prod = xor_multiply(i, j)
-                if prod >= limit:
-                    break
-                sieve[prod] = 0
-                j += 1
+        cnt += 1
+        if cnt == tgt_idx:
+            return num
 
-    raise ValueError("Failed to locate the requested XOR-prime")
+        mul = num
+        while True:
+            prod = xor_multiply(mul, num)
+            if prod >= lim:
+                break
+            is_prime[prod] = False
+            mul += 2
 
+    raise ValueError("Could not compute the XOR-prime.")
 
-def solution(limit: int = 5000001) -> int:
-    """
-    Wrapper for Project Euler-style solution function.
 
+def solution(nth: int = 5000000) -> int:
+    """
+    Compute the Nth XOR prime and print timing.
     >>> solution(10)
     41
     """
-    result = find_xor_prime(limit)
+    result = xor_prime_func(nth)
     return result
 
 
 if __name__ == "__main__":
-    print(f"{solution(5000000) = }")
+    print(f"{solution() = }")