"""Finite-mission arithmetic and intact-cell sampling, independent of SSA solvers."""
from dataclasses import dataclass,field
from fractions import Fraction as Q
from functools import lru_cache
from math import comb,ceil
from itertools import combinations

D0=512*10**18
N0=140*10**18
EPS=Q(1,50)
RHO=Q(49,204)
RHOS=Q(49,136)

def log_ratio(x,terms=180):
    """Exact enclosure of log(x) for x>0, using binary range reduction."""
    x=Q(x)
    if x<=0:raise ValueError('Positive logarithm argument required.')
    k=0
    while x>=2:x/=2;k+=1
    while x<1:x*=2;k-=1
    def core(z):
        s=sum((z**(2*j+1)/Q(2*j+1) for j in range(terms)),Q(0))*2
        return s,s+2*z**(2*terms+1)/((2*terms+1)*(1-z*z))
    a,b=core((x-1)/(x+1));l,h=core(Q(1,3))
    return a+min(k*l,k*h),b+max(k*l,k*h)

LOG2=log_ratio(2)
LOG5149=log_ratio(Q(51,49))
G=(Q(6,5)*LOG2[0]-Q(19,500)-LOG5149[1],Q(6,5)*LOG2[1]-Q(19,500)-LOG5149[0])

def upper_grid(x,digits=30):return Q(ceil(x*10**digits),10**digits)
def lower_grid(x,digits=30):return Q(x*10**digits//1,10**digits)

@lru_cache(maxsize=4096)
def exp_negative_upper(x):
    """Outward rational bound: e^-x <= [sum_0^79 (min(x,512)/16)^j/j!]^-16.
    Capping x weakens the bound; it never silently drops a tiny term.
    """
    x=Q(x)
    if x<0:raise ValueError('Nonnegative exponent required.')
    y=min(x,Q(512))/16;t=Q(1);s=t
    for j in range(1,80):t*=y/j;s+=t
    return s**-16

def tail(mu):
    mu=Q(mu)
    if mu<=0:raise ValueError('Positive conditional mean required.')
    return upper_grid(min(1/(EPS*EPS*mu),exp_negative_upper(EPS*EPS*mu/2)+exp_negative_upper(EPS*EPS*mu/(2+EPS))))

@dataclass(frozen=True)
class MissionCertificate:
    N:int=65536*10**18
    M:int=4*10**9
    K:int=10
    gamma:Q=Q(1,10**11)
    batch_allowance:Q=Q(1,10000)
    recovery_allowance:Q=Q(1,10000)
    refined:bool=True
    def __post_init__(self):
        for key in ('gamma','batch_allowance','recovery_allowance'):object.__setattr__(self,key,Q(getattr(self,key)))
        if type(self.N)is not int or self.N<N0 or type(self.M)is not int or self.M<2 or self.M%2 or type(self.K)is not int or self.K<1:raise ValueError('Admitted molecular scale, even M>=2, and K>=1 required.')
        if not 0<self.gamma<=Q(1,10**11) or not 0<self.batch_allowance<1 or not 0<self.recovery_allowance<1:raise ValueError('The certificate requires 0<gamma<=1e-11 and service allowances in (0,1).')
    @property
    def T(self):return 8/self.gamma
    @property
    def rho(self):return RHOS if self.refined else RHO
    def chemistry(self):
        u=Q(self.N,D0);M=self.M;N=self.N;T=self.T
        # coefficient, positive exponential argument. These are full formulas,
        # not the paper's deliberately weakened binary display envelope.
        terms={'outer_initial':(M,7*u),'daughter_initial':(14*M,4*u),'outer_drift':(M*T*u/42,15*u/2),'division_initial':(M,u),'division_drift':(M*T*u/42,3*u/2),'partition':(56*M,Q(N,35*10**12)),'deadline':(1,Q(N,2500)),'size_gain':(1,Q(19*N,500000)),'recovery_initial':(M,u),'recovery_return':(2*M,u/2),'recovery_exit':(M,8*u),'recovery_drift':(16*M*u,31*u/2)}
        if self.refined:terms['minority_growth']=(1,Q(19*N,500000))
        return {k:upper_grid(a*exp_negative_upper(x)) for k,(a,x) in terms.items()}
    def transfers(self):
        if not self.refined:return [Q(80000,self.M)*RHO**-j for j in range(self.K)]
        return [tail(Q(3*self.M,32)*RHOS**j)+tail(Q(3*self.M,32)) for j in range(self.K)]
    def binary_replay(self):
        if self.N!=65536*10**18:raise ValueError('This display replay requires u=128.')
        return self.K*(self.M*(74+16*128+self.T*128/21)+2+int(self.refined))/2**128+Q(2*self.K*self.M,2**64)
    def resources(self):
        N,M,K=self.N,self.M,self.K;NM=N*M;qB=280000*NM;qR=34000*N
        JB=ceil(self.T*qB/self.batch_allowance)+1;JR=ceil(M*5376*qR/self.recovery_allowance)+1
        return dict(operating_precursor_upper=(8*K-4)*NM,growth_consumption_upper=(6*K-3)*NM,residual_discard_upper=(2*K-1)*NM,discarded_cell_size_upper=(7*K-4)*NM,terminal_refill_upper=8*NM,total_with_terminal_refill_upper=(8*K+4)*NM,elapsed_time=str(K*(self.T+5376)),batch_clock=qB,recovery_clock_per_cell=qR,JB=JB,JR=JR,gross_each_resident_reservoir_upper=K*(JB+M*JR),actual_batch_allowance=str(self.T*qB/JB),actual_recovery_allowance=str(Q(M*5376*qR,JR)))
    def evaluate(self):
        chem=self.K*sum(self.chemistry().values());transfer=sum(self.transfers());service=self.K*(self.batch_allowance+self.recovery_allowance);failure=upper_grid(chem+transfer+service)
        return dict(N=self.N,M=self.M,K=self.K,refined=self.refined,chemical_failure_upper=str(chem),transfer_failure_upper=str(upper_grid(transfer)),service_failure_upper=str(service),failure_upper=str(failure),joint_success_lower=str(max(Q(0),1-failure)),joint_success_decimal=float(max(Q(0),1-failure)),last_census_size_logodds_lower=str(lower_grid(self.K*G[0])),last_census_count_logodds_lower=str(lower_grid(self.K*G[0]-LOG2[1])),minority_count_floor=ceil(Q(self.M,4)*self.rho**self.K),scope='Uniform bounds apply jointly to all K marked censuses from balanced ready newborns. The stochastic source theorem is imported; no Lean compilation is performed.')

def population_design(K,delta,refined=True,max_M=10**16):
    """Smallest even M satisfying this monotone, outward transfer bound only.
    A budget exhaustion is unknown, not evidence that no design exists.
    """
    if type(K)is not int or K<1 or not 0<Q(delta)<1 or type(max_M)is not int or max_M<2:raise ValueError('Invalid design inputs.')
    def bound(M):return sum(MissionCertificate(M=M,K=K,refined=refined).transfers())
    hi=max_M//2
    if bound(2*hi)>Q(delta):return dict(status='unknown_population_cap',max_M=max_M)
    lo=0
    while hi-lo>1:
        mid=(lo+hi)//2
        if bound(2*mid)<=Q(delta):hi=mid
        else:lo=mid
    M=2*hi
    return dict(status='certified_transfer_design',M=M,bound=str(upper_grid(bound(M))),previous_even_bound=str(upper_grid(bound(M-2))) if M>2 else None,scope='Minimal for this implemented sufficient transfer bound, not minimal for the physical source. Chemical and service terms must still be added.')

def event_ceiling(M,K,reserve=1):
    if type(M)is not int or type(K)is not int or type(reserve)is not int or M<2 or K<1 or not 1<=reserve<M:raise ValueError('Invalid retained population, horizon or minority reserve.')
    l,h=log_ratio(Q(2*(M-reserve),reserve));a,b=K*G[0],K*G[1]
    return dict(status='impossible_target_event' if a>=h else 'not_excluded' if b<l else 'unresolved',ceiling_log_interval=[str(lower_grid(l)),str(upper_grid(h))],scope='Necessary ceiling for the specified simultaneous gain and minority-reserve event, not a theorem of physical extinction.')

class WeightedTransfer:
    @staticmethod
    def moments(weights,M):
        w=tuple(map(Q,weights));n=len(w)
        if not n or type(M)is not int or not 0<=M<=n:raise ValueError('Nonempty finite population and admissible sample count required.')
        mean=Q(M,n)*sum(w)
        variance=Q(M*(n-M),n*(n-1))*(sum(x*x for x in w)-sum(w)**2/n) if n>1 else Q(0)
        return mean,variance
    @staticmethod
    def minority_loss(n,M,c):
        if any(type(x)is not int for x in (n,M,c)) or n<0 or not 0<=M<=n or not 0<=c<=n:raise ValueError('Invalid hypergeometric population.')
        return Q(comb(n-c,M),comb(n,M)) if n-c>=M else Q(0)
    @staticmethod
    def enumerate(sizes,tags,M,epsilon):
        sizes=tuple(map(Q,sizes));tags=tuple(tags);n=len(sizes);eps=Q(epsilon)
        if len(tags)!=n or not n or min(sizes)<=0 or set(tags)!={'H','L'} or not 0<=M<=n or not 0<eps<1 or comb(n,M)>100000:raise ValueError('Invalid population or enumeration budget exceeded.')
        weights={tag:tuple(s if t==tag else Q(0) for s,t in zip(sizes,tags)) for tag in ('H','L')}
        moments={tag:WeightedTransfer.moments(w,M) for tag,w in weights.items()};rows=[]
        for sample in combinations(range(n),M):
            values={tag:sum(w[i] for i in sample) for tag,w in weights.items()}
            rows.append(dict(sample=sample,**values,failed=any(abs(values[t]-moments[t][0])>eps*moments[t][0] for t in weights)))
        return dict(failure=Q(sum(r['failed'] for r in rows),len(rows)),moments=moments,rows=rows)

@dataclass
class AbsorbingAudit:
    """Record supplied mark decisions without retrying or dropping failed histories.
    This transport is reusable; it does not itself verify every source-theorem mark.
    """
    status:str='active'
    records:list=field(default_factory=list)
    def observe(self,name,decision,evidence=None):
        if self.status!='active':return self.status
        if decision not in (True,False,None):raise ValueError('A mark must pass, fail or remain unresolved.')
        self.records.append(dict(mark=name,decision=decision,evidence=evidence))
        if decision is False:self.status='failed'
        elif decision is None:self.status='unresolved'
        return self.status
