#!/usr/bin/env python3
"""COMPLEXA C171–C180: reproducible scientific SVGs and numerical source CSVs.

Uses only Python 3.10+ standard library. No stochastic draws, online data,
plotting packages, or opaque external dependencies. Original manuscripts remain
unaltered.  `--check REPO_ROOT` regenerates and compares exact bytes against
content/supplementary/{figures,repro/c171-180}.  `--out PATH` generates them.
"""
from __future__ import annotations
import argparse, csv, io, math, tempfile
from pathlib import Path
from html import escape

BASE = 'complexa-c{}-fig01.svg'
COLORS = {'blue':'#315f88','orange':'#b7783d','teal':'#247e78','purple':'#795b9a','red':'#a4524c','gray':'#647580','mint':'#dbefe9'}
INK='#223643'; LIGHT='#eff3f1'; WHITE='#ffffff'
W,H=1080,640

def txt(x,y,t,size=16,color=INK,anchor='start',weight='normal',rotate=None):
    a=f' transform="rotate({rotate} {x} {y})"' if rotate else ''
    return f'<text x="{x:.2f}" y="{y:.2f}" text-anchor="{anchor}" font-family="Arial, Helvetica, sans-serif" font-size="{size}" font-weight="{weight}" fill="{color}"{a}>{escape(str(t))}</text>'

def polyline(points,color,width=3,dash=None):
    if not points:return ''
    s=' '.join(f'{x:.2f},{y:.2f}' for x,y in points)
    d=f' stroke-dasharray="{dash}"' if dash else ''
    return f'<polyline points="{s}" fill="none" stroke="{color}" stroke-width="{width}" stroke-linecap="round" stroke-linejoin="round"{d}/>'

def line(x1,y1,x2,y2,color='#aab6b9',width=1,dash=None):
    d=f' stroke-dasharray="{dash}"' if dash else ''
    return f'<line x1="{x1:.2f}" y1="{y1:.2f}" x2="{x2:.2f}" y2="{y2:.2f}" stroke="{color}" stroke-width="{width}"{d}/>'

def rect(x,y,w,h,fill,stroke='none',rad=0,opacity=1):
    return f'<rect x="{x:.2f}" y="{y:.2f}" width="{w:.2f}" height="{h:.2f}" rx="{rad}" fill="{fill}" stroke="{stroke}" opacity="{opacity}"/>'

def circle(x,y,r,fill,stroke='none',sw=1):
    return f'<circle cx="{x:.2f}" cy="{y:.2f}" r="{r:.2f}" fill="{fill}" stroke="{stroke}" stroke-width="{sw}"/>'

def figure(title,subtitle,desc,parts):
    return ('<svg xmlns="http://www.w3.org/2000/svg" width="1080" height="640" viewBox="0 0 1080 640" '
            'role="img" aria-labelledby="chart-title chart-desc">\n'
            f'<title id="chart-title">{escape(title)}</title>\n<desc id="chart-desc">{escape(desc)}</desc>\n'
            +rect(0,0,W,H,WHITE)+rect(0,0,12,H,COLORS['teal'])
            +txt(64,59,title,27,INK,weight='bold')+txt(64,87,subtitle,15,COLORS['gray'])
            +''.join(parts)+txt(1018,615,'COMPLEXA · Figura didáctica reproducible',12,COLORS['gray'],anchor='end')+'\n</svg>\n')

def chart(xlab,ylab,xmin,xmax,ymin,ymax,xticks,yticks,parts):
    left,top,width,height=116,145,846,390
    def xf(v):return left+(v-xmin)/(xmax-xmin)*width
    def yf(v):return top+height-(v-ymin)/(ymax-ymin)*height
    p=[rect(left,top,width,height,WHITE)]
    for a,label in yticks:
        y=yf(a);p.append(line(left,y,left+width,y,LIGHT,1))
        p.append(txt(left-12,y+5,label,13,COLORS['gray'],anchor='end'))
    for a,label in xticks:
        x=xf(a);p.append(line(x,top+height,x,top+height+5,'#829397'))
        p.append(txt(x,top+height+24,label,13,COLORS['gray'],anchor='middle'))
    p.append(line(left,top,left,top+height,'#8b9ba0',1.4))
    p.append(line(left,top+height,left+width,top+height,'#8b9ba0',1.4))
    p.extend(parts(xf,yf))
    p.append(txt(left+width/2,589,xlab,15,INK,anchor='middle'))
    p.append(txt(39,top+height/2,ylab,15,INK,anchor='middle',rotate=-90))
    return p

