#!/usr/bin/env python3
"""COMPLEXA C281-C290: source-derived, deterministic, standard-library figures.

Usage:
    python make_figures.py --out content/supplementary
    python make_figures.py --check content/supplementary

Data are exact evaluations of explicit didactic models, not empirical samples.
SVG is generated directly so graphical elements remain vector and searchable.
"""
from __future__ import annotations
import argparse
import base64
import csv
from hashlib import sha256
from html import escape
from io import BytesIO, StringIO
from math import exp, log, sqrt, isclose
from xml.etree import ElementTree
from pathlib import Path
import zipfile

ROOT = Path(__file__).resolve().parents[2]
NAMES = {
    281: ("Picos y rutas en un paisaje de ocho estados", "Tres posiciones binarias; aptitud y búsqueda por una mutación"),
    282: ("Una variante modifica su propia ventaja", "Tasas dependientes de la frecuencia; equilibrio interior"),
    283: ("La memoria ambiental cambia el destino", "Misma frecuencia inicial y distinto entorno heredado"),
    284: ("Selección de modelos y transformación cultural", "Frecuencias esperadas entre fuentes y productos transmitidos"),
    285: ("Una mejora global pese al declive local", "Composición celular y diferente contribución de grupos"),
    286: ("Especializarse no siempre rinde más", "Rendimientos productivos, costos de intercambio y umbrales"),
    287: ("Genética, cultura e infraestructura coevolucionan", "Tres variables, historia ambiental y retroalimentación"),
    288: ("Retener mejor hace crecer la acumulación", "Ganancia constante y pérdidas durante la transmisión"),
    289: ("Diversidad, encuentros y recombinación", "Probabilidad exacta de una oportunidad en diez encuentros"),
    290: ("Más solicitudes cambian la carga del servicio", "Balance de pedidos sintéticos; capacidad y trabajo equivalente"),
}
NAVY="#183448"; TEAL="#087F87"; OCHRE="#C56D35"; PURPLE="#7465AD"
INK="#263B4A"; MUTED="#66808B"; GRID="#E2EAEF"; BG="#F5F8FA"; GREEN="#4A9870"


def v(x):
    if isinstance(x,str): return x
    if isinstance(x,int): return str(x)
    return f"{x:.10f}"


def csv_bytes(header, rows):
    f=StringIO(newline="")
    w=csv.writer(f,lineterminator="\n")
    w.writerow(header)
    for row in rows: w.writerow([v(x) for x in row])
    return f.getvalue().encode("utf-8")


def X(text): return escape(str(text),quote=True)


class SVG:
    def __init__(self,n,title,desc):
        self.n=n; self.a=[]
        self.a.append(f'<svg xmlns="http://www.w3.org/2000/svg" width="1200" height="760" viewBox="0 0 1200 760" role="img" aria-labelledby="ttl{n} des{n}">')
        self.a.append(f'<title id="ttl{n}">Figura {n}.1. {X(title)}</title>')
        self.a.append(f'<desc id="des{n}">{X(desc)}. Valores de modelos definidos en el manuscrito, no resultados observacionales.</desc>')
        self.rect(0,0,1200,760,"#FFFFFF")
        self.rect(0,0,13,760,TEAL)
        self.rect(43,28,1114,99,BG,15)
        self.text(65,56,"C O M P L E X A   /   LABORATORIO REPRODUCIBLE",13,TEAL,700)
        self.text(65,93,title,29,NAVY,700)
        self.text(1133,62,f"C{n}",24,NAVY,700,anchor="end")
        self.text(68,150,desc,17,MUTED,400)
        self.line(55,691,1153,691,GRID,2)
        self.text(64,718,"MODELO MATEMÁTICO / DATOS DIDÁCTICOS",12,MUTED,700)
        self.text(1138,718,f"complexa.science  ·  {n}.1",13,MUTED,500,anchor="end")
    def rect(self,x,y,w,h,color,r=0,stroke=None,sw=1):
        self.a.append(f'<rect x="{x}" y="{y}" width="{w}" height="{h}" rx="{r}" fill="{color}"'+(f' stroke="{stroke}" stroke-width="{sw}"' if stroke else "")+'/>')
    def text(self,x,y,s,size=16,color=INK,weight=400,anchor="start",rotate=None):
        rot=f' transform="rotate({rotate} {x} {y})"' if rotate is not None else ""
        self.a.append(f'<text x="{x}" y="{y}" text-anchor="{anchor}" fill="{color}" font-family="DejaVu Sans, Arial, sans-serif" font-size="{size}" font-weight="{weight}"{rot}>{X(s)}</text>')
    def line(self,x1,y1,x2,y2,color,sw=1,dash=None):
        self.a.append(f'<line x1="{x1:.2f}" y1="{y1:.2f}" x2="{x2:.2f}" y2="{y2:.2f}" stroke="{color}" stroke-width="{sw}" stroke-linecap="round"'+(f' stroke-dasharray="{dash}"' if dash else "")+'/>')
    def circle(self,x,y,r,color,stroke=None):
        self.a.append(f'<circle cx="{x:.2f}" cy="{y:.2f}" r="{r}" fill="{color}"'+(f' stroke="{stroke}"' if stroke else "")+'/>')
    def path(self,p,color,sw=3,fill="none",dash=None):
        self.a.append(f'<path d="{p}" stroke="{color}" stroke-width="{sw}" stroke-linejoin="round" stroke-linecap="round" fill="{fill}"'+(f' stroke-dasharray="{dash}"' if dash else "")+'/>')
    def note(self,s):
        self.rect(65,641,1070,33,BG,8)
        self.text(80,663,s,14,INK,500)
    def svg(self):return ("\n".join(self.a)+"\n</svg>\n").encode("utf-8")


