"""Независимые точные решения кольцевого и плоского слоистых зазоров."""
import math,json
MU0=4e-7*math.pi
source=json.load(open('ring-source.json'))['records'];base=next(r for r in source if r['id']=='Кольцо-Угол-90')
native_ring=json.load(open('ring-precision-result.json'))['records'];native_plane=json.load(open('planar-precision-result.json'))['records']
def linear(M,Y):
    n=len(Y);A=[list(M[k])+[Y[k]]for k in range(n)]
    for k in range(n):
        pivot=max(range(k,n),key=lambda j:abs(A[j][k]));A[k],A[pivot]=A[pivot],A[k]
        v=A[k][k];A[k]=[z/v for z in A[k]]
        for j in range(n):
            if j!=k:
                v=A[j][k];A[j]=[A[j][l]-v*A[k][l]for l in range(n+1)]
    return [A[k][-1]for k in range(n)]
def dot(a,b):return sum(a[k]*b[k]for k in range(len(a)))
def reference(gap,curved=True):
    p=3;R=(100-gap/2)*.001;k=p/R;ri=100-gap
    bounds=[v*.001 for v in ([ri,ri+.3,ri+.7,99.3,99.7,100]if curved else[0,.3,.7,gap-.7,gap-.3,gap])]
    J=[0,-1j*base['inner_amp_MApm2']*1e6,0,base['outer_amp_MApm2']*1e6,0]
    def fun(r,d):
        if curved:return [p*r**(p-1),-p*r**(-p-1)]if d else[r**p,r**(-p)]
        return [k*math.exp(k*r),-k*math.exp(-k*r)]if d else[math.exp(k*r),math.exp(-k*r)]
    def part(r,j,d):return MU0*j*(2*r if d else r*r)/(p*p-4)if curved else(0 if d else MU0*j/(k*k))
    M=[[0]*10 for _ in range(10)];Y=[0j]*10
    M[0][0:2]=fun(bounds[0],1);M[1][8:10]=fun(bounds[-1],1);row=2
    for i,r in enumerate(bounds[1:-1]):
        for d in[0,1]:
            M[row][2*i:2*i+2]=fun(r,d);M[row][2*i+2:2*i+4]=[-v for v in fun(r,d)];Y[row]=part(r,J[i+1],d)-part(r,J[i],d);row+=1
    c=linear(M,Y);r=R if curved else gap*.0005;az=dot(fun(r,0),c[4:6]);der=dot(fun(r,1),c[4:6]);br=1j*p/R*az;bt=-der
    return math.pi*.1*R*R/MU0*(br*bt.conjugate()).real
out=[]
for gap in[10,8]:
    exact_ring=reference(gap,True);exact_plane=reference(gap,False)
    nr=next(r['T_stress_Nm']for r in native_ring if r['gap_mm']==gap);np=next(r['T_equivalent_Nm']for r in native_plane if r['gap_mm']==gap)
    out.append({'gap_mm':gap,'g_over_R':gap/(100-gap/2),'exact_ring_Nm':exact_ring,'exact_plane_Nm':exact_plane,'native_ring_Nm':nr,'native_plane_Nm':np,'exact_curvature_relative_error':exact_plane/exact_ring-1,'native_curvature_relative_difference':np/nr-1,'native_ring_error':nr/exact_ring-1,'native_plane_error':np/exact_plane-1})
with open('comparison.json','w')as f:json.dump({'comparisons':out},f)
print('RING_REFERENCE_COMPLETE')