def legend(series,x=740,y=124,spacing=148):
    p=[]
    for i,(label,color,dash) in enumerate(series):
        xx=x+i*spacing;p.append(line(xx,y-5,xx+26,y-5,color,3,dash))
        p.append(txt(xx+33,y,label,13,INK))
    return p

def to_csv(header,rows):
    out=io.StringIO(newline='');w=csv.writer(out,lineterminator='\n');w.writerow(header)
    for r in rows:w.writerow([('%0.10g'%a if isinstance(a,float) else a) for a in r])
    return out.getvalue()

def save(n,header,rows,title,subtitle,desc,plotted):
    return (to_csv(header,rows),figure(title,subtitle,desc,plotted))

def rk4_step(fun,t,y,h):
    k1=fun(t,y);k2=fun(t+h/2,[a+h*b/2 for a,b in zip(y,k1)])
    k3=fun(t+h/2,[a+h*b/2 for a,b in zip(y,k2)]);k4=fun(t+h,[a+h*b for a,b in zip(y,k3)])
    return [a+h*(b+2*c+2*d+e)/6 for a,b,c,d,e in zip(y,k1,k2,k3,k4)]

def c171():
    # Two timescales for reactor: reactant x and product p in mmol/L.
    rows=[]
    for i in range(121):
        t=i/12;x=5*(1-math.exp(-t));z=10*(1-math.exp(-.5*t));p=z-x
        rows.append((t,x,p,z))
    assert abs(rows[-1][3]-10)<.07
    def p(xf,yf):
        a=[]
        for key,color in [(1,COLORS['blue']),(2,COLORS['orange']),(3,COLORS['teal'])]:
            a.append(polyline([(xf(r[0]),yf(r[key])) for r in rows],color,3))
        a.append(line(xf(0),yf(5),xf(10),yf(5),'#9eabb0',1,'4 4'))
        a.extend(legend([('Reactivo x',COLORS['blue'],None),('Producto p',COLORS['orange'],None),('Inventario x + p',COLORS['teal'],None)],x=157,y=135,spacing=234))
        return a
    return save(171,['tiempo_h','reactivo_mmol_L','producto_mmol_L','inventario_mmol_L'],rows,
        'Reactor abierto: dos escalas de relajación','D = k = 0,5 h⁻¹; cₑₙₜ = 10 mmol/L; x(0) = p(0) = 0',
        'Reactivo y producto convergen a cinco milimoles por litro; la suma converge a diez más lentamente, con una escala de dos horas.',
        chart('Tiempo (h)','Concentración (mmol/L)',0,10,0,11,[(i,str(i)) for i in range(0,11,2)],[(i,str(i)) for i in range(0,11,2)],p))

def c172():
    rows=[]
    for regime,a,b in [('Irreversible (sesgado)',3.,1.),('Reversible (simétrico)',2.,2.)]:
        rows.append((regime,a,b,1/3,a/(a+b),(a-b)/3,1/(a+b)))
    assert rows[0][5]==2/3 and rows[1][5]==0
    # two simple comparisons, same occupancy and same event rate
    p=[]
    p.append(txt(110,150,'Ocupación estacionaria',18,INK,weight='bold'))
    p.append(txt(606,150,'Corriente neta por enlace',18,INK,weight='bold'))
    for i,lab in enumerate(['Anillo sesgado 3 : 1','Anillo simétrico 2 : 2']):
        yy=216+180*i
        p.append(txt(110,yy-13,lab,16))
        for z in range(3):
            xx=110+z*130;p.append(rect(xx,yy,90,35,COLORS['blue'] if i==0 else COLORS['teal'],rad=3));p.append(txt(xx+45,yy+24,'⅓',18,WHITE,anchor='middle'))
        p.append(txt(605,yy-13,lab,16))
        p.append(rect(605,yy,295,34,'#e9efee',rad=3))
        width=295*rows[i][5]/(2/3)
        if width>0:p.append(rect(605,yy,width,34,COLORS['orange'],rad=3))
        p.append(txt(913,yy+24,('2/3' if i==0 else '0')+' s⁻¹',14,INK))
    p.append(txt(112,561,'π = (⅓, ⅓, ⅓) en ambos regímenes',15,COLORS['gray']))
    p.append(txt(606,561,'Actividad de salida idéntica: 4 s⁻¹',15,COLORS['gray']))
    return save(172,['modelo','tasa_horaria_s_inv','tasa_antihoraria_s_inv','pi_por_estado','probabilidad_horaria','corriente_neta_s_inv','espera_media_s'],rows,
        'Mismas ocupaciones, distintas corrientes','Anillo de tres estados · tasas (3, 1) frente a (2, 2) s⁻¹',
        'Los tres estados tienen probabilidad un tercio en ambos anillos. Solo el anillo asimétrico tiene circulación estacionaria neta de dos tercios por segundo.',p)