def plot(svg,area,xmin,xmax,ymin,ymax,xticks,yticks,xlabel,ylabel,series,tag=None,highlight=None,legend=True):
    """area=(x,y,w,h); series list (label,color,[(x,y)...])."""
    x,y,w,h=area
    svg.rect(x,y,w,h,"#FFFFFF",12,GRID)
    svg.text(x+20,y+33,tag or "MODELO Y COMPARACIÓN",17,NAVY,700)
    left=x+72; right=x+w-25; top=y+(91 if legend else 65); bottom=y+h-64
    def xx(t):return left+(t-xmin)*(right-left)/(xmax-xmin)
    def yy(t):return bottom-(t-ymin)*(bottom-top)/(ymax-ymin)
    for val in yticks:
        cy=yy(val);svg.line(left,cy,right,cy,GRID,1)
        label=f"{val:.2f}" if abs(val)<2 and val!=int(val) else f"{val:g}"
        svg.text(left-12,cy+5,label,12,MUTED,400,anchor="end")
    for val in xticks:
        cx=xx(val);svg.line(cx,bottom,cx,bottom+5,MUTED,1.1)
        label=f"{val:.2f}" if abs(val)<2 and val!=int(val) else f"{val:g}"
        svg.text(cx,bottom+23,label,12,MUTED,400,anchor="middle")
    svg.line(left,bottom,right,bottom,NAVY,1.7);svg.line(left,top,left,bottom,NAVY,1.7)
    svg.text((left+right)/2,y+h-17,xlabel,14,MUTED,500,anchor="middle")
    svg.text(x+23,(top+bottom)/2,ylabel,13,MUTED,500,anchor="middle",rotate=-90)
    for lab,color,vals in series:
        points=[(xx(a),yy(b)) for a,b in vals]
        d=" ".join(("M" if j==0 else "L")+f"{a:.2f},{b:.2f}" for j,(a,b) in enumerate(points))
        svg.path(d,color,3.8)
    if highlight:
        for item in highlight:
            hx,hy,itemtext,color=item
            cx=xx(hx);cy=yy(hy)
            svg.circle(cx,cy,5.3,"#FFFFFF",stroke=color)
            svg.circle(cx,cy,2.4,color)
            svg.text(cx+8,cy-10,itemtext,12,color,700)
    if legend:
        # Draw compact legend below panel heading, above plot.
        lx=x+w-22
        for label,color,vals in reversed(series):
            length=min(280,len(label)*8.3+28);lx-=length
            svg.line(lx,y+60,lx+18,y+60,color,3)
            svg.text(lx+24,y+65,label,13,INK,600)
            lx-=14
    return xx,yy


