"""Owned unpreconditioned exact-rational CG for symmetric positive-definite 2x2."""
from fractions import Fraction as F
import json
def dot(a,b):return sum(x*y for x,y in zip(a,b))
def mv(A,x):return [dot(row,x) for row in A]
def cg(A,b,x,maxiter):
    if len(A)!=2 or any(len(r)!=2 for r in A) or len(b)!=2 or len(x)!=2:raise ValueError('2-by-2 subset')
    if A[0][1]!=A[1][0] or A[0][0]<=0 or A[0][0]*A[1][1]-A[0][1]**2<=0:raise ValueError('2-by-2 symmetric positive definite required')
    if type(maxiter)is not int or maxiter<0:raise ValueError('nonnegative iteration cap')
    x=x[:];r=[a-c for a,c in zip(b,mv(A,x))];p=r[:];trace=[]
    for k in range(maxiter):
        rr=dot(r,r)
        if rr==0:break
        Ap=mv(A,p);den=dot(p,Ap)
        if den<=0:raise ValueError('nonpositive direction curvature')
        alpha=rr/den;xnew=[a+alpha*d for a,d in zip(x,p)];rnew=[a-alpha*d for a,d in zip(r,Ap)];beta=dot(rnew,rnew)/rr
        trace.append({'step':k+1,'x_before':x,'r_before':r,'p':p,'alpha':alpha,'x_after':xnew,'r_after':rnew,'beta':beta})
        p=[a+beta*d for a,d in zip(rnew,p)];x,r=xnew,rnew
    actual=[a-c for a,c in zip(b,mv(A,x))];assert actual==r
    return {'x':x,'residual':r,'residual_squared_norm':dot(r,r),'converged_exact_zero_residual':dot(r,r)==0,'iteration_cap':maxiter,'steps':len(trace),'trace':trace}
def demo():
    A=[[F(4),F(1)],[F(1),F(3)]];b=[F(1),F(2)];zero=[F(0),F(0)];full=cg(A,b,zero,2);short=cg(A,b,zero,1)
    direct=[(A[1][1]*b[0]-A[0][1]*b[1])/11,(A[0][0]*b[1]-A[1][0]*b[0])/11]
    assert full['x']==direct==[F(1,11),F(7,11)] and full['converged_exact_zero_residual']
    assert not short['converged_exact_zero_residual'] and short['residual_squared_norm']==F(5,16)
    t=full['trace'];r_orth=dot(t[0]['r_before'],t[0]['r_after']);a_conj=dot(t[0]['p'],mv(A,t[1]['p']));assert r_orth==a_conj==0
    already=cg(A,b,direct,2);assert already['steps']==0
    try:cg([[F(1),F(0)],[F(0),F(-1)]],b,zero,2)
    except ValueError:rejected=True
    else:raise AssertionError('indefinite accepted')
    return {'A':A,'b':b,'two_step_run':full,'one_step_cap_run':short,'direct_solution':direct,'successive_residual_dot':r_orth,'first_two_directions_A_inner_product':a_conj,'already_solved_run':already,'indefinite_rejected_by_full_2by2_SPD_check':rejected,'arithmetic':'exact Fraction, not SciPy tolerance behavior','SciPy_preconditioning_or_large_problem_executed':False}
if __name__=='__main__':print(json.dumps(demo(),default=str,indent=2))