"""Repeated chemical-state selection: full physical source and finite-mission design."""
# EDITABLE INPUTS. These are theorem-scale witnesses, not calibrated biological cells.
NEWBORN_SIZE = 65536 * 10**18
RETAINED_CELLS = 4 * 10**9
CYCLES = 10
GROWTH_COUPLING = '1/100000000000'
BATCH_SERVICE_ALLOWANCE = '1/10000'
RECOVERY_SERVICE_ALLOWANCE = '1/10000'
TRANSFER_FAILURE_BUDGET = '1/250'
DESIGN_POPULATION_CAP = 10**16
MINORITY_RESERVE = 1
REFERENCE_MOLAR = '1/1000'     # Optional concentration convention, not a fit.
TIME_UNIT_SECONDS = '1'       # Optional time convention, not a measured rate.
# Small physical diagnostic: deliberately outside the proved source regime.
SMALL_SIZE = 10
SMALL_CELLS = 4
SMALL_CYCLES = 3
SMALL_GROWTH = '2'
SMALL_RECOVERY = '.05'
SMALL_EVENT_BUDGET = 100000
RANDOM_SEED = 77092026

from pathlib import Path
from dataclasses import asdict
from fractions import Fraction as Q
import argparse,csv,hashlib,json,platform
import numpy as np
import scipy
from resident import Cell,ChemicalRegions,A
from population import ChemicalModel,PopulationSource,SerialProtocol,UniformTransfer
from selection import MissionCertificate,WeightedTransfer,AbsorbingAudit,population_design,event_ceiling,RHO,RHOS,G,LOG2
import paper_checks

MANUSCRIPT_SHA256='89018b0f5f84c348611b59ed44e378c4f360972dbf601feabc5f1c0e9f0ee912'