def bar_panel(svg,area,title,bars,ymax,unit="puntos",baseline=0):
    x,y,w,h=area
    svg.rect(x,y,w,h,"#FFFFFF",12,GRID)
    svg.text(x+20,y+33,title,17,NAVY,700)
    x0=x+70; x1=x+w-25; y0=y+66; y1=y+h-64
    def scale(z):return y1-(z-baseline)*(y1-y0)/(ymax-baseline)
    for q in range(5):
        z=baseline+(ymax-baseline)*q/4
        cy=scale(z);svg.line(x0,cy,x1,cy,GRID,1)
        svg.text(x0-9,cy+4,f"{z:.2f}" if ymax<3 else f"{z:.1f}",12,MUTED,400,"end")
    count=len(bars);block=(x1-x0)/count
    for i,(label,val,color) in enumerate(bars):
        bw=min(67,block*.61);bx=x0+i*block+block/2-bw/2
        by=scale(val);svg.rect(bx,by,bw,y1-by,color,4)
        svg.text(bx+bw/2,by-10,f"{val:.3f}" if ymax<=3 else f"{val:.2f}",13,color,700,anchor="middle")
        svg.text(bx+bw/2,y1+26,label,12,INK,600,anchor="middle")
    svg.text((x0+x1)/2,y+h-17,unit,13,MUTED,500,anchor="middle")



# Every numerical record is calculated from an explicitly stated model in the original chapter.
# No table or diagram below represents an empirical dataset.
def rk4_niche(p0, E0, dt=.01, steps=1000):
    a,decay,b,c=.6,.3,1,.7
    def grad(z):
        p,E=z
        return (p*(1-p)*(b*E-c),a*p-decay*E)
    z=(p0,E0);out=[z]
    for _ in range(steps):
        k1=grad(z)
        k2=grad(tuple(z[i]+dt*k1[i]/2 for i in range(2)))
        k3=grad(tuple(z[i]+dt*k2[i]/2 for i in range(2)))
        k4=grad(tuple(z[i]+dt*k3[i] for i in range(2)))
        z=tuple(z[i]+dt*(k1[i]+2*k2[i]+2*k3[i]+k4[i])/6 for i in range(2))
        if not 0<=z[0]<=1 or z[1]<0: raise AssertionError('Niche ODE out of domain')
        out.append(z)
    return out


