Composite Numbers: Factorization, Primality Testing, and Cryptographic Implications

import math
from typing import List, Tuple, Dict, Optional
from dataclasses import dataclass


@dataclass
class Factorization:
    """Result of prime factorization."""
    number: int
    is_composite: bool
    prime_factors: Dict[int, int]  # prime -> exponent
    smallest_factor: Optional[int]
    divisor_count: int


class CompositeDetector:
    """Detect composite numbers and perform prime factorization."""

    def is_composite(self, n: int) -> bool:
        """
        Check if n is composite using trial division.
        Time: O(sqrt(n))
        Space: O(1)
        """
        if n <= 1:
            return False  # 1 is neither prime nor composite
        if n <= 3:
            return False  # 2 and 3 are prime
        if n % 2 == 0:
            return True   # all even numbers > 2 are composite
        if n % 3 == 0:
            return True   # all multiples of 3 > 3 are composite

        # Check divisors of form 6k +/- 1 up to sqrt(n)
        # All primes > 3 are of this form — skips 2/3 of candidates
        i = 5
        while i * i <= n:
            if n % i == 0 or n % (i + 2) == 0:
                return True
            i += 6

        return False

    def factorize(self, n: int) -> Factorization:
        """
        Complete prime factorization using trial division.
        Returns all prime factors and their exponents.
        """
        if n <= 1:
            return Factorization(n, False, {}, None, 1)

        factors: Dict[int, int] = {}
        temp = n

        # Factor out all 2s
        while temp % 2 == 0:
            factors[2] = factors.get(2, 0) + 1
            temp //= 2

        # Factor out odd primes from 3 onwards
        d = 3
        while d * d <= temp:
            while temp % d == 0:
                factors[d] = factors.get(d, 0) + 1
                temp //= d
            d += 2

        # If temp > 1, it is a prime factor
        if temp > 1:
            factors[temp] = 1

        is_composite = len(factors) > 1 or (len(factors) == 1 and list(factors.values())[0] > 1)
        smallest = min(factors.keys()) if factors else None

        # Calculate divisor count: if n = p1^a1 * p2^a2 * ...
        # then divisor_count = (a1+1) * (a2+1) * ...
        divisor_count = 1
        for exp in factors.values():
            divisor_count *= (exp + 1)

        return Factorization(
            number=n,
            is_composite=is_composite,
            prime_factors=factors,
            smallest_factor=smallest,
            divisor_count=divisor_count,
        )

    def get_divisors(self, n: int) -> List[int]:
        """Return all positive divisors of n."""
        if n <= 0:
            return []
        divisors = []
        for i in range(1, int(math.isqrt(n)) + 1):
            if n % i == 0:
                divisors.append(i)
                if i != n // i:
                    divisors.append(n // i)
        return sorted(divisors)

    def sieve_of_eratosthenes(self, limit: int) -> Tuple[List[bool], List[int]]:
        """
        Generate all primes up to limit using Sieve of Eratosthenes.
        is_composite[i] = True means i is composite.
        Time: O(n log log n)
        Space: O(n)
        """
        is_composite = [False] * (limit + 1)
        is_composite[0] = True  # 0 is not prime
        is_composite[1] = True  # 1 is not prime

        for i in range(2, int(math.isqrt(limit)) + 1):
            if not is_composite[i]:
                # Mark all multiples of i as composite
                for j in range(i * i, limit + 1, i):
                    is_composite[j] = True

        primes = [i for i in range(2, limit + 1) if not is_composite[i]]
        return is_composite, primes

    def segmented_sieve(self, low: int, high: int, segment_size: int = 65536) -> List[int]:
        """
        Find all primes in range [low, high] using segmented sieve.
        Memory: O(sqrt(high) + segment_size) instead of O(high).
        Use when range is large but you only need primes in a specific window.
        """
        sqrt_high = int(math.isqrt(high)) + 1
        _, base_primes = self.sieve_of_eratosthenes(sqrt_high)

        primes_in_range = []

        for seg_start in range(low, high + 1, segment_size):
            seg_end = min(seg_start + segment_size - 1, high)
            seg_size = seg_end - seg_start + 1
            is_composite_seg = [False] * seg_size

            for p in base_primes:
                # Find first multiple of p in [seg_start, seg_end]
                start = max(p * p, ((seg_start + p - 1) // p) * p)
                for j in range(start, seg_end + 1, p):
                    is_composite_seg[j - seg_start] = True

            for i in range(seg_size):
                num = seg_start + i
                if num >= 2 and not is_composite_seg[i]:
                    primes_in_range.append(num)

        return primes_in_range

    def count_composites_up_to(self, n: int) -> int:
        """
        Count composite numbers from 4 to n inclusive.
        Uses prime counting: composites = (n - 1) - prime_count
        (subtract 1 for the number 1 which is neither)
        """
        if n < 4:
            return 0
        _, primes = self.sieve_of_eratosthenes(n)
        return (n - 1) - len(primes)

import random
from typing import List


class PrimalityTester:
    """Production-grade primality testing with Miller-Rabin and deterministic witnesses."""

    # Deterministic witnesses for n < 2^64
    DETERMINISTIC_WITNESSES_64BIT = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]

    # Small primes for quick divisibility check
    SMALL_PRIMES = [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]

    def is_prime_trial_division(self, n: int) -> bool:
        """Deterministic primality test via trial division. O(sqrt(n))."""
        if n < 2:
            return False
        if n in (2, 3):
            return True
        if n % 2 == 0 or n % 3 == 0:
            return False
        i = 5
        while i * i <= n:
            if n % i == 0 or n % (i + 2) == 0:
                return False
            i += 6
        return True

    def _miller_rabin_round(self, n: int, witness: int) -> bool:
        """
        Single round of Miller-Rabin test.
        Returns True if n passes this round (possibly prime).
        Returns False if n is definitely composite.
        """
        # Write n-1 as 2^r * d where d is odd
        r, d = 0, n - 1
        while d % 2 == 0:
            r += 1
            d //= 2

        # Compute witness^d mod n
        x = pow(witness, d, n)
        if x == 1 or x == n - 1:
            return True  # passes this round

        for _ in range(r - 1):
            x = pow(x, 2, n)
            if x == n - 1:
                return True  # passes this round

        return False  # definitely composite

    def is_prime_miller_rabin(self, n: int, num_rounds: int = 20) -> bool:
        """
        Probabilistic primality test using Miller-Rabin.
        Error probability <= 4^(-num_rounds).
        For n < 2^64 with 12 fixed witnesses: deterministic.
        """
        if n < 2:
            return False
        if n in (2, 3):
            return True
        if n % 2 == 0:
            return False

        # Quick check against small primes
        for p in self.SMALL_PRIMES:
            if n == p:
                return True
            if n % p == 0:
                return False

        # For n < 2^64, use deterministic witnesses
        if n < (1 << 64):
            witnesses = self.DETERMINISTIC_WITNESSES_64BIT
        else:
            # Random witnesses for numbers beyond 64-bit
            witnesses = [random.randint(2, n - 2) for _ in range(num_rounds)]

        for w in witnesses:
            if not self._miller_rabin_round(n, w):
                return False  # definitely composite

        return True  # probably prime (or definitely prime for n < 2^64)

    def next_prime(self, n: int) -> int:
        """Find the smallest prime >= n."""
        if n <= 2:
            return 2
        candidate = n if n % 2 == 1 else n + 1
        while not self.is_prime_miller_rabin(candidate):
            candidate += 2
        return candidate

    def generate_rsa_primes(self, bit_length: int) -> tuple:
        """
        Generate two primes of specified bit length for RSA key generation.
        Uses Miller-Rabin with 40 rounds for cryptographic assurance.
        """
        while True:
            p = random.getrandbits(bit_length) | (1 << (bit_length - 1)) | 1
            if self.is_prime_miller_rabin(p, num_rounds=40):
                break

        while True:
            q = random.getrandbits(bit_length) | (1 << (bit_length - 1)) | 1
            if q != p and self.is_prime_miller_rabin(q, num_rounds=40):
                break

        return p, q

    def verify_rsa_modulus_uniqueness(self, moduli: List[int]) -> List[tuple]:
        """
        Check for shared prime factors across RSA moduli.
        Returns list of (index_i, index_j, shared_factor) for any pair with GCD > 1.
        """
        from math import gcd
        collisions = []
        for i in range(len(moduli)):
            for j in range(i + 1, len(moduli)):
                g = gcd(moduli[i], moduli[j])
                if g > 1:
                    collisions.append((i, j, g))
        return collisions

import random
from math import gcd
from typing import List, Tuple, Optional


class RSACompositeAnalyzer:
    """Analyze RSA moduli for composite number vulnerabilities."""

    def generate_rsa_keypair(self, bit_length: int = 1024) -> dict:
        """
        Generate RSA keypair with post-generation verification.
        Ensures p != q and both are verified prime.
        """
        from io.thecodeforge.math.primality_tester import PrimalityTester
        tester = PrimalityTester()

        p, q = tester.generate_rsa_primes(bit_length)
        n = p * q
        phi_n = (p - 1) * (q - 1)
        e = 65537  # standard public exponent

        # Compute modular inverse: d = e^(-1) mod phi(n)
        d = pow(e, -1, phi_n)

        return {
            'public_key': (n, e),
            'private_key': (n, d),
            'modulus': n,
            'modulus_bit_length': n.bit_length(),
            'p': p,
            'q': q,
        }

    def detect_shared_factors(self, moduli: List[int]) -> List[dict]:
        """
        Detect RSA moduli that share prime factors.
        GCD(n1, n2) > 1 means both keys are compromised.
        """
        findings = []
        for i in range(len(moduli)):
            for j in range(i + 1, len(moduli)):
                g = gcd(moduli[i], moduli[j])
                if g > 1 and g != moduli[i] and g != moduli[j]:
                    p = g
                    q1 = moduli[i] // g
                    q2 = moduli[j] // g
                    findings.append({
                        'key_pair': (i, j),
                        'shared_factor': g,
                        'key_i_factors': (g, q1),
                        'key_j_factors': (g, q2),
                        'severity': 'CRITICAL',
                        'action': 'Revoke both keys immediately. Re-generate with verified entropy.',
                    })
        return findings

    def factor_via_shared_prime(self, n1: int, n2: int) -> Optional[Tuple[Tuple[int, int], Tuple[int, int]]]:
        """
        If two RSA moduli share a prime factor, factor both instantly.
        This is why shared primes are catastrophic.
        """
        g = gcd(n1, n2)
        if g > 1 and g != n1 and g != n2:
            p = g
            q1 = n1 // g
            q2 = n2 // g
            return ((p, q1), (p, q2))
        return None

    def check_smooth_number(self, n: int, smoothness_bound: int = 1000) -> dict:
        """
        Check if n is B-smooth (all prime factors <= B).
        Smooth numbers are easier to factor using algorithms like Pollard's p-1.
        """
        from io.thecodeforge.math.composite_detector import CompositeDetector
        detector = CompositeDetector()

        temp = n
        factors = {}

        # Simplified prime check; in production use a proper prime list
        def is_prime_simple(x: int) -> bool:
            if x < 2:
                return False
            for i in range(2, int(x**0.5) + 1):
                if x % i == 0:
                    return False
            return True

        for p in range(2, smoothness_bound + 1):
            if is_prime_simple(p):
                while temp % p == 0:
                    factors[p] = factors.get(p, 0) + 1
                    temp //= p

        is_smooth = temp == 1

        return {
            'number': n,
            'smoothness_bound': smoothness_bound,
            'is_smooth': is_smooth,
            'factors_found': factors,
            'remaining_cofactor': temp if temp > 1 else None,
            'risk': 'HIGH' if is_smooth else 'LOW',
            'note': (
                'B-smooth numbers are vulnerable to Pollard p-1 factorization. '
                'RSA primes should be generated to avoid smoothness.'
                if is_smooth else
                'Number is not B-smooth at the given bound. Standard factoring difficulty applies.'
            ),
        }

    def estimate_factoring_difficulty(self, n: int) -> dict:
        """
        Estimate the difficulty of factoring a composite number
        based on its bit length and structure.
        """
        bit_length = n.bit_length()

        if bit_length <= 64:
            difficulty = 'TRIVIAL'
            method = 'Trial division or Pollard Rho. Seconds on modern hardware.'
            time_estimate = '< 1 second'
        elif bit_length <= 128:
            difficulty = 'EASY'
            method = 'Pollard Rho or ECM. Minutes on modern hardware.'
            time_estimate = '1-10 minutes'
        elif bit_length <= 256:
            difficulty = 'MODERATE'
            method = 'ECM or Quadratic Sieve. Hours to days.'
            time_estimate = '1-48 hours'
        elif bit_length <= 512:
            difficulty = 'HARD'
            method = 'Quadratic Sieve or GNFS. Days to weeks on a cluster.'
            time_estimate = '1-30 days'
        elif bit_length <= 1024:
            difficulty = 'VERY HARD'
            method = 'GNFS. Months to years on a large cluster.'
            time_estimate = '1-5 years'
        elif bit_length <= 2048:
            difficulty = 'EXTREME'
            method = 'GNFS. Requires massive computational resources.'
            time_estimate = '10-100 years on current hardware'
        else:
            difficulty = 'IMPOSSIBLE (classical)'
            method = 'GNFS'
            time_estimate = '> age of universe'

        return {
            'bit_length': bit_length,
            'difficulty': difficulty,
            'method': method,
            'time_estimate': time_estimate,
            'quantum_note': (
                'Shor algorithm on a sufficiently large quantum computer '
                'could factor this in polynomial time. Post-quantum migration recommended.'
                if bit_length <= 2048 else
                'Quantum threat applies to all RSA key sizes.'
            ),
        }

import math
from typing import List, Iterator
import array


class OptimizedSieve:
    """Production-grade sieve implementations for prime generation and composite detection."""

    def basic_sieve(self, limit: int) -> List[int]:
        """
        Standard Sieve of Eratosthenes.
        Time: O(n log log n), Space: O(n)
        """
        if limit < 2:
            return []

        is_composite = bytearray(limit + 1)

        for i in range(2, int(math.isqrt(limit)) + 1):
            if not is_composite[i]:
                for j in range(i * i, limit + 1, i):
                    is_composite[j] = 1

        return [i for i in range(2, limit + 1) if not is_composite[i]]

    def bit_packed_sieve(self, limit: int) -> List[int]:
        """
        Bit-packed sieve — 8x memory reduction over boolean array.
        Each number uses 1 bit instead of 1 byte.
        For limit=10^9: 125MB instead of 1GB.
        """
        if limit < 2:
            return []

        # Each byte stores 8 flags (odd numbers only)
        size = (limit + 1) // 2
        sieve = bytearray((size + 7) // 8)

        def is_set(idx: int) -> bool:
            return bool(sieve[idx // 8] & (1 << (idx % 8)))

        def clear_bit(idx: int):
            sieve[idx // 8] |= (1 << (idx % 8))

        # Only sieve odd numbers: index i represents number 2*i + 1
        sqrt_limit = int(math.isqrt(limit))
        for i in range(3, sqrt_limit + 1, 2):
            if not is_set(i // 2):
                # Mark odd multiples of i starting from i*i
                step = i  # step in terms of odd numbers
                start = (i * i) // 2
                for j in range(start, size, step):
                    clear_bit(j)

        # Collect primes
        primes = [2]  # 2 is the only even prime
        for i in range(1, size):
            if not is_set(i):
                primes.append(2 * i + 1)

        return [p for p in primes if p <= limit]

    def segmented_sieve(self, low: int, high: int, segment_size: int = 65536) -> List[int]:
        """
        Segmented sieve for finding primes in [low, high].
        Memory: O(sqrt(high) + segment_size)
        Use when high > 10^7 and you cannot allocate O(high) memory.
        """
        sqrt_high = int(math.isqrt(high)) + 1

        # Generate base primes up to sqrt(high)
        base_is_composite = bytearray(sqrt_high + 1)
        for i in range(2, int(math.isqrt(sqrt_high)) + 1):
            if not base_is_composite[i]:
                for j in range(i * i, sqrt_high + 1, i):
                    base_is_composite[j] = 1
        base_primes = [i for i in range(2, sqrt_high + 1) if not base_is_composite[i]]

        primes = []

        for seg_start in range(max(2, low), high + 1, segment_size):
            seg_end = min(seg_start + segment_size - 1, high)
            seg_is_composite = bytearray(seg_end - seg_start + 1)

            for p in base_primes:
                # First multiple of p in [seg_start, seg_end]
                start = max(p * p, ((seg_start + p - 1) // p) * p)
                for j in range(start, seg_end + 1, p):
                    seg_is_composite[j - seg_start] = 1

            for i in range(len(seg_is_composite)):
                num = seg_start + i
                if num >= 2 and not seg_is_composite[i]:
                    primes.append(num)

        return primes

    def wheel_sieve(self, limit: int) -> List[int]:
        """
        Wheel factorization sieve — skips multiples of 2, 3, 5.
        Checks only numbers of form 30k + {1,7,11,13,17,19,23,29}.
        3.75x fewer candidates than basic sieve.
        """
        if limit < 2:
            return []

        # Base primes up to the wheel size
        small_primes = [2, 3, 5]
        if limit < 7:
            return [p for p in small_primes if p <= limit]

        # Wheel offsets: numbers mod 30 that are coprime to 2, 3, 5
        wheel = [1, 7, 11, 13, 17, 19, 23, 29]
        wheel_gaps = [6, 4, 2, 4, 2, 4, 6, 2]  # gaps between wheel positions

        # Sieve array for numbers in the wheel pattern
        size = (limit + 29) // 30
        sieve = bytearray(size * 8)  # 8 flags per 30-number block

        sqrt_limit = int(math.isqrt(limit))

        # Sieve using wheel offsets
        for i in range(1, size):
            for j, offset in enumerate(wheel):
                num = 30 * i + offset
                if num > sqrt_limit:
                    break
                if not sieve[i * 8 + j]:
                    # Mark multiples of num
                    for k in range(i, size):
                        for l, inner_offset in enumerate(wheel):
                            multiple = 30 * k + inner_offset
                            if multiple > limit:
                                break
                            if multiple % num == 0:
                                sieve[k * 8 + l] = 1

        # Collect primes
        primes = list(small_primes)
        for i in range(size):
            for j, offset in enumerate(wheel):
                num = 30 * i + offset
                if num > limit:
                    break
                if num >= 7 and not sieve[i * 8 + j]:
                    primes.append(num)

        return primes

    def prime_iterator(self, limit: int) -> Iterator[int]:
        """
        Memory-efficient prime iterator using segmented sieve.
        Yields primes one at a time without storing all in memory.
        """
        segment_size = 65536
        for seg_start in range(2, limit + 1, segment_size):
            seg_end = min(seg_start + segment_size - 1, limit)
            for p in self.segmented_sieve(seg_start, seg_end, segment_size):
                yield p

import math
from typing import List, Dict, Tuple


class CompositeApplications:
    """Composite number applications in hash tables, PRNGs, and modular arithmetic."""

    def analyze_hash_distribution(self, keys: List[int], bucket_count: int) -> dict:
        """
        Analyze key distribution across buckets to detect composite-related clustering.
        """
        buckets = [0] * bucket_count
        for key in keys:
            buckets[key % bucket_count] += 1

        non_empty = sum(1 for b in buckets if b > 0)
        max_bucket = max(buckets)
        min_bucket = min(b for b in buckets if b > 0) if non_empty > 0 else 0
        avg_bucket = len(keys) / non_empty if non_empty > 0 else 0

        # Calculate variance from ideal distribution
        ideal = len(keys) / bucket_count
        variance = sum((b - ideal) ** 2 for b in buckets) / bucket_count

        return {
            'bucket_count': bucket_count,
            'is_prime_bucket_count': self._is_prime(bucket_count),
            'total_keys': len(keys),
            'non_empty_buckets': non_empty,
            'utilization': round(non_empty / bucket_count * 100, 1),
            'max_bucket_size': max_bucket,
            'min_bucket_size': min_bucket,
            'avg_bucket_size': round(avg_bucket, 2),
            'variance': round(variance, 2),
            'recommendation': (
                'Prime bucket count — good distribution expected.'
                if self._is_prime(bucket_count) else
                'Composite bucket count — consider switching to nearest prime for better distribution.'
            ),
        }

    def _is_prime(self, n: int) -> bool:
        """Quick primality check for small numbers."""
        if n < 2:
            return False
        if n < 4:
            return True
        if n % 2 == 0 or n % 3 == 0:
            return False
        i = 5
        while i * i <= n:
            if n % i == 0 or n % (i + 2) == 0:
                return False
            i += 6
        return True

    def euler_totient(self, n: int) -> int:
        """
        Compute Euler's totient function phi(n).
        phi(n) counts integers 1..n that are coprime to n.
        For prime p: phi(p) = p - 1
        For composite n = p1^a1 * p2^a2 * ...: phi(n) = n * prod(1 - 1/pi)
        """
        result = n
        temp = n
        p = 2
        while p * p <= temp:
            if temp % p == 0:
                while temp % p == 0:
                    temp //= p
                result -= result // p  # result *= (1 - 1/p)
            p += 1
        if temp > 1:
            result -= result // temp
        return result

    def find_carmichael_numbers(self, limit: int) -> List[int]:
        """
        Find Carmichael numbers up to limit.
        Carmichael numbers are composite numbers n where a^(n-1) = 1 mod n
        for all a coprime to n. They fool Fermat's primality test.
        """
        carmichael = []
        for n in range(2, limit + 1):
            if self._is_prime(n):
                continue  # Carmichael numbers must be composite

            # Must be square-free (no repeated prime factors)
            temp = n
            is_square_free = True
            p = 2
            while p * p <= temp:
                count = 0
                while temp % p == 0:
                    count += 1
                    temp //= p
                if count > 1:
                    is_square_free = False
                    break
                p += 1
            if not is_square_free:
                continue

            # Must have at least 3 prime factors
            phi = self.euler_totient(n)
            if (n - 1) % phi == 0:
                # Verify Fermat condition for a few bases
                passes_all = True
                for a in [2, 3, 5, 7, 11]:
                    if math.gcd(a, n) == 1:
                        if pow(a, n - 1, n) != 1:
                            passes_all = False
                            break
                if passes_all:
                    carmichael.append(n)

        return carmichael

    def lcg_period_analysis(self, a: int, c: int, m: int, num_samples: int = 100000) -> dict:
        """
        Analyze the period of a Linear Congruential Generator.
        X_{n+1} = (a * X_n + c) mod m
        Maximum period = m if: gcd(c,m)=1, and for every prime p|m: (a-1)=0 mod p,
        and if 4|m: (a-1)=0 mod 4.
        """
        x = 1  # seed
        seen = {}
        period = 0

        for i in range(num_samples):
            if x in seen:
                period = i - seen[x]
                break
            seen[x] = i
            x = (a * x + c) % m

        max_possible_period = m
        achieves_max = period == max_possible_period

        # Check conditions for maximum period
        conditions = {
            'gcd_c_m_is_1': math.gcd(c, m) == 1,
        }

        # Check: for every prime p dividing m, (a-1) = 0 mod p
        temp = m
        p = 2
        prime_conditions = True
        while p * p <= temp:
            if temp % p == 0:
                if (a - 1) % p != 0:
                    prime_conditions = False
                while temp % p == 0:
                    temp //= p
            p += 1
        if temp > 1 and (a - 1) % temp != 0:
            prime_conditions = False
        conditions['a_minus_1_divisible_by_all_prime_factors_of_m'] = prime_conditions

        # If 4 divides m, check (a-1) = 0 mod 4
        if m % 4 == 0:
            conditions['a_minus_1_divisible_by_4'] = (a - 1) % 4 == 0

        return {
            'parameters': {'a': a, 'c': c, 'm': m},
            'm_is_prime': self._is_prime(m),
            'm_is_composite': not self._is_prime(m) and m > 1,
            'observed_period': period,
            'max_possible_period': max_possible_period,
            'achieves_max_period': achieves_max,
            'conditions': conditions,
        }