def c173():
    rows=[]
    for i in range(161):
        th=300+i/4;sigma=10*(1/300-1/th);rows.append((th,sigma,sigma*1000))
    target=10*(1/300-1/310)*1000
    assert abs(target-1.075268817)<1e-7
    def p(xf,yf):
        parts=[polyline([(xf(r[0]),yf(r[2])) for r in rows],COLORS['teal'],3)]
        parts.append(line(xf(310),yf(0),xf(310),yf(target),'#829398',1.5,'5 5'))
        parts.append(circle(xf(310),yf(target),5,COLORS['orange'],WHITE,1.5))
        parts.append(txt(xf(310)+12,yf(target)-12,'310 K → 1,075 mW/K',14,COLORS['orange'],weight='bold'))
        return parts
    return save(173,['reservorio_caliente_K','produccion_W_K','produccion_mW_K'],rows,
        'Producción entrópica en dos reservorios','Potencia térmica fija: Q̇ = 10 W · Reservorio frío: 300 K',
        'Producción de entropía 10 por diferencia de temperaturas recíprocas; se anula si ambas temperaturas son iguales y aumenta al calentar el reservorio caliente.',
        chart('Temperatura del reservorio caliente (K)','Producción de entropía (mW/K)',300,340,0,4,[(v,str(v)) for v in range(300,341,10)],[(v,str(v)) for v in [0,1,2,3,4]],p))

def c174():
    rows=[];state=[0.,0.];h=.01
    def f(t,z):
        c,d=z;a=1-c;b=1-c-d
        return [a*b-2*c,c-d]
    for i in range(1001):
        t=i*h;c,d=state;a,b=1-c,1-c-d
        rows.append((t,a,b,c,d,a+c,b+c+d))
        state=rk4_step(f,t,state,h)
    assert max(abs(r[5]-1)+abs(r[6]-1) for r in rows)<1e-12
    expected=(5-math.sqrt(17))/4
    assert abs(rows[-1][3]-expected)<0.0002
    def p(xf,yf):
        out=[]
        for idx,color in [(1,COLORS['blue']),(2,COLORS['orange']),(3,COLORS['teal']),(4,COLORS['purple'])]:
            out.append(polyline([(xf(r[0]),yf(r[idx])) for r in rows[::3]],color,3))
        out.extend(legend([('A',COLORS['blue'],None),('B',COLORS['orange'],None),('C',COLORS['teal'],None),('D',COLORS['purple'],None)],x=380,y=134,spacing=115))
        out.append(txt(650,196,'A + C = B + C + D = 1',14,COLORS['gray']))
        return out
    return save(174,['tiempo_adimensional','A','B','C','D','conservacion_A_C','conservacion_B_C_D'],rows,
        'Red de cuatro especies: dos inventarios constantes','A(0) = B(0) = 1 · C(0) = D(0) = 0 · k₁ = k₂ = k₃ = k₄ = 1',
        'Cuatro concentraciones cambian con el tiempo, mientras A más C y B más C más D permanecen exactamente iguales a uno.',
        chart('Tiempo adimensional','Concentración normalizada',0,10,0,1.05,[(v,str(v)) for v in range(0,11,2)],[(v,str(v/10)) for v in [0,2,4,6,8,10]],p))

