"""Source-based rational checks. General validity rests on the separate algebraic audit."""
from fractions import Fraction as F
from pathlib import Path
import random,json,hashlib
def mm(a,b):return [[sum(x*y for x,y in zip(row,col)) for col in zip(*b)] for row in a]
def avg(rows):return [sum(col)/len(rows) for col in zip(*rows)]
def normalize(row):return [F(x,sum(row)) for x in row]
def check(w,beta,gamma,t,q):
 n=len(w);a=[1-x for x in beta]+[gamma]
 assert n>=1 and len(t)==n+1 and all(len(row)==n for row in w)
 assert all(sum(row)==1 and min(row)>0 for row in w)
 assert len(q)==n and sum(q)==1 and min(q)>0
 assert all(0<x<1 for x in beta) and 0<gamma<1
 assert all(len(row)==n+1 and min(row)>=0 and sum(row)==1 for row in t)
 j=[[(1-beta[i])*x for x in row]+[beta[i]] for i,row in enumerate(w)]+[[gamma*x for x in q]+[1-gamma]]
 z=mm(t,j);u=[[x/sum(row[:n]) for x in row[:n]] for row in z[:n]]
 d=[sum(t[i][r]*a[r] for r in range(n+1)) for i in range(n)]
 theta=[a[r]*sum(t[i][r]/d[i] for i in range(n))/n for r in range(n+1)]
 b=[avg(w)[k]-sum(theta[r]*w[r][k] for r in range(n)) for k in range(n)]
 qstar=[x/theta[-1] for x in b] if theta[-1] else None
 expected=[sum(theta[r]*w[r][k] for r in range(n))+theta[-1]*q[k] for k in range(n)]
 assert avg(u)==expected and sum(theta)==1 and min(theta)>=0
 assert all(0<x<1 for x in d) and all(min(row)>0 and sum(row)==1 for row in u)
 if theta[-1]:
  assert sum(qstar)==1
  assert [avg(u)[k]-avg(w)[k] for k in range(n)]==[theta[-1]*(q[k]-qstar[k]) for k in range(n)]
  assert (avg(u)==avg(w))==(q==qstar)
 else:assert (avg(u)==avg(w))==all(x==0 for x in b)
 return {'theta':theta,'b':b,'qstar':qstar,'U':u,'mean':avg(u),'error':[avg(u)[k]-avg(w)[k] for k in range(n)]}
def plain(x):
 if isinstance(x,F):return str(x)
 if isinstance(x,list):return [plain(v) for v in x]
 if isinstance(x,dict):return {k:plain(v) for k,v in x.items()}
 return x
rng=random.Random(410420260915);cases=[]
for n in range(1,6):
 for iteration in range(40):
  w=[normalize([rng.randint(1,9) for _ in range(n)]) for _ in range(n)]
  beta=[F(rng.randint(1,8),9) for _ in range(n)];gamma=F(rng.randint(1,8),9)
  t=[]
  for _ in range(n+1):
   raw=[rng.randint(0,4) for _ in range(n+1)]
   if not sum(raw):raw[0]=1
   t.append(normalize(raw))
  q=normalize([rng.randint(1,9) for _ in range(n)]);out=check(w,beta,gamma,t,q)
  if out['qstar'] and min(out['qstar'])>0:assert check(w,beta,gamma,t,out['qstar'])['mean']==avg(w)
  cases.append({'n':n,'iteration':iteration,'W':w,'beta':beta,'gamma':gamma,'T':t,'q':q,'result':out})
# Original equations9-12.
w=[[F(9,10),F(1,10)],[F(1,10),F(9,10)]]
p=[[F(1,2) if i==j else F(1,4) for j in range(3)] for i in range(3)]
eq12=check(w,[F(1,3)]*2,F(7,10),p,[F(6,7),F(1,7)])
assert eq12['U']==[[F(56,81),F(25,81)],[F(40,81),F(41,81)]] and eq12['mean']==[F(16,27),F(11,27)]
# Original equation19 and its one/two-step variants.
w=[[F(3,4),F(1,4)],[F(1,4),F(3,4)]];p=[[F(3,4),F(1,4),0],[F(1,4),F(1,2),F(1,4)],[0,F(1,4),F(3,4)]]
one=check(w,[F(1,3)]*2,F(2,3),p,[F(1,4),F(3,4)])
two=check(w,[F(1,3)]*2,F(2,3),mm(p,p),[F(1,4),F(3,4)])
assert one['U']==[[F(5,8),F(3,8)],[F(3,8),F(5,8)]]
assert two['mean']==[F(31,64),F(33,64)] and two['qstar']==[F(1,3),F(2,3)]
assert check(w,[F(1,3)]*2,F(2,3),mm(p,p),two['qstar'])['mean']==avg(w)
# Theta_new=0 does NOT itself imply mean preservation for arbitrary T.
zero_shift=check(w,[F(1,2)]*2,F(1,2),[[1,0,0],[1,0,0],[0,0,1]],[F(1,2)]*2)
assert zero_shift['theta'][-1]==0 and zero_shift['error']!=[0,0]
# Nonmonotone stay even inside strict-positive-diagonal Laplacian specialization.
p=[[F(1,5) if i==j else F(2,5) for j in range(3)] for i in range(3)]
t=[[F(i==j) for j in range(3)] for i in range(3)];stays=[]
for ell in [1,2,3]:
 t=mm(t,p);stays.append(check(w,[F(1,2)]*2,F(1,2),t,[F(3,4),F(1,4)]))
assert [x['theta'][-1] for x in stays]==[F(2,5),F(8,25),F(42,125)]
assert max(map(abs,stays[1]['error']))<max(map(abs,stays[0]['error']))
# Equal-retention corollary across arbitrary convex combinations of permutations.
corollary_cases=[]
for n in range(1,6):
 for it in range(12):
  w=[normalize([rng.randint(1,9) for _ in range(n)]) for _ in range(n)]
  t=[[F(0) for _ in range(n+1)] for _ in range(n+1)]
  for _ in range(3):
   permutation=list(range(n+1));rng.shuffle(permutation)
   for i,j in enumerate(permutation):t[i][j]+=F(1,3)
  b=F(2,5);q=normalize([rng.randint(1,9) for _ in range(n)]);out=check(w,[b]*n,1-b,t,q);s=1-t[-1][-1]
  if s:
   expected=[sum(t[-1][j]*w[j][k] for j in range(n))/s for k in range(n)]
   assert out['qstar']==expected and min(expected)>0
   assert check(w,[b]*n,1-b,t,expected)['mean']==avg(w)
  else:assert out['mean']==avg(w)
  corollary_cases.append({'n':n,'iteration':it,'W':w,'T':t,'s':s,'result':out})
result={'version':'1.1','status':'PASS_EXACT_CHECKS','seed':410420260915,'general_cases':len(cases),'corollary_cases':len(corollary_cases),'n_range':[1,5],'original_eq12':eq12,'original_eq19_one':one,'original_eq19_two':two,'zero_influence_mean_shift':zero_shift,'strict_laplacian_nonmonotone_stays':stays,'cases':cases,'corollary_records':corollary_cases,'limits':'Finite checks, not a proof of general validity, novelty, physical implementation or applicability to human/LLM auditors. No independent human review.'}
Path(__file__).with_name('THEOREM_SOURCE_CHECKS_v1.1.json').write_bytes((json.dumps(plain(result),indent=2)+'\n').encode())
print(json.dumps({'status':'PASS','general_cases':len(cases),'corollary_cases':len(corollary_cases),'original_examples':'PASS','strict_laplacian_counterexample':'PASS'}))