Optimized TNFR Primality Testing Implementation
This module provides an enhanced version of the TNFR primality test that leverages the full infrastructure of the TNFR repository:
Author: TNFR Research Team Date: 2025-11-29 Status: OPTIMIZED IMPLEMENTATION
"""
Optimized TNFR Primality Testing Implementation
This module provides an enhanced version of the TNFR primality test
that leverages the full infrastructure of the TNFR repository:
- Centralized caching system
- Vectorized operations
- Structural field computations
- GPU backends when available
- Mathematical optimization techniques
Author: TNFR Research Team
Date: 2025-11-29
Status: OPTIMIZED IMPLEMENTATION
"""
from __future__ import annotations
import logging
import math
import time
from dataclasses import dataclass
from typing import Any
from ..backends.optimized_numpy import OptimizedNumpyBackend
from ..metrics.coherence import compute_coherence
from ..physics.fields import compute_structural_potential
# Core TNFR infrastructure
from .unified_cache import CacheLevel, cache_tnfr_computation
from .unified_numerical import np
# GPU acceleration if available
try:
from ..backends.torch_backend import TorchBackend
HAS_TORCH = True
except ImportError:
HAS_TORCH = False
try:
from ..backends.jax_backend import JAXBackend
HAS_JAX = True
except ImportError:
HAS_JAX = False
# Mathematical libraries
try:
HAS_NUMBA = True
except ImportError:
HAS_NUMBA = False
logger = logging.getLogger(__name__)
# TNFR primality parameters: ΔNFR pressure coefficients.
# Prime ⟺ ΔNFR = 0 with canonical unit-scale coefficients (TNFR_NUMBER_THEORY).
ZETA_CANONICAL = 1.0 # Factorization pressure coefficient
ETA_CANONICAL = 0.8 # Divisor pressure coefficient
THETA_CANONICAL = 0.6 # Sigma pressure coefficient
# Primality decision cut on |ΔNFR|: primes give ΔNFR = 0 exactly while
# composites give |ΔNFR| ≳ 2 with the coefficients above, so any threshold in
# (0, 2) is robust. 0.5 is the plain canonical separator.
PRIME_THRESHOLD_HP = 0.5
@dataclass
class PrimalityResult:
"""Enhanced result structure for optimized primality testing."""
n: int
is_prime: bool
delta_nfr: float
computation_time_ms: float
method: str
confidence: float = 1.0
structural_metrics: dict[str, float] | None = None
cache_hit: bool = False
class OptimizedTNFRPrimality:
"""
Optimized TNFR primality testing with multiple acceleration strategies.
Features:
- Multi-tier caching (LRU + persistent)
- Vectorized batch operations
- GPU acceleration when available
- Sieve-based preprocessing
- Structural field integration
"""
def __init__(
self,
*,
backend: str = "auto",
cache_size: int = 10000,
sieve_limit: int = 1000000,
enable_gpu: bool = True,
precision_mode: str = "standard",
):
self.backend_name = backend
self.cache_size = cache_size
self.sieve_limit = sieve_limit
self.enable_gpu = enable_gpu
self.precision_mode = precision_mode
# Initialize backend
self.backend = self._init_backend()
# Precompute sieve for fast factorization
self.sieve_data = self._build_sieve(sieve_limit)
# Initialize caches
self._init_caches()
logger.info(
f"OptimizedTNFRPrimality initialized: backend={self.backend_name}, "
f"cache_size={cache_size}, sieve_limit={sieve_limit}"
)
def _init_backend(self):
"""Initialize the computational backend."""
if self.backend_name == "auto":
if HAS_JAX and self.enable_gpu:
return JAXBackend()
elif HAS_TORCH and self.enable_gpu:
return TorchBackend()
else:
return OptimizedNumpyBackend()
elif self.backend_name == "jax" and HAS_JAX:
return JAXBackend()
elif self.backend_name == "torch" and HAS_TORCH:
return TorchBackend()
elif self.backend_name == "numpy":
return OptimizedNumpyBackend()
else:
logger.warning(f"Backend {self.backend_name} not available, using numpy")
return OptimizedNumpyBackend()
def _build_sieve(self, limit: int) -> dict[str, np.ndarray]:
"""Build optimized sieve for fast factorization."""
start_time = time.perf_counter()
# Sieve of Eratosthenes for primes
is_prime = np.ones(limit + 1, dtype=bool)
is_prime[0] = is_prime[1] = False
for i in range(2, int(math.sqrt(limit)) + 1):
if is_prime[i]:
is_prime[i * i : limit + 1 : i] = False
primes = np.where(is_prime)[0]
# Minimum prime factor for each number
min_prime_factor = np.arange(limit + 1, dtype=np.int32)
for p in primes:
if p * p <= limit:
for i in range(p * p, limit + 1, p):
if min_prime_factor[i] == i: # First time we see this number
min_prime_factor[i] = p
elapsed = (time.perf_counter() - start_time) * 1000
logger.info(
f"Sieve built in {elapsed:.2f}ms: {len(primes)} primes up to {limit}"
)
return {
"is_prime": is_prime,
"primes": primes,
"min_prime_factor": min_prime_factor,
"limit": limit,
}
def _init_caches(self):
"""Initialize multi-tier caching system."""
self.result_cache = {} # Simple dict cache for results
self.arithmetic_cache = {} # Cache for arithmetic functions
# Statistics
self.cache_hits = 0
self.cache_misses = 0
@cache_tnfr_computation(level=CacheLevel.DERIVED_METRICS)
def _fast_divisor_count(self, n: int) -> int:
"""Optimized divisor count using sieve when possible."""
if n <= self.sieve_data["limit"] and n >= 1:
return self._sieve_divisor_count(n)
return self._trial_divisor_count(n)
def _sieve_divisor_count(self, n: int) -> int:
"""Ultra-fast divisor count using precomputed sieve."""
if n == 1:
return 1
count = 1 # Start with divisor 1
temp = n
min_pf = self.sieve_data["min_prime_factor"]
while temp > 1:
p = min_pf[temp]
exp = 0
while temp % p == 0:
exp += 1
temp //= p
count *= exp + 1
return count
def _trial_divisor_count(self, n: int) -> int:
"""Fallback divisor count for large numbers."""
count = 0
i = 1
sqrt_n = int(math.sqrt(n))
while i <= sqrt_n:
if n % i == 0:
count += 1
if i != n // i:
count += 1
i += 1
return count
@cache_tnfr_computation(level=CacheLevel.DERIVED_METRICS)
def _fast_divisor_sum(self, n: int) -> int:
"""Optimized divisor sum using sieve when possible."""
if n <= self.sieve_data["limit"] and n >= 1:
return self._sieve_divisor_sum(n)
return self._trial_divisor_sum(n)
def _sieve_divisor_sum(self, n: int) -> int:
"""Ultra-fast divisor sum using precomputed sieve."""
if n == 1:
return 1
total = 1 # Start with divisor 1
temp = n
min_pf = self.sieve_data["min_prime_factor"]
while temp > 1:
p = min_pf[temp]
exp = 0
p_power = 1
while temp % p == 0:
exp += 1
p_power *= p
temp //= p
# Sum of geometric series: (p^(exp+1) - 1) / (p - 1)
total *= (p_power * p - 1) // (p - 1)
return total
def _trial_divisor_sum(self, n: int) -> int:
"""Fallback divisor sum for large numbers."""
total = 0
i = 1
sqrt_n = int(math.sqrt(n))
while i <= sqrt_n:
if n % i == 0:
total += i
if i != n // i:
total += n // i
i += 1
return total
@cache_tnfr_computation(level=CacheLevel.DERIVED_METRICS)
def _fast_omega(self, n: int) -> int:
"""Optimized prime factor count (ω function)."""
if n <= self.sieve_data["limit"] and n >= 1:
return self._sieve_omega(n)
return self._trial_omega(n)
def _sieve_omega(self, n: int) -> int:
"""Ultra-fast ω(n) using precomputed sieve."""
if n <= 1:
return 0
count = 0
temp = n
min_pf = self.sieve_data["min_prime_factor"]
last_p = 0
while temp > 1:
p = min_pf[temp]
if p != last_p:
count += 1
last_p = p
temp //= p
return count
def _trial_omega(self, n: int) -> int:
"""Fallback ω(n) for large numbers."""
if n <= 1:
return 0
count = 0
d = 2
while d * d <= n:
if n % d == 0:
count += 1
while n % d == 0:
n //= d
d += 1
if n > 1:
count += 1
return count
def compute_delta_nfr(
self,
n: int,
*,
zeta: float = ZETA_CANONICAL,
eta: float = ETA_CANONICAL,
theta: float = THETA_CANONICAL,
) -> float:
"""
Optimized TNFR ΔNFR computation with caching and vectorization.
ΔNFR(n) = ζ·(ω(n)−1) + η·(τ(n)−2) + θ·(σ(n)/n − (1+1/n))
"""
if n < 2:
return float("inf")
# Check cache first
cache_key = (n, zeta, eta, theta)
if cache_key in self.arithmetic_cache:
self.cache_hits += 1
return self.arithmetic_cache[cache_key]
self.cache_misses += 1
# Compute arithmetic functions
tau_n = self._fast_divisor_count(n)
sigma_n = self._fast_divisor_sum(n)
omega_n = self._fast_omega(n)
# TNFR pressure components
factorization_pressure = zeta * (omega_n - 1)
divisor_pressure = eta * (tau_n - 2)
sigma_pressure = theta * (sigma_n / n - (1 + 1 / n))
delta_nfr = factorization_pressure + divisor_pressure + sigma_pressure
# Cache result
self.arithmetic_cache[cache_key] = delta_nfr
return delta_nfr
def is_prime_optimized(
self,
n: int,
*,
threshold: float = PRIME_THRESHOLD_HP,
include_metrics: bool = False,
) -> PrimalityResult:
"""
Optimized TNFR primality test with comprehensive result structure.
Args:
n: Integer to test for primality
threshold: ΔNFR threshold for primality (default: theoretical optimum)
include_metrics: Whether to compute structural field metrics
Returns:
PrimalityResult with comprehensive information
"""
start_time = time.perf_counter()
# Check result cache
cache_key = (n, threshold, include_metrics)
if cache_key in self.result_cache:
result = self.result_cache[cache_key]
result.cache_hit = True
return result
# Fast path for small numbers using sieve
if n <= self.sieve_data["limit"] and n >= 2:
is_prime_sieve = self.sieve_data["is_prime"][n]
delta_nfr = self.compute_delta_nfr(n)
method = "sieve_lookup"
else:
# TNFR computation for large numbers
delta_nfr = self.compute_delta_nfr(n)
is_prime_sieve = abs(delta_nfr) < threshold
method = "tnfr_computation"
# Compute structural metrics if requested
structural_metrics = None
if include_metrics and n >= 2:
structural_metrics = self._compute_structural_metrics(n)
# Calculate confidence based on ΔNFR distance from threshold
confidence = min(1.0, abs(delta_nfr - threshold) / threshold + 0.5)
elapsed_ms = (time.perf_counter() - start_time) * 1000
result = PrimalityResult(
n=n,
is_prime=bool(is_prime_sieve),
delta_nfr=float(delta_nfr),
computation_time_ms=elapsed_ms,
method=method,
confidence=confidence,
structural_metrics=structural_metrics,
cache_hit=False,
)
# Cache result
self.result_cache[cache_key] = result
return result
def _compute_structural_metrics(self, n: int) -> dict[str, float]:
"""Compute TNFR structural field metrics for the number."""
try:
# Create minimal graph for structural computations
import networkx as nx
G = nx.Graph()
G.add_node(n)
# Add some context nodes for field computation
for i in range(max(2, n - 2), min(n + 3, n + 10)):
if i != n:
G.add_node(i)
if abs(i - n) <= 2: # Connect nearby numbers
G.add_edge(n, i)
# set phases based on logarithmic scaling
for node in G.nodes():
G.nodes[node]["phase"] = math.log(node) if node > 1 else 0.0
G.nodes[node]["nu_f"] = 1.0 # Base frequency
# Compute structural potential
phi_s = compute_structural_potential(G)
# Compute coherence if graph has edges
coherence = 0.0
if G.number_of_edges() > 0:
coherence = compute_coherence(G)
return {
"structural_potential": float(phi_s.get(n, 0.0)),
"coherence": float(coherence),
"node_degree": G.degree(n),
"graph_size": G.number_of_nodes(),
}
except Exception as e:
logger.debug(f"Structural metrics computation failed for n={n}: {e}")
return {}
def batch_test(
self,
numbers: list[int],
*,
threshold: float = PRIME_THRESHOLD_HP,
include_metrics: bool = False,
) -> list[PrimalityResult]:
"""
Batch primality testing with vectorized optimizations.
Args:
numbers: list of integers to test
threshold: ΔNFR threshold for primality
include_metrics: Whether to compute structural metrics
Returns:
list of PrimalityResult objects
"""
results = []
# Sort numbers for better cache locality
sorted_numbers = sorted(set(numbers))
start_time = time.perf_counter()
for n in sorted_numbers:
result = self.is_prime_optimized(
n, threshold=threshold, include_metrics=include_metrics
)
results.append(result)
elapsed_ms = (time.perf_counter() - start_time) * 1000
# Log batch statistics
cache_hit_rate = sum(1 for r in results if r.cache_hit) / len(results)
avg_time = sum(r.computation_time_ms for r in results) / len(results)
logger.info(
f"Batch tested {len(results)} numbers in {elapsed_ms:.2f}ms "
f"(cache hit rate: {cache_hit_rate:.1%}, avg: {avg_time:.3f}ms/number)"
)
return results
def get_statistics(self) -> dict[str, Any]:
"""Get performance and cache statistics."""
total_requests = self.cache_hits + self.cache_misses
hit_rate = self.cache_hits / total_requests if total_requests > 0 else 0
return {
"backend": self.backend_name,
"sieve_limit": self.sieve_limit,
"cache_size": len(self.result_cache),
"arithmetic_cache_size": len(self.arithmetic_cache),
"cache_hits": self.cache_hits,
"cache_misses": self.cache_misses,
"hit_rate": hit_rate,
"primes_in_sieve": len(self.sieve_data["primes"]),
"sieve_coverage": self.sieve_data["limit"],
}
def clear_caches(self):
"""Clear all caches and reset statistics."""
self.result_cache.clear()
self.arithmetic_cache.clear()
self.cache_hits = 0
self.cache_misses = 0
logger.info("All caches cleared")
# Convenience functions for backward compatibility
def tnfr_is_prime_optimized(
n: int, *, threshold: float = PRIME_THRESHOLD_HP
) -> tuple[bool, float]:
"""
Optimized TNFR primality test (backward compatible interface).
Returns:
tuple of (is_prime, delta_nfr)
"""
# Global instance for stateless usage
global _global_optimizer
if "_global_optimizer" not in globals():
_global_optimizer = OptimizedTNFRPrimality()
result = _global_optimizer.is_prime_optimized(n, threshold=threshold)
return result.is_prime, result.delta_nfr
def benchmark_optimization(
max_n: int = 100000, sample_size: int = 1000
) -> dict[str, Any]:
"""
Benchmark the optimized implementation against baseline.
Args:
max_n: Maximum number to test
sample_size: Number of random samples to test
Returns:
Benchmark results dictionary
"""
import random
# Generate test numbers
test_numbers = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29] # Small primes
test_numbers.extend([4, 6, 8, 9, 10, 12, 14, 15, 16, 18]) # Small composites
# Add random samples
random.seed(42) # Reproducible
test_numbers.extend(random.sample(range(100, max_n), min(sample_size, max_n - 100)))
# Initialize optimizer
optimizer = OptimizedTNFRPrimality(sieve_limit=max_n)
# Benchmark
start_time = time.perf_counter()
results = optimizer.batch_test(test_numbers)
elapsed_time = time.perf_counter() - start_time
# Analyze results
prime_count = sum(1 for r in results if r.is_prime)
avg_time_per_number = elapsed_time * 1000 / len(results) # ms
max_time = max(r.computation_time_ms for r in results)
min_time = min(r.computation_time_ms for r in results)
stats = optimizer.get_statistics()
return {
"total_numbers_tested": len(results),
"primes_found": prime_count,
"prime_ratio": prime_count / len(results),
"total_time_ms": elapsed_time * 1000,
"avg_time_per_number_ms": avg_time_per_number,
"max_time_ms": max_time,
"min_time_ms": min_time,
"throughput_numbers_per_sec": len(results) / elapsed_time,
"cache_statistics": stats,
"largest_number_tested": max(test_numbers),
"backend_used": optimizer.backend_name,
}
if __name__ == "__main__":
# Quick demonstration
optimizer = OptimizedTNFRPrimality()
# Test some numbers
test_cases = [2, 3, 17, 97, 1009, 10007, 982451653]
print("Optimized TNFR Primality Testing Demo")
print("=" * 50)
for n in test_cases:
result = optimizer.is_prime_optimized(n, include_metrics=True)
print(
f"n={n:>10}: prime={result.is_prime}, "
f"ΔNFR={result.delta_nfr:8.6f}, "
f"time={result.computation_time_ms:6.3f}ms, "
f"method={result.method}"
)
print("\nPerformance Statistics:")
print("-" * 30)
stats = optimizer.get_statistics()
for key, value in stats.items():
print(f"{key}: {value}")