def c175():
    initials=[.1,.19,.21,.4,.9]
    states=initials[:];h=.02;rows=[]
    def f(t,z):return [x*x*(1-x)-.16*x for x in z]
    for i in range(4001):
        if i%20==0:rows.append((i*h,*states))
        states=rk4_step(f,i*h,states,h)
        assert all(-1e-10<=x<=1.000001 for x in states)
    assert rows[-1][1]<.001 and rows[-1][4]>.79
    def p(xf,yf):
        out=[line(xf(0),yf(.2),xf(80),yf(.2),'#818f94',1.5,'6 5'),line(xf(0),yf(.8),xf(80),yf(.8),'#818f94',1.5,'6 5')]
        for i,col in enumerate([COLORS['blue'],COLORS['purple'],COLORS['teal'],COLORS['orange'],COLORS['red']]):
            out.append(polyline([(xf(r[0]),yf(r[i+1])) for r in rows],col,2.7))
        out.append(txt(xf(79),yf(.8)-8,'Estable: X = 0,8',14,COLORS['gray'],anchor='end'))
        out.append(txt(xf(79),yf(.2)-8,'Umbral: X = 0,2',14,COLORS['gray'],anchor='end'))
        out.extend(legend([(str(a).replace('.',','),col,None) for a,col in zip(initials,[COLORS['blue'],COLORS['purple'],COLORS['teal'],COLORS['orange'],COLORS['red']])],x=215,y=134,spacing=135))
        return out
    return save(175,['tiempo_adimensional']+[f'X_inicial_{a:.2f}' for a in initials],rows,
        'Autocatálisis cooperativa: un umbral separa destinos','Ẋ = X²(1 − X) − 0,16X · Inventario total y k normalizados a 1',
        'Soluciones deterministas desde cinco semillas iniciales. Por debajo de 0,2 se extinguen; por encima convergen a 0,8.',
        chart('Tiempo adimensional','Concentración X',0,80,0,1.0,[(v,str(v)) for v in range(0,81,20)],[(v/10,str(v/10)) for v in range(0,11,2)],p))

def c176():
    rows=[]
    for i in range(161):
        b=1.2+i*.01;trace=b-2;re=trace/2;im=math.sqrt(max(0,1-(trace/2)**2));rows.append((b,trace,re,im))
    assert abs(rows[80][2])<1e-12
    def p(xf,yf):
        parts=[line(xf(1.2),yf(0),xf(2.8),yf(0),'#849699',1.5),line(xf(2),yf(-.5),xf(2),yf(.5),'#819396',1,'5 4')]
        parts.append(polyline([(xf(r[0]),yf(r[2])) for r in rows],COLORS['teal'],3))
        parts.extend([circle(xf(b),yf((b-2)/2),5,COLORS['orange'],WHITE,1.3) for b in [1.5,2,2.2]])
        parts.append(txt(xf(2)+14,yf(.5)+20,'Hopf: b = 2',14,COLORS['orange'],weight='bold'))
        parts.append(txt(xf(1.4),yf(-.25),'Foco estable',15,COLORS['blue']))
        parts.append(txt(xf(2.22),yf(.27),'Foco inestable',15,COLORS['red']))
        return parts
    return save(176,['parametro_b','traza','real_autovalor','frecuencia_angular'],rows,
        'Brusselator: cruce de autovalores en Hopf','Parámetro a = 1 · b crítico = 1 + a² = 2',
        'Parte real de los autovalores complejos del punto fijo frente al parámetro b; cambia de negativa a positiva en dos.',
        chart('Parámetro b','Re(λ) · tiempo⁻¹',1.2,2.8,-.5,.5,[(v/10,str(v/10)) for v in [12,16,20,24,28]],[(v/10,str(v/10)) for v in [-5,-25,0,25,5] if abs(v)!=25]+[(-.25,'−0,25'),(.25,'0,25')],p))