def make_data():
    D={}
    # C281 -- eight exactly specified genotype-fitness pairs; one-bit neighbors only.
    fitness={'000':1.0,'100':1.35,'010':1.10,'001':1.05,
             '110':1.20,'101':1.15,'011':1.40,'111':1.60}
    rows=[]
    for genotype in ['000','100','010','001','110','101','011','111']:
        neigh=[''.join(('1' if c=='0' else '0') if j==i else c for j,c in enumerate(genotype)) for i in range(3)]
        upgrades=[v for v in neigh if fitness[v]>fitness[genotype]]
        rows.append((genotype,fitness[genotype],len(upgrades),int(not upgrades),int(genotype=='100'),int(genotype=='111')))
    D[281]=(['genotype','fitness','improving_single_step_neighbors','is_local_maximum','is_local_trap_100','is_global_optimum_111'],rows)
    # C282 -- composition-dependent rates, selection field and total-rate derivative.
    rows=[]
    for i in range(101):
        p=i/100;ma=.6-.4*p;mb=.2+.2*p;adv=ma-mb
        dp=p*(1-p)*adv;mm=.2+.6*p-.6*p*p;dmm=(.6-1.2*p)*dp
        rows.append((p,ma,mb,adv,dp,mm,dmm))
    D[282]=(['frequency_A','growth_A','growth_B','growth_advantage_A','frequency_derivative','mean_growth','mean_growth_time_derivative'],rows)
    # C283 -- exact model solved using deterministic fourth-order Runge--Kutta (dt=.01).
    low=rk4_niche(.5,0);high=rk4_niche(.5,1.8)
    D[283]=(['time','frequency_A_E0_0','environment_E0_0','frequency_A_E0_1_8','environment_E0_1_8'],
            [(i/100,low[i][0],low[i][1],high[i][0],high[i][1]) for i in range(1001)])
    # C284 -- cultural source bias followed by binary transmission error.
    rows=[];x=.2;z=.2
    for t in range(31):
        q=1.5*x/(1+.5*x)
        qplain=z
        rows.append((t,x,q,z,qplain))
        x=.9*q+.05*(1-q)
        z=.9*z+.05*(1-z)
    D[284]=(['generation','transmitted_A_selection_1_5','selected_sources_A','transmitted_A_no_bias','sources_A_no_bias'],rows)
    # C285 -- declines within both groups but rises in weighted total.
    D[285]=(['class','initial_producer_fraction','end_within_fraction','offspring_weight','weighted_contribution','unweighted_mean_after_within'],[
        ('Group A',.8,.7,4,4*.7,.4),('Group B',.2,.1,1,.1,.4),('Population',.5,.58,5,2.9,.4)])
    # C286 -- a convex technology can favor specialists only if transfer is sufficiently efficient.
    rows=[]
    for i in range(111):
        alpha=.25+.025*i
        general=2*(.5**alpha);specialist=.7
        rows.append((alpha,general,specialist,specialist-general,int(specialist>general)))
    D[286]=(['production_exponent_alpha','generalist_output','specialist_output_eta_0_7','specialist_advantage','specialist_is_better'],rows)
    # C287 -- simultaneous updates (never use newly updated values within a generation).
    p,c,e=.2,.4,.3;rows=[]
    for t in range(101):
        rows.append((t,p,c,e,.5*c*e,.3*p+.2*e-.15))
        newp=p*(1+.5*c*e)/(1+.5*p*c*e)
        newc=c+.4*c*(1-c)*(.3*p+.2*e-.15)
        newe=.8*e+.2*c
        p,c,e=newp,newc,newe
        if not all(0<=q<=1 for q in (p,c,e)):raise AssertionError('C287 state out of domain')
    D[287]=(['generation','allele_frequency','practice_users','infrastructure','allele_fitness_bonus','culture_net_attraction'],rows)
    # C288 -- sum of new gains and retained previous achievement, not human measurement.
    D[288]=(['generation','retention_0_80','retention_0_95','equilibrium_0_80','equilibrium_0_95'],
            [(t,2*(1-.8**t)/.2,2*(1-.95**t)/.05,10,40) for t in range(41)])
    # C289 -- ten independent meetings, mixed type probability and 10% interface success.
    D[289]=(['frequency_A','mixed_meeting_probability','at_least_one_mixed_in_10','at_least_one_success_in_10_at_10pct'],[
            (i/100,2*(i/100)*(1-i/100),1-(1-2*(i/100)*(1-i/100))**10,
             1-(1-.1*2*(i/100)*(1-i/100))**10) for i in range(101)])
    # C290 -- model-defined hypothetical workload, no duration estimates or patient observations.
    arrivals=[20,24,28,30,30,30]
    queue=30;reference=30;qh=[30];rh=[30]
    for a in arrivals:
        queue=max(0,queue+a-24)
        reference=max(0,reference+18-24)
        qh.append(queue);rh.append(reference)
    rows=[]
    for i in range(101):
        f=i/100
        wk=i if i<=6 else ''
        arrival=(arrivals[i-1] if 1<=i<=6 else '')
        backlog=(qh[i] if i<=6 else '')
        stable=(rh[i] if i<=6 else '')
        rows.append((f,30*(1-.6*f),26-30*(1-.6*f),wk,arrival,backlog,stable))
    D[290]=(['eligible_share_f','weekly_workload_units','margin_with_capacity_26','week_if_le_6','new_requests_if_week','backlog_variable_requests','backlog_18_requests'],rows)
    return D


