BioSight/python/phi/test_verified_units.py

import os
import json
import sys

# Add package directory to path so we can import from phi
sys.path.insert(0, os.path.dirname(os.path.dirname(os.path.abspath(__file__))))

from phi.charclass import classify_char
from phi.consistency import check_consistency

# Helper functions for Q16_16 math in Python matching Lean's FixedPoint behavior.
#
# V4 MITIGATION — Lean 4 Int.div uses EUCLIDEAN division:
#   - The remainder is always non-negative: a = q*b + r, r >= 0
#   - For POSITIVE divisor: ediv == Python's // (floor division)
#   - For NEGATIVE divisor: ediv != Python's //
#     Example: Lean:   65536 / (-3) = -21845  (rem 1)
#              Python: 65536 // (-3) = -21846  (rem -2, floor)
#
# Since q16_mul always divides by 65536 (positive), Python's // is safe there.
# q16_div divides by an arbitrary b which may be negative, so we use _ediv.

def _ediv(a, b):
    """Euclidean division matching Lean 4's Int.div.
    Guarantees: a = _ediv(a,b)*b + _emod(a,b) and _emod(a,b) >= 0."""
    if b == 0:
        return 0
    # Python's divmod gives floor division (q, r) with r same sign as b.
    # Euclidean requires r >= 0 always.
    q, r = divmod(a, b)
    if r < 0:
        # This happens when b < 0 (Python's divmod gives r with sign of b)
        q += 1
        # r -= b  (but we don't need r for the return)
    return q

def q16_mul(a, b):
    # Divisor is always 65536 (positive), so Python // matches Lean ediv.
    return (a * b) // 65536

def q16_div(a, b):
    if b == 0:
        return 2147483647  # Q16_16 infinity sentinel
    # b may be negative, so we must use Euclidean division.
    return _ediv(a * 65536, b)

def int_sqrt(n):
    """Integer square root via Newton's method (floor(√n)).
    Uses floor division intentionally here because we want floor(√n),
    not truncation — the Newton iteration converges correctly with
    floor division on non-negative integers."""
    if n < 0:
        raise ValueError("Square root of negative number")
    if n == 0:
        return 0
    x = n
    y = (x + 1) // 2  # floor div is correct here (n ≥ 0)
    while y < x:
        x = y
        y = (x + n // x) // 2  # floor div is correct here (both operands ≥ 0)
    return x

def q16_sqrt(q_raw):
    if q_raw <= 0:
        return 0
    return int_sqrt(q_raw * 65536)

def parabola_conjugate_pair(m_raw):
    m_sq = q16_mul(m_raw, m_raw)
    sqrt_term = q16_sqrt(m_sq + 65536)
    s1 = m_raw + sqrt_term
    s2 = q16_div(-65536, s1)
    return s1, s2

def run_verified_tests():
    fixtures_path = os.path.join(os.path.dirname(__file__), "fixtures.json")
    with open(fixtures_path, "r", encoding="utf-8") as f:
        fixtures = json.load(f)

    failures = 0

    print("=== Running Verified Test Units ===")

    # 1. Character Classification Verification
    print("Running Character Classification Tests...")
    for item in fixtures["character_classification"]:
        char = item["char"]
        expected = item["class"]
        actual = classify_char(char)
        if actual == expected:
            print(f"  ✓ Character '{char}': expected class {expected}, got {actual}")
        else:
            print(f"  ✗ Character '{char}': expected class {expected}, got {actual}")
            failures += 1

    # 2. Consistency Checks Verification
    print("Running Consistency Checks Tests...")
    for item in fixtures["consistency_checks"]:
        expr = item["expr"]
        expected_dna = item["expected_dna"]
        fail_rule = item["fail_rule"]
        
        # Runs the rules in consistency.py and encodes as bases
        from phi.consistency import check_consistency, CONSISTENCY_RULES
        results_dict = check_consistency(expr)
        actual_dna = "".join(["G" if results_dict[rule] else "T" for rule in CONSISTENCY_RULES])
        
        if actual_dna == expected_dna:
            print(f"  ✓ Expr '{expr}': expected DNA '{expected_dna}' (fail rule: {fail_rule}), got '{actual_dna}'")
        else:
            print(f"  ✗ Expr '{expr}': expected DNA '{expected_dna}', got '{actual_dna}'")
            failures += 1

    # 3. Fixed-Point Arithmetic & Slopes Verification
    print("Running Fixed-Point Slope Tests...")
    for name, item in fixtures["fixed_point_slopes"].items():
        input_scale = item["input_scale"]
        expected_large = item["expected_slope_large"]
        expected_small = item["expected_slope_small"]
        tolerance_lsb = item["tolerance_lsb"]
        
        actual_large, actual_small = parabola_conjugate_pair(input_scale)
        
        # Verify large slope match
        large_match = (actual_large == expected_large)
        # Verify small slope match
        small_match = (actual_small == expected_small)
        
        # Verify orthogonality within tolerance
        prod = q16_mul(actual_large, actual_small)
        target = -65536
        diff = abs(prod - target)
        tol_match = (diff <= tolerance_lsb)
        
        if large_match and small_match and tol_match:
            print(f"  ✓ slopes for {name}: expected ({expected_large}, {expected_small}), got ({actual_large}, {actual_small}) within {tolerance_lsb} LSB tolerance")
        else:
            print(f"  ✗ slopes for {name}: expected ({expected_large}, {expected_small}) with tol {tolerance_lsb}, got ({actual_large}, {actual_small}) diff {diff}")
            failures += 1

    print("===================================")
    if failures == 0:
        print("ALL VERIFIED TEST UNITS PASSED")
        return True
    else:
        print(f"{failures} VERIFIED TEST UNIT(S) FAILED")
        return False

if __name__ == "__main__":
    success = run_verified_tests()
    sys.exit(0 if success else 1)