Upload problem_solvers/ktuple_constants.py

problem_solvers/ktuple_constants.py

@@ -0,0 +1,97 @@
+"""
+PROBLEM 5: Hardy-Littlewood k-Tuple Constants
+===============================================
+Verify the density of prime k-tuples against Hardy-Littlewood predictions.
+"""
+
+import numpy as np
+from typing import Dict, List
+
+
+class KTupleAnalyzer:
+    def __init__(self):
+        self.results = {}
+
+    def sieve_primes(self, limit: int) -> List[int]:
+        sieve = np.ones(limit + 1, dtype=bool)
+        sieve[:2] = False
+        for i in range(2, int(limit ** 0.5) + 1):
+            if sieve[i]:
+                sieve[i * i::i] = False
+        return np.where(sieve)[0].tolist()
+
+    def count_patterns(self, primes: List[int], pattern: List[int]) -> int:
+        """Count occurrences of a k-tuple pattern."""
+        prime_set = set(primes)
+        count = 0
+        for p in primes:
+            if all((p + d) in prime_set for d in pattern):
+                count += 1
+        return count
+
+    def hardy_littlewood_constant(self, pattern: List[int]) -> float:
+        """
+        Compute the Hardy-Littlewood constant 𝔖(H) for pattern H.
+        𝔖(H) = ∏_p (1 - 1/p)^{-k} * (1 - ν_H(p)/p)
+        where ν_H(p) is number of distinct residues mod p occupied by H.
+        Approximate for small p.
+        """
+        from sympy import primerange, Rational
+        k = len(pattern)
+        product = 1.0
+        for p in primerange(2, 100):
+            residues = set((d % p) for d in pattern)
+            nu_H = len(residues)
+            if nu_H == p:
+                continue  # Pattern impossible mod p
+            term = (1.0 - 1.0 / p) ** (-k) * (1.0 - nu_H / p)
+            product *= term
+        return product
+
+    def analyze(self, limit: int = 2_000_000) -> Dict:
+        primes = self.sieve_primes(limit)
+        prime_set = set(primes)
+        n_primes = len(primes)
+        x = limit
+
+        patterns = {
+            'twin': [0, 2],
+            'cousin': [0, 4],
+            'sexy': [0, 6],
+            'triplet1': [0, 2, 6],
+            'triplet2': [0, 4, 6],
+            'quadruplet': [0, 2, 6, 8],
+        }
+
+        results = {}
+        for name, pat in patterns.items():
+            count = self.count_patterns(primes, pat)
+            k = len(pat)
+            hl_const = self.hardy_littlewood_constant(pat)
+            # Predicted: count ~ 𝔖(H) · x / (log x)^k
+            predicted = hl_const * x / (np.log(x) ** k)
+            error = abs(count - predicted) / max(predicted, 1)
+
+            results[name] = {
+                'pattern': pat,
+                'count': count,
+                'predicted': float(predicted),
+                'hardy_littlewood_constant': float(hl_const),
+                'relative_error': float(error),
+                'k': k,
+            }
+
+        self.results = {
+            'limit': limit,
+            'n_primes': n_primes,
+            'patterns': results,
+        }
+        return self.results
+
+    def summary(self) -> str:
+        r = self.results
+        s = f"Hardy-Littlewood k-Tuple Analysis (up to {r['limit']:,})\n{'='*50}\n"
+        for name, pat in r['patterns'].items():
+            s += f"  {name:12s}: count={pat['count']:,}, predicted={pat['predicted']:.1f}, "
+            s += f"𝔖={pat['hardy_littlewood_constant']:.6f}, rel_err={pat['relative_error']:.4f}\n"
+        return s