File size: 4,725 Bytes
fa128fd | 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 | -- True Entropy — Haskell verification layer.
-- H(P) = log2(N) - (1/N) * sum(c_i * log2(c_i))
-- Symbolic representation: exact rational coefficients * log2(integer base).
-- No Double ever in the core computation.
{-# LANGUAGE DataKinds #-}
{-# LANGUAGE ScopedTypeVariables #-}
module Verified.Entropy where
import Data.Ratio ((%))
import qualified Data.Map.Strict as Map
import Data.List (foldl', toList)
-- ─── Symbolic entropy ────────────────────────────────────────────────────────
-- | coeff * log2(base)
data SymLog2 = SymLog2
{ coeff :: Rational
, base :: Integer
} deriving (Eq, Show)
-- | Symbolic entropy: sum of SymLog2 terms
-- H = log2(N) + sum(-c_i/N * log2(c_i))
newtype SymbolicEntropy = SymbolicEntropy { terms :: [SymLog2] }
deriving Show
-- | Exact symbolic entropy from integer counts.
-- No evaluation, no rounding, no Double.
entropy :: (Foldable f, Integral a) => f a -> SymbolicEntropy
entropy counts = SymbolicEntropy (positiveTerm : negativeTerms)
where
cs = map fromIntegral $ filter (> 0) $ toList counts
n = sum cs :: Integer
-- H = log2(N) term
positiveTerm = SymLog2 { coeff = 1 % n, base = n }
-- -c_i/N * log2(c_i) terms
negativeTerms = map (\c -> SymLog2 { coeff = -(c % n), base = c }) cs
-- ─── Interval evaluation ─────────────────────────────────────────────────────
-- | Rational approximation of log2 for interval bounds.
-- Uses convergents of continued fraction for log2(n).
-- Returns (lower, upper) as Rational pair at given precision.
log2Interval :: Integer -> Int -> (Rational, Rational)
log2Interval n prec
| n <= 0 = error "log2Interval: non-positive argument"
| n == 1 = (0, 0)
| otherwise =
-- log2(n) = log(n) / log(2)
-- Use rational approximation: ln(2) ≈ 6931471805599453/10000000000000000
let ln2_lo = 6931471805599453 % 10000000000000000
ln2_hi = 6931471805599454 % 10000000000000000
-- ln(n) approximated via Taylor series for small n, else recursion
lnn = lnRational n prec
lo = fst lnn / ln2_hi -- divide by larger denominator → smaller result
hi = snd lnn / ln2_lo
in (lo, hi)
-- | Rational interval for ln(n), precision as number of terms
lnRational :: Integer -> Int -> (Rational, Rational)
lnRational n prec
| n == 1 = (0, 0)
| n == 2 = (6931471805599453 % 10000000000000000,
6931471805599454 % 10000000000000000)
| even n = let (l, h) = lnRational (n `div` 2) prec
(l2, h2) = lnRational 2 prec
in (l + l2, h + h2)
| otherwise = -- ln(n) ≈ ln(n-1) + 2/(2n-1) + ... (first term bound)
let (l, h) = lnRational (n - 1) prec
delta_lo = 2 % (2 * n - 1 + 1)
delta_hi = 2 % (2 * n - 1 - 0)
in (l + delta_lo, h + delta_hi)
-- ─── Connection to Bifrost WORM ──────────────────────────────────────────────
-- | Serialize entropy for WORM audit log
serializeEntropy :: SymbolicEntropy -> String
serializeEntropy (SymbolicEntropy ts) =
"H = " ++ unwords (map showTerm ts)
where
showTerm (SymLog2 c b) =
"(" ++ show (numerator c) ++ "/" ++ show (denominator c) ++
")*log2(" ++ show b ++ ")"
numerator r = let (n, _) = (floor (r * 10^15), ()) in n
denominator _ = 10^15 :: Integer
-- ─── Ramanujan connection ─────────────────────────────────────────────────────
-- | Partition entropy: entropy of the partition number sequence p(0)..p(n)
-- This is the entropy of the complexity invariant distribution.
partitionEntropy :: Int -> SymbolicEntropy
partitionEntropy n = entropy (partitionNumbers n)
-- | Euler pentagonal partition numbers (exact)
partitionNumbers :: Int -> [Integer]
partitionNumbers n = take (n + 1) ps
where
ps = 1 : [compute k | k <- [1..n]]
compute i = sum
[ sign k * safeIdx (i - penta k)
| k <- [1..i]
, penta k <= i
]
+
sum
[ sign (-k) * safeIdx (i - penta (-k))
| k <- [1..i]
, penta (-k) <= i
]
penta k = k * (3 * k - 1) `div` 2
sign k = if odd k then 1 else -1
safeIdx j = if j < 0 then 0 else ps !! j
|