def c177():
    rows=[]
    for i in range(241):
        z=i*.005;trace=-.2-11*z;det=1-7*z+10*z*z;dis=trace*trace-4*det
        growth=(trace+math.sqrt(dis))/2 if dis>=0 else trace/2
        rows.append((z,trace,det,growth))
    z5=(math.pi/5)**2;z3=(math.pi/3)**2
    def growth(z):
        t=-.2-11*z;d=1-7*z+10*z*z;dis=t*t-4*d
        return (t+math.sqrt(dis))/2 if dis>=0 else t/2
    assert abs(growth(z5)-.0447)<.001
    assert growth(z3)<0
    def p(xf,yf):
        a=[rect(xf(.2),yf(.14),xf(.5)-xf(.2),yf(-.7)-yf(.14),'#e4f3ef',opacity=.9)]
        a.append(line(xf(0),yf(0),xf(1.2),yf(0),'#81999a',1.4))
        a.append(polyline([(xf(r[0]),yf(r[3])) for r in rows],COLORS['teal'],3))
        for z,lab,col in [(z5,'L = 5 · n = 1',COLORS['orange']),(z3,'L = 3 · n = 1',COLORS['blue'])]:
            a.append(circle(xf(z),yf(growth(z)),5.5,col,WHITE,1.5));a.append(txt(xf(z)+10,yf(growth(z))-13,lab,13,col,weight='bold'))
        a.append(txt(xf(.35),yf(.125),'Banda de Turing',14,COLORS['teal'],anchor='middle',weight='bold'))
        return a
    # Negative eigenvalue may reach ~ -0.4; clip to chosen chart if needed at right.
    return save(177,['numero_onda_cuadrado','traza_modal','determinante_modal','crecimiento_maximo'],rows,
        'Turing: banda continua y modos espaciales permitidos','Brusselator a = 1, b = 1,8 · Dₓ = 1, Dᵧ = 10 · Neumann',
        'Tasa de crecimiento modal positiva solo para q al cuadrado entre 0,2 y 0,5; el primer modo del dominio de longitud cinco entra en la banda y el de longitud tres no.',
        chart('Número de onda al cuadrado (q²)','Re(λₘₐₓ) · tiempo⁻¹',0,1.2,-.7,.14,[(v/10,str(v/10)) for v in range(0,13,2)],[(-.6,'−0,6'),(-.4,'−0,4'),(-.2,'−0,2'),(0,'0'),(.1,'0,1')],p))

def c178():
    norm=lambda x,mu: math.exp(-(x-mu)**2/4)/(math.sqrt(4*math.pi))
    rows=[]
    for i in range(201):
        w=-3+i*.045;pF=norm(w,2);pR=norm(w,0);weighted=pF*math.exp(-w+1)
        assert abs(weighted-pR)<1e-12
        rows.append((w,pF,pR,weighted))
    def p(xf,yf):
        a=[polyline([(xf(r[0]),yf(r[1])) for r in rows],COLORS['blue'],3),polyline([(xf(r[0]),yf(r[2])) for r in rows],COLORS['orange'],3)]
        a.append(line(xf(1),yf(0),xf(1),yf(norm(1,2)),'#84959a',1.4,'5 4'))
        a.append(circle(xf(1),yf(norm(1,2)),5,COLORS['teal'],WHITE,1.5))
        a.append(txt(xf(1)+12,yf(.22),'Cruce: ΔF = 1 E₀',15,COLORS['teal'],weight='bold'))
        a.extend(legend([('Trabajo directo',COLORS['blue'],None),('Inverso reflejado',COLORS['orange'],None)],x=440,y=128,spacing=245))
        return a
    return save(178,['trabajo_en_E0','densidad_directa','densidad_inversa_reflejada','directa_ponderada_normalizada'],rows,
        'Crooks y Jarzynski: el papel de la cola de bajo trabajo','Modelo gaussiano: Wᴳ ~ N(2, 2) · β = 1 / E₀ · ΔF = 1 E₀',
        'Densidades gaussianas de trabajo directo y trabajo inverso con signo invertido; cruzan en uno y la distribución directa reponderada exponencialmente coincide con la inversa reflejada.',
        chart('Trabajo W / E₀','Densidad por E₀',-3,6,0,.32,[(v,str(v)) for v in [-3,-1,1,3,5]],[(v,str(v)) for v in [0,.08,.16,.24,.32]],p))

