#!/usr/bin/env python3
"""Exact arithmetic check of the sharp 2x2 similarity-metric control."""

from fractions import Fraction as F


def matmul(a, b):
    return [[sum(a[i][k] * b[k][j] for k in range(2)) for j in range(2)] for i in range(2)]


def transpose(a):
    return [list(row) for row in zip(*a)]


N = [[F(0), F(2)], [F(0), F(0)]]
P = [[F(1), F(0)], [F(0), F(4)]]
stein = matmul(transpose(N), matmul(P, N))
defect = [[P[i][j] - stein[i][j] for j in range(2)] for i in range(2)]

assert stein == [[F(0), F(0)], [F(0), F(4)]]
assert defect == [[F(1), F(0)], [F(0), F(0)]]

print("N^* P N =", stein)
print("P - N^* P N =", defect)
print("kappa(P) = 4; hence sqrt(kappa(P)) = 2")
print("Sharpness: the (2,2) Stein entry forces p22 >= 4 p11,")
print("so lambda_max(P) >= 4 lambda_min(P) for every feasible P.")