def charts(n,data):
    title,subtitle=NAMES[n];s=SVG(n,title,subtitle)
    if n==281:
        x,y,w,h=(55,176,1090,450)
        s.rect(x,y,w,h,'#FFFFFF',12,GRID)
        s.text(75,212,'APTITUDES EN OCHO GENOTIPOS',17,NAVY,700)
        labels=[r[0] for r in data];vals=[r[1] for r in data]
        L=134;R=1104;top=250;bottom=544
        for z in [0,.4,.8,1.2,1.6]:
            cy=bottom-(z/1.75)*(bottom-top)
            s.line(L,cy,R,cy,GRID);s.text(L-13,cy+5,f'{z:.1f}',13,MUTED,400,'end')
        for i,(lab,val) in enumerate(zip(labels,vals)):
            cx=L+(i+.5)*(R-L)/8;cy=bottom-(val/1.75)*(bottom-top)
            color=OCHRE if lab=='100' else TEAL if lab=='111' else PURPLE
            s.rect(cx-38,cy,76,bottom-cy,color,6)
            s.text(cx,cy-13,f'{val:.2f}',15,color,700,'middle')
            s.text(cx,bottom+27,lab,14,INK,700,'middle')
            if lab in ('100','111'):
                s.text(cx,bottom+53,'local' if lab=='100' else 'global',13,color,700,'middle')
        s.note('Desde 000, el ascenso más pronunciado termina en 100; el máximo global es 111.')
    elif n==282:
        plot(s,(55,176,539,450),0,1,.1,.7,[0,.25,.5,.75,1],[.1,.2,.3,.4,.5,.6,.7],
             'Frecuencia de A','Crecimiento',[
             ('mA',TEAL,[(r[0],r[1]) for r in data]),
             ('mB',OCHRE,[(r[0],r[2]) for r in data]),
             ('Media',PURPLE,[(r[0],r[5]) for r in data])],
             'CRECIMIENTO Y COMPOSICIÓN')
        plot(s,(606,176,539,450),0,1,-.02,.045,[0,.25,.5,.75,1],[-.02,-.01,0,.01,.02,.03,.04],
             'Frecuencia de A','Cambio por tiempo',[
             ('dp/dt',TEAL,[(r[0],r[4]) for r in data]),
             ('dmedia/dt',OCHRE,[(r[0],r[6]) for r in data])],
             'SELECCIÓN Y CAMBIO DEL PROMEDIO')
        s.note('p* = 2/3 es estable; entre p=0,5 y p*, A se expande aunque la media disminuya.')
    elif n==283:
        plot(s,(55,176,539,450),0,10,0,1.04,[0,2,4,6,8,10],[0,.25,.5,.75,1],
             'Tiempo','Frecuencia de A',[
             ('E inicial 0',OCHRE,[(r[0],r[1]) for r in data]),
             ('E inicial 1,8',TEAL,[(r[0],r[3]) for r in data])],
             'MISMA POBLACIÓN, DOS DESTINOS')
        plot(s,(606,176,539,450),0,10,0,2.1,[0,2,4,6,8,10],[0,.5,1,1.5,2],
             'Tiempo','Estado ambiental E',[
             ('E inicial 0',OCHRE,[(r[0],r[2]) for r in data]),
             ('E inicial 1,8',TEAL,[(r[0],r[4]) for r in data])],
             'PERSISTENCIA DEL CONTEXTO')
        s.note('RK4, paso 0,01; p(0)=0,5. La memoria E(0) invierte la dirección selectiva.')
    elif n==284:
        plot(s,(55,176,1090,450),0,30,0,.85,[0,5,10,15,20,25,30],[0,.2,.4,.6,.8],
             'Rondas de transmisión','Proporción esperada de A',[
             ('Selección de fuentes',TEAL,[(r[0],r[1]) for r in data]),
             ('Sin ventaja de elección',OCHRE,[(r[0],r[3]) for r in data])],
             'COPIA TRANSFORMADA TRAS ELEGIR UN MODELO',[(5,data[5][1],'x5 ≈ 0,57',TEAL)])
        s.note('Pesos 1,5:1; tasas A→B 0,10 y B→A 0,05. Modelo cultural, no genética.')
    elif n==285:
        x,y,w,h=(55,176,545,450);s.rect(x,y,w,h,'#FFFFFF',12,GRID)
        s.text(x+20,y+33,'COMPOSICIÓN ANTES Y DESPUÉS',16,NAVY,700)
        l=x+70;r=x+w-25;top=y+73;bot=y+h-64
        for z in [0,.2,.4,.6,.8,1]:
            cy=bot-z*(bot-top);s.line(l,cy,r,cy,GRID);s.text(l-12,cy+5,f'{z:.1f}',12,MUTED,400,'end')
        B=[('Grupo A',.8,.7),('Grupo B',.2,.1),('Global',.5,.58)]
        step=(r-l)/3
        for i,(lab,prev,nxt) in enumerate(B):
            cx=l+(i+.5)*step
            for dx,z,col in [(-23,prev,PURPLE),(23,nxt,TEAL)]:
                cy=bot-z*(bot-top);s.rect(cx+dx-16,cy,32,bot-cy,col,3)
                s.text(cx+dx,cy-11,f'{z:.2f}',12,col,700,'middle')
            s.text(cx,bot+26,lab,12,INK,600,'middle')
        s.line(x+29,y+52,x+47,y+52,PURPLE,4);s.text(x+53,y+57,'Inicial',12,INK,500)
        s.line(x+119,y+52,x+137,y+52,TEAL,4);s.text(x+143,y+57,'Final',12,INK,500)
        x,y,w,h=(606,176,539,450);s.rect(x,y,w,h,'#FFFFFF',12,GRID)
        s.text(x+20,y+33,'COMPONENTES DE PRICE',16,NAVY,700)
        l=x+78;r=x+w-25;top=y+74;bot=y+h-64;mi=-.13;ma=.22
        def yy(z):return bot-(z-mi)/(ma-mi)*(bot-top)
        for z in [-.1,0,.1,.2]:
            cy=yy(z);s.line(l,cy,r,cy,NAVY if z==0 else GRID,1.7 if z==0 else 1);s.text(l-12,cy+4,f'{z:+.2f}',12,MUTED,400,'end')
        for i,(lab,z,col) in enumerate([('Entre',.18,TEAL),('Dentro',-.1,OCHRE),('Neto',.08,PURPLE)]):
            cx=l+(i+.5)*(r-l)/3;cy=yy(z);zero=yy(0)
            s.rect(cx-35,min(cy,zero),70,abs(cy-zero),col,4)
            s.text(cx,cy-11 if z>=0 else cy+22,f'{z:+.2f}',14,col,700,'middle')
            s.text(cx,bot+29,lab,12,INK,600,'middle')
        s.note('Los productores disminuyen dentro de ambos grupos, pero la frecuencia global sube a 0,58.')
    elif n==286:
        plot(s,(55,176,539,450),.25,3,0,1.85,[.25,.75,1.25,1.75,2.25,2.75],[0,.4,.8,1.2,1.6],
             'Exponente productivo α','Producto limitante Q',[
             ('Generalistas',OCHRE,[(r[0],r[1]) for r in data]),
             ('Especialistas η=0,7',TEAL,[(r[0],r[2]) for r in data])],
             'VENTAJA CON RENDIMIENTOS CRECIENTES')
        rows=[(m,-.2+.6*m/(10+m)) for m in [i/10 for i in range(301)]]
        plot(s,(606,176,539,450),0,30,-.22,.25,[0,5,10,15,20,25,30],[-.2,-.1,0,.1,.2],
             'Densidad del socio M','Crecimiento de N escaso',[
             ('rN + beneficio',TEAL,rows)],'UMBRAL DE ESTABLECIMIENTO',[(5,0,'M crítico=5',TEAL)])
        s.note('La especialización necesita intercambio eficiente y un socio presente; Q no mide aptitud empírica.')
    elif n==287:
        plot(s,(55,176,1090,450),0,100,0,1.05,[0,20,40,60,80,100],[0,.2,.4,.6,.8,1],
             'Generaciones','Frecuencia / disponibilidad',[
             ('Alelo p',TEAL,[(r[0],r[1]) for r in data]),
             ('Práctica c',OCHRE,[(r[0],r[2]) for r in data]),
             ('Infraestructura e',PURPLE,[(r[0],r[3]) for r in data])],
             'COEVOLUCIÓN: ACTUALIZACIÓN SIMULTÁNEA')
        s.note('La práctica primero disminuye y después se extiende; el entorno conserva memoria de uso previo.')
    elif n==288:
        plot(s,(55,176,1090,450),0,40,0,42,[0,10,20,30,40],[0,10,20,30,40],
             'Etapas de transmisión','Desempeño didáctico acumulado',[
             ('Retención 0,80',OCHRE,[(r[0],r[1]) for r in data]),
             ('Retención 0,95',TEAL,[(r[0],r[2]) for r in data])],
             'GANANCIA CONSTANTE g=2')
        s.note('Equilibrios exactos: 10 y 40. No es una escala general de conocimiento cultural.')
    elif n==289:
        plot(s,(55,176,1090,450),0,1,0,1.07,[0,.25,.5,.75,1],[0,.2,.4,.6,.8,1],
             'Fracción que dispone de la variante A','Probabilidad acumulada',[
             ('Al menos 1 encuentro mixto',TEAL,[(r[0],r[2]) for r in data]),
             ('Al menos 1 combinación exitosa',OCHRE,[(r[0],r[3]) for r in data])],
             '10 ENCUENTROS, ÉXITO CONDICIONAL 10 %')
        s.note('Preservar variantes es distinto de conectarlas: encuentros y capacidad de integración se multiplican.')
    elif n==290:
        rows=[r for r in data if r[3]!='']
        plot(s,(55,176,539,450),0,6,0,53,[0,1,2,3,4,5,6],[0,10,20,30,40,50],
             'Semanas','Solicitudes pendientes',[
             ('Demanda creciente',OCHRE,[(r[3],r[5]) for r in rows]),
             ('18 por semana',TEAL,[(r[3],r[6]) for r in rows])],
             'LISTA INICIAL DE 30 SOLICITUDES')
        plot(s,(606,176,539,450),0,1,10,31,[0,.2,.4,.6,.8,1],[10,15,20,25,30],
             'Proporción elegible f','Unidades de trabajo por semana',[
             ('L(f)',PURPLE,[(r[0],r[1]) for r in data]),
             ('Capacidad 26',TEAL,[(0,26),(1,26)])],
             'EQUIVALENCIA DIDÁCTICA DE TRABAJO',[(2/9,26,'f=2/9',PURPLE)])
        s.note('Son solicitudes y equivalencias contables sintéticas, no duraciones clínicas ni criterios de derivación.')
    return s.svg()