def c179():
    states=[('Estado estable',(1,1,1),0,0,0),('Adición al centro',(1,2,1),0,0,0),('Vuelco 1: centro',(2,0,2),1,0,1),('Vuelco 2: extremos',(0,2,0),2,2,3),('Vuelco 3: centro',(1,0,1),1,0,4)]
    assert states[-1][1]==(1,0,1) and sum(s[3] for s in states)==2
    rows=[(i,label,*st,topples,lost,cumulative) for i,(label,st,topples,lost,cumulative) in enumerate(states)]
    p=[txt(122,147,'Actualización paralela (sitios 1 – 2 – 3)',17,INK,weight='bold')]
    heights={0:'#e8eeeb',1:'#76a9a3',2:'#b8723d'}
    for i,(label,st,topples,lost,cum) in enumerate(states):
        y=181+i*75;p.append(txt(116,y+31,label,14,INK))
        for k,height in enumerate(st):
            x=548+k*125;p.append(rect(x,y,70,48,heights[height],rad=5))
            p.append(txt(x+35,y+32,str(height),22,WHITE if height>0 else INK,anchor='middle',weight='bold'))
        p.append(txt(960,y+30,(f'{topples} vuelcos' if i>=2 else '—'),13,COLORS['gray'],anchor='end'))
    p.append(txt(125,594,'Total: 4 vuelcos · área 3 sitios · 3 generaciones · 2 granos perdidos',14,COLORS['gray']))
    return save(179,['paso','etiqueta','altura_1','altura_2','altura_3','vuelcos_generacion','perdida_generacion','vuelcos_acumulados'],rows,
        'Pila de arena: la secuencia de una avalancha','Cadena de tres sitios · umbral 2 · extremos abiertos · un grano al centro',
        'Diagrama discreto de cinco configuraciones; la adición al centro causa cuatro vuelcos repartidos en tres generaciones, concluye en alturas uno cero uno y pierde dos granos.',p)

def c180():
    rows=[]
    for i in range(241):
        t=i*.5;p=.1;m=math.sqrt(.01/(2*math.exp(.02*t)-1))
        rows.append((t,p,m,p*math.cos(t),m*math.cos(t)))
    assert abs(rows[20][2]-.0832523)<.000001
    assert abs(rows[200][2]-.0269405)<.000001
    def p(xf,yf):
        a=[line(xf(0),yf(.1),xf(120),yf(.1),'#91a1a5',1,'4 4')]
        a.append(polyline([(xf(r[0]),yf(r[1])) for r in rows],COLORS['teal'],3))
        a.append(polyline([(xf(r[0]),yf(r[2])) for r in rows],COLORS['orange'],3))
        for t in [10,100]:
            val=math.sqrt(.01/(2*math.exp(.02*t)-1));a.append(circle(xf(t),yf(val),5,COLORS['orange'],WHITE,1))
            a.append(txt(xf(t)+10,yf(val)-14,f't = {t}: {val:.4f}',13,COLORS['orange']))
        a.extend(legend([('μ = +0,01',COLORS['teal'],None),('μ = −0,01',COLORS['orange'],None)],x=520,y=130,spacing=218))
        return a
    return save(180,['tiempo_adimensional','amplitud_regimen_positivo','amplitud_regimen_negativo','observable_positivo','observable_negativo'],rows,
        'Oscilación sostenida y transitorio prolongado','Forma normal radial · r(0) = 0,1 · ω = 1 · μ = ±0,01',
        'Amplitud constante 0,1 en régimen positivo y decaimiento gradual a cero en régimen negativo; al tiempo diez, la negativa todavía vale 0,08325.',
        chart('Tiempo adimensional','Amplitud radial r(t)',0,120,0,.11,[(i,str(i)) for i in range(0,121,20)],[(i/100,str(i/100)) for i in [0,2,4,6,8,10]],p))

FUNCTIONS=[c171,c172,c173,c174,c175,c176,c177,c178,c179,c180]

def files():
    for i,fn in enumerate(FUNCTIONS):
        n=171+i; csv_data, svg=fn()
        yield f'repro/c171-180/c{n}.csv',csv_data.encode('utf-8')
        yield f'figures/complexa-c{n}-fig01.svg',svg.encode('utf-8')

def main():
    parser=argparse.ArgumentParser()
    group=parser.add_mutually_exclusive_group(required=True)
    group.add_argument('--out',type=Path);group.add_argument('--check',type=Path)
    args=parser.parse_args();root=args.out or args.check
    errors=[]
    for rel,content in files():
        dest=root/rel if args.out else root/'content/supplementary'/rel
        if args.out:
            dest.parent.mkdir(parents=True,exist_ok=True);dest.write_bytes(content)
        else:
            if not dest.exists():errors.append(f'MISSING {dest}')
            elif dest.read_bytes()!=content:errors.append(f'DIFF {dest}')
    if errors:raise RuntimeError('\n'.join(errors))
    print('PASS: 10 deterministic SVG figures + 10 numeric CSVs; analytical checks and source invariants verified.')

if __name__=='__main__':main()