def main():
    parser=argparse.ArgumentParser(description=__doc__);parser.add_argument('--output',type=Path,default=Path('outputs'));out=parser.parse_args().output;out.mkdir(parents=True,exist_ok=True)
    def dump(name,value):(out/name).write_text(json.dumps(value,indent=2,default=str)+'\n')
    def table(name,rows):
        with (out/name).open('w',newline='') as f:w=csv.DictWriter(f,fieldnames=list(rows[0]));w.writeheader();w.writerows(rows)
    cert=MissionCertificate(NEWBORN_SIZE,RETAINED_CELLS,CYCLES,Q(GROWTH_COUPLING),Q(BATCH_SERVICE_ALLOWANCE),Q(RECOVERY_SERVICE_ALLOWANCE))
    dump('configured_mission.json',cert.evaluate());dump('resource_accounts.json',cert.resources());dump('chemical_terms.json',{k:str(v) for k,v in cert.chemistry().items()})
    dump('population_design.json',population_design(CYCLES,TRANSFER_FAILURE_BUDGET,max_M=DESIGN_POPULATION_CAP))
    dump('necessary_event_ceiling.json',event_ceiling(RETAINED_CELLS,CYCLES,MINORITY_RESERVE))
    comparisons=[]
    for name,M,refined,s in [('original',10**13,False,Q(1,10000)),('refined',4*10**9,True,Q(1,10000)),('refined_finer_service',4*10**9,True,Q(1,10**6))]:
        c=MissionCertificate(M=M,refined=refined,batch_allowance=s,recovery_allowance=s);row=c.evaluate();row.update(case=name,binary_chemical_replay=str(c.binary_replay()),resources=c.resources());comparisons.append(row)
    dump('canonical_comparison.json',comparisons)
    dump('source_arithmetic_replay.json',dict(results=paper_checks.OUT,checks=paper_checks.CHECKS,scope='Exact enclosures and explicitly labeled numerical diagnostics from the source. Neither the process theorem nor Lean is reproved.'))
    rows=[]
    for k in range(1,19):
        c=MissionCertificate(N=NEWBORN_SIZE,M=RETAINED_CELLS,K=k,gamma=Q(GROWTH_COUPLING),batch_allowance=Q(BATCH_SERVICE_ALLOWANCE),recovery_allowance=Q(RECOVERY_SERVICE_ALLOWANCE));r=c.evaluate()
        rows.append(dict(cycle=k,success_lower=r['joint_success_decimal'],count_logodds_lower=float(k*G[0]-LOG2[1]),high_fraction_lower=1/(1+np.exp(-float(k*G[0]-LOG2[1]))),minority_floor=r['minority_count_floor'],transfer_failure_upper=float(sum(c.transfers())),original_transfer_failure_upper=float(Q(80000,RETAINED_CELLS)*sum(RHO**-j for j in range(k)))))
    table('horizon_sweep.csv',rows)
    exact=[];cases=[]
    for name,sizes in [('equal',[10]*8),('mixed',[10,10,19,19]*2)]:
        for tolerance in [Q(1,50),Q(1,2)]:
            result=WeightedTransfer.enumerate(sizes,['H']*4+['L']*4,4,tolerance)
            exact.append(dict(case=name,tolerance=str(tolerance),failure=str(result['failure']),moments=result['moments'],scope='Conditional endpoint population; not a probability of reaching that population.'))
            if tolerance==Q(1,50):cases.append((name,result))
            table(f'transfer_{name}_{tolerance.denominator}.csv',[{k:str(v) for k,v in r.items()} for r in result['rows']])
    dump('weighted_transfer_comparison.json',exact)
    table('conditional_minority_loss.csv',[dict(endpoint_count=400,retained=100,minority=c,loss_exact=str(WeightedTransfer.minority_loss(400,100,c)),loss=float(WeightedTransfer.minority_loss(400,100,c))) for c in [1,4,20]])
    regions=ChemicalRegions();chem=ChemicalModel();dump('chemical_identities.json',chem.certificates())
    prepared=[regions.newborn(t,NEWBORN_SIZE) for t in ('H','L')]
    dump('ready_representatives.json',[dict(cell=asdict(c),energy_interval=regions.energy(c).json(),readout=c.readout,scope='Two representatives, not allocation of the full theorem-scale population.') for c in prepared])
    if SMALL_CELLS<2 or SMALL_CELLS%2:raise ValueError('Use an even diagnostic population.')
    small=[]
    for tag in ['H']*(SMALL_CELLS//2)+['L']*(SMALL_CELLS//2):
        counts=tuple(int(float((x.lo+x.hi)/2)*SMALL_SIZE) for x in regions.centers[tag]);small.append(Cell(tag,SMALL_SIZE,counts))
    physical=SerialProtocol(PopulationSource(SMALL_SIZE,SMALL_GROWTH),SMALL_CELLS,float(SMALL_RECOVERY)).run(small,SMALL_CYCLES,RANDOM_SEED,SMALL_EVENT_BUDGET)
    dump('small_physical_history.json',physical)
    audit=AbsorbingAudit();ready=all(regions.admits(c,A,closed=True) for c in small);audit.observe('initial_readiness',ready,'Exact energy comparison for every diagnostic founder.');audit.observe('illustrative_later_pass',True,'A failed or unresolved history stays absorbed.')
    dump('small_marked_audit.json',dict(**asdict(audit),scope='Readiness mark demonstration only. Physical history is saved separately; the full paper marked law additionally tests intermediate energies, division, growth, transfer, deadline and services. Do not classify physical completion as theorem success.'))
    table('physical_censuses.csv',[dict(cycle=h['cycle']+1,start_size=h['start_size'],selected_H=sum(c['size'] for c in h['selected'] if c['tag']=='H'),selected_L=sum(c['size'] for c in h['selected'] if c['tag']=='L'),residual_precursor=h['residual_discarded'],discarded_size=h['discarded_size']) for h in physical['history']] or [dict(cycle=0,start_size=sum(c.size for c in small),selected_H=0,selected_L=0,residual_precursor=0,discarded_size=0)])
    NA=Q(602214076000000000000000);concentration=Q(REFERENCE_MOLAR);timeunit=Q(TIME_UNIT_SECONDS)
    if min(concentration,timeunit)<=0:raise ValueError('Positive reference units required.')
    dump('illustrative_units.json',dict(newborn_litres=float(Q(NEWBORN_SIZE)/NA/concentration),elapsed_years=float(CYCLES*(cert.T+5376)*timeunit/Q(31557600)),retained_cells=RETAINED_CELLS,scope='N is a size/copy convention. This concentration/time interpretation is illustrative. M counts whole cells and is never divided by Avogadro.'))
    plot(out,rows,cases)
    print('Full resident source: 13 chemical labels plus growth; complementary division; intact uniform transfer; inherited size and counts.',flush=True)
    print('Configured joint success lower:',cert.evaluate()['joint_success_decimal'],'; minority floor:',cert.evaluate()['minority_count_floor'],flush=True)
    print('Physical small run completed:',physical['complete'],'; readiness-mark status:',audit.status,flush=True)
    print('Source checks:',len(paper_checks.CHECKS),'; bounds are sufficient and imported process assumptions remain explicit. Lean not rerun.',flush=True)
    here=Path(__file__).parent;digest=lambda p:hashlib.sha256(p.read_bytes()).hexdigest()
    dump('run_metadata.json',dict(manuscript_sha256=MANUSCRIPT_SHA256,source_sha256=digest(Path(__file__)),python=platform.python_version(),numpy=np.__version__,scipy=scipy.__version__,module_sha256={p.name:digest(p) for p in sorted(here.glob('*.py')) if p.name not in ['example.py','test_example.py']},input_sha256={p.name:digest(p) for p in sorted(here.glob('*.json')) if p.name!='release.json'},output_sha256={p.name:digest(p) for p in sorted(out.iterdir()) if p.is_file() and p.name!='run_metadata.json'}))

def plot(out,rows,cases):
    import matplotlib
    matplotlib.use('Agg')
    import matplotlib.pyplot as plt
    fig,axs=plt.subplots(1,2,figsize=(11,4),layout='constrained')
    axs[0].plot([r['cycle'] for r in rows],[r['high_fraction_lower'] for r in rows],label='High fraction on success')
    axs[0].plot([r['cycle'] for r in rows],[r['success_lower'] for r in rows],label='Joint mission confidence')
    axs[0].set(xlabel='Number of cycles',ylabel='Lower bound',ylim=(-.03,1.04),title='Composition and joint-success bounds over\nrepeated cycles');axs[0].legend(fontsize=8)
    axs[1].semilogy([r['cycle'] for r in rows],[r['transfer_failure_upper'] for r in rows],label='Minority growth + exponential tail');axs[1].semilogy([r['cycle'] for r in rows],[r['original_transfer_failure_upper'] for r in rows],label='Original bound at same M');axs[1].axhline(1,color='gray',ls=':');axs[1].set(xlabel='Number of cycles',ylabel='Sufficient transfer error (uncapped)',title='Both comparisons retain four billion cells');axs[1].legend(fontsize=8)
    for a in axs:a.grid(alpha=.2)
    fig.savefig(out/'mission_horizon.png',dpi=180);fig.savefig(out/'mission_horizon.svg');plt.close(fig)
    fig,axs=plt.subplots(1,2,figsize=(11,4),layout='constrained')
    for name,result in cases:
        xs=sorted(set(float(r['H']) for r in result['rows']));ps=[sum(float(r['H'])==x for r in result['rows'])/len(result['rows']) for x in xs]
        axs[0].plot(xs,ps,'o-',label=name+' sizes')
    axs[0].set(xlabel='Retained high-type size',ylabel='Exact subset probability',title='Retained high-type size under two cell-size\ndistributions');axs[0].legend(fontsize=8)
    cs=np.arange(1,21);axs[1].semilogy(cs,[float(WeightedTransfer.minority_loss(400,100,int(c))) for c in cs],'o-');axs[1].set(xlabel='Minority cells among 400 endpoint cells',ylabel='Probability all are lost',title='Uniform transfer retains 100 whole cells')
    for a in axs:a.grid(alpha=.2)
    fig.savefig(out/'intact_cell_transfer.png',dpi=180);fig.savefig(out/'intact_cell_transfer.svg');plt.close(fig)

if __name__=='__main__':main()