def assert_near(a,b,name,tol=2e-5):
    if not isclose(a,b,rel_tol=tol,abs_tol=tol):raise AssertionError(f'{name}: {a} versus {b}')


def audit_model(d):
    f={str(row[0]):row[1] for row in d[281][1]}
    assert_near(f['100'],1.35,'C281 local peak');assert_near(f['111'],1.6,'C281 global peak')
    assert sorted(r[0] for r in d[281][1] if r[3])==['100','111']
    # Probability of reaching the global maximum under uniform choice of improving neighbors.
    from functools import lru_cache
    @lru_cache(None)
    def hit(g):
        neigh=[''.join(('1' if c=='0' else '0') if j==i else c for j,c in enumerate(g)) for i in range(3)]
        ups=[v for v in neigh if f[v]>f[g]]
        return (1 if g=='111' else 0) if not ups else sum(hit(v) for v in ups)/len(ups)
    assert_near(hit('000'),.5,'C281 uniform uphill search hits global peak')
    r=d[282][1][60]
    assert_near(r[3],.04,'C282 advantage at p=.6');assert_near(r[4],.0096,'C282 dpdt')
    assert_near(r[6],-.001152,'C282 declining mean despite selection')
    for idx,tgt1,tgt2 in [(200,.269,.877),(500,.114,.994),(1000,.012,.99998)]:
        row=d[283][1][idx]
        assert_near(row[1],tgt1,'C283 low-environment basin',tol=.015)
        assert_near(row[3],tgt2,'C283 high-environment basin',tol=.015)
    assert_near(d[284][1][1][1],.28181818181818186,'C284 frequency after one transmission')
    assert_near(d[284][1][5][1],.570,'C284 after five transmissions',tol=.004)
    assert_near(d[285][1][-1][2],.58,'C285 increasing global production frequency')
    assert_near(.18-.1,.08,'C285 Price total')
    assert_near(d[286][1][70][1],.5,'C286 alpha2 generalist output')
    assert_near(.6*5/(10+5)-.2,0,'C286 establishment threshold')
    assert_near(d[287][1][1][1],.2094861660079,'C287 initial genotype change')
    assert_near(d[287][1][1][2],.39712,'C287 initial culture change')
    assert_near(d[287][1][1][3],.32,'C287 initial environmental change')
    for i,reference in [(10,(.3284,.3963,.3833)),(50,(.9965,.919,.849)),(100,(1,1,1))]:
        for k,val in enumerate(reference,1):assert_near(d[287][1][i][k],val,'C287 coupling over time',tol=.001)
    assert_near(d[288][1][10][1],8.926258176,'C288 10 rounds f=.8')
    assert_near(d[288][1][10][2],16.050522449,'C288 10 rounds f=.95')
    assert_near(d[289][1][90][2],.862552397,'C289 at least one encounter with p=.9',tol=.002)
    assert_near(d[289][1][99][3],.019623657,'C289 interface success with p=.99',tol=.001)
    assert [row[5] for row in d[290][1][:7]]==[30,26,26,30,36,42,48]
    assert [row[6] for row in d[290][1][:7]]==[30,24,18,12,6,0,0]
    assert_near(d[290][1][25][1],25.5,'C290 25 percent eligible work')
    assert_near(30*(1-.6*(2/9)),26,'C290 workload threshold')
    for n,(header,rows) in d.items():
        assert len(header)>=4 and rows and all(len(r)==len(header) for r in rows),f'Wrong CSV shape C{n}'
    assert len(d)==10


def byte_outputs():
    data=make_data();audit_model(data);files={}
    for n in range(281,291):
        header,rows=data[n]
        files[f'repro/c281-290/c{n}.csv']=csv_bytes(header,rows)
        fig=charts(n,rows)
        parsed=ElementTree.fromstring(fig);ns='{http://www.w3.org/2000/svg}'
        if parsed.find(ns+'title') is None or parsed.find(ns+'desc') is None:raise AssertionError(f'C{n} accessibility')
        files[f'figures/complexa-c{n}-fig01.svg']=fig
    return files


def make_archive(files,readme,script):
    entries=dict(files)
    entries['repro/c281-290/README.md']=readme
    entries['repro/c281-290/make_figures.py']=script
    stream=BytesIO()
    with zipfile.ZipFile(stream,'w',compression=zipfile.ZIP_STORED,allowZip64=False) as archive:
        for name in sorted(entries):
            info=zipfile.ZipInfo(name,date_time=(2026,10,26,12,0,0))
            info.compress_type=zipfile.ZIP_STORED;info.create_system=3;info.external_attr=0o100644<<16
            archive.writestr(info,entries[name])
    return stream.getvalue()


def outputs():
    d=byte_outputs();here=Path(__file__)
    readme=here.with_name('README.md').read_bytes()
    script=here.read_bytes()
    archive=make_archive(d,readme,script)
    d['zip-data/complexa-c281-c290-reproducible.zip.b64']=base64.b64encode(archive)+b'\n'
    return d,archive


def run():
    ap=argparse.ArgumentParser(description=__doc__)
    group=ap.add_mutually_exclusive_group()
    group.add_argument('--out',type=Path,help='Directory content/supplementary to write')
    group.add_argument('--check',type=Path,help='Directory to verify generated bytes')
    args=ap.parse_args();dest=args.out or args.check or ROOT
    d,archive=outputs()
    for rel,blob in sorted(d.items()):
        target=dest/rel
        if args.check:
            if not target.is_file() or target.read_bytes()!=blob:raise SystemExit(f'FAIL reproducibility: {rel}')
        else:
            target.parent.mkdir(parents=True,exist_ok=True);target.write_bytes(blob)
    if len(d)!=21:raise AssertionError('Expected 10 CSV, 10 SVG and a base64 ZIP')
    with zipfile.ZipFile(BytesIO(archive)) as z:
        if len(z.namelist())!=22 or z.testzip():raise AssertionError('ZIP integrity failed')
    print(f"{'CHECK' if args.check else 'GENERATED'}: 10 SVG, 10 CSV, deterministic ZIP SHA-256 {sha256(archive).hexdigest()}")


if __name__=='__main__':run()