#!/usr/bin/env python3
"""Exact static bitstream examples: RRR classes and partitioned Elias--Fano.
Standard library only; stdout only. Stored bytes include headers and directories.
No constant-time rank/select machine-word implementation is claimed for Python.
"""
from math import comb
from bisect import bisect_left
import json

def need(ok, message):
    if not ok: raise ValueError(message)
def integer(x,lo,hi=None):
    need(type(x) is int and x>=lo and (hi is None or x<hi),'integer range')
    return x

def width(states):
    integer(states,1)
    return (states-1).bit_length()
class Writer:
    def __init__(self):self.data=bytearray();self.n=0
    def put(self,x,w):
        integer(w,0);integer(x,0,1<<w)
        for j in range(w-1,-1,-1):
            if self.n%8==0:self.data.append(0)
            self.data[-1]|=((x>>j)&1)<<(7-self.n%8);self.n+=1
    def gamma(self,x):
        integer(x,1);w=x.bit_length();self.put(0,w-1);self.put(x,w)
class Reader:
    def __init__(self,data):
        need(type(data) is bytes,'bytes required');self.data=data;self.pos=0
    def at(self,pos,w):
        need(0<=pos and 0<=w and pos+w<=8*len(self.data),'truncated field')
        x=0
        for j in range(pos,pos+w):x=(x<<1)|((self.data[j//8]>>(7-j%8))&1)
        return x
    def get(self,w):x=self.at(self.pos,w);self.pos+=w;return x
    def gamma(self):
        z=0
        while self.get(1)==0:z+=1
        return (1<<z)|self.get(z)
    def finish(self,end):
        need((end+7)//8==len(self.data),'extra bytes')
        need(self.at(end,8*len(self.data)-end)==0,'nonzero byte padding')

def class_rank(bits):
    need(all(type(x) is int and x in (0,1) for x in bits),'bits required')
    positions=[i for i,x in enumerate(bits) if x]
    return len(positions),sum(comb(p,j+1) for j,p in enumerate(positions))
def class_unrank(s,c,z):
    integer(s,0);integer(c,0,s+1);integer(z,0,comb(s,c))
    bits=[0]*s;p=s-1
    for j in range(c,0,-1):
        while comb(p,j)>z:p-=1
        bits[p]=1;z-=comb(p,j);p-=1
    need(z==0,'rank remainder')
    return bits

def pack_rrr(bits,b=8,t=2):
    bits=tuple(bits);integer(b,1);integer(t,1)
    need(all(type(x) is int and x in (0,1) for x in bits),'bits required')
    n=len(bits);q=(n+b-1)//b;gw=n.bit_length();lw=(b*t).bit_length();cw=b.bit_length()
    classes=[];offsets=[];globals_=[];locals_=[];ones=payload=base1=basep=0
    for j in range(q):
        s=min(b,n-j*b);c,z=class_rank(bits[j*b:j*b+s]);w=width(comb(s,c))
        if j%t==0:base1,basep=ones,payload;globals_.append((ones,payload))
        locals_.append((ones-base1,payload-basep));classes.append(c);offsets.append((z,w))
        ones+=c;payload+=w
    out=Writer();out.gamma(n+1);out.gamma(b);out.gamma(t);out.put(ones,gw);out.put(payload,gw)
    for c in classes:out.put(c,cw)
    for x,y in globals_:out.put(x,gw);out.put(y,gw)
    for x,y in locals_:out.put(x,lw);out.put(y,lw)
    for z,w in offsets:out.put(z,w)
    return bytes(out.data)

class RRR:
    def __init__(self,data):
        r=Reader(data);self.r=r;self.n=r.gamma()-1;self.b=r.gamma();self.t=r.gamma()
        n,b,t=self.n,self.b,self.t;self.gw=n.bit_length();self.lw=(b*t).bit_length();self.cw=b.bit_length()
        self.ones=r.get(self.gw);self.payload=r.get(self.gw);self.q=(n+b-1)//b;ns=(self.q+t-1)//t
        self.classes=r.pos;self.globals=self.classes+self.q*self.cw
        self.locals=self.globals+ns*2*self.gw;self.offsets=self.locals+self.q*2*self.lw
        self.bit_length=self.offsets+self.payload;r.finish(self.bit_length)
        ones=payload=base1=basep=0
        for j in range(self.q):
            s=min(b,n-j*b);c=self.cls(j);need(c<=s,'invalid class');w=width(comb(s,c))
            if j%t==0:base1,basep=ones,payload
            gp,lp=self.parts(j)
            need(gp==(base1,basep) and lp==(ones-base1,payload-basep),'directory mismatch')
            need(r.at(self.offsets+payload,w)<comb(s,c),'unused class code')
            ones+=c;payload+=w
        need((ones,payload)==(self.ones,self.payload),'header totals')
    def cls(self,j):return self.r.at(self.classes+j*self.cw,self.cw)
    def parts(self,j):
        k=j//self.t;g=self.globals+2*k*self.gw;l=self.locals+2*j*self.lw;r=self.r
        return (r.at(g,self.gw),r.at(g+self.gw,self.gw)),(r.at(l,self.lw),r.at(l+self.lw,self.lw))
    def directory(self,j):
        gp,lp=self.parts(j)
        return gp[0]+lp[0],gp[1]+lp[1]
    def block(self,j):
        integer(j,0,self.q);s=min(self.b,self.n-j*self.b);c=self.cls(j);_,off=self.directory(j)
        return class_unrank(s,c,self.r.at(self.offsets+off,width(comb(s,c))))
    def rank(self,i):
        integer(i,0,self.n+1)
        if i==self.n:return self.ones
        j,h=divmod(i,self.b);prefix,_=self.directory(j)
        return prefix+sum(self.block(j)[:h])
    def access(self,i):integer(i,0,self.n);j,h=divmod(i,self.b);return self.block(j)[h]
    def select(self,k):
        integer(k,1)
        if k>self.ones:return None
        lo,hi=0,self.n
        while lo<hi:
            mid=(lo+hi)//2
            if self.rank(mid+1)<k:lo=mid+1
            else:hi=mid
        return lo
    def decode(self):return [x for j in range(self.q) for x in self.block(j)]

def ef_params(k,u):
    integer(k,1);integer(u,k)
    ell=(u//k).bit_length()-1;h=(u+(1<<ell)-1)>>ell
    return ell,h,k*ell+k+h

def choice(k,u):
    if k==u:return 0,0
    ef=ef_params(k,u)[2]
    return min((u,1),(ef,2)) # payload length, then deterministic mode

def pef_widths(n,u):
    W=(u-1).bit_length();C=n.bit_length();P=(n*(W+3)).bit_length()
    return W,C,P,2+W+C+P

def validate_sequence(xs,u):
    integer(u,1);xs=tuple(xs);last=-1
    for x in xs:integer(x,0,u);need(x>last,'strictly increasing');last=x
    return xs

def optimal_partition(xs,u):
    xs=validate_sequence(xs,u);n=len(xs);F=pef_widths(n,u)[3];dp=[0]+[None]*n;prev=[None]*(n+1)
    for j in range(1,n+1):
        for i in range(j):
            base=0 if i==0 else xs[i-1]+1;local=xs[j-1]-base+1
            cost=dp[i]+F+choice(j-i,local)[0]
            if dp[j] is None or cost<dp[j]:dp[j],prev[j]=cost,i
    cuts=[];j=n
    while j:cuts.append(j);j=prev[j]
    return list(reversed(cuts)),dp

def block_payload(values,base,mode,out):
    k=len(values);u=values[-1]-base+1;ys=[x-base for x in values]
    if mode==0:return
    if mode==1:
        cursor=0
        for y in ys:
            while cursor<y:out.put(0,1);cursor+=1
            out.put(1,1);cursor+=1
        return
    ell,h,_=ef_params(k,u)
    for y in ys:out.put(y% (1<<ell),ell)
    cursor=0
    for i,y in enumerate(ys):
        target=(y>>ell)+i
        while cursor<target:out.put(0,1);cursor+=1
        out.put(1,1);cursor+=1
    while cursor<k+h:out.put(0,1);cursor+=1

def pack_pef(xs,u,cuts=None):
    xs=validate_sequence(xs,u);n=len(xs);W,C,P,F=pef_widths(n,u)
    if cuts is None:cuts,_=optimal_partition(xs,u)
    else:cuts=list(cuts)
    prev=0
    for end in cuts:integer(end,prev+1,n+1);prev=end
    need((prev==n and (n>0 or not cuts)),'partition coverage')
    rows=[];payload=Writer();start=0
    for end in cuts:
        base=0 if start==0 else xs[start-1]+1;size=end-start;local=xs[end-1]-base+1;cost,mode=choice(size,local)
        rows.append((mode,xs[end-1],end,payload.n))
        before=payload.n;block_payload(xs[start:end],base,mode,payload);need(payload.n-before==cost,'payload size');start=end
    out=Writer();out.gamma(n+1);out.gamma(u);out.put(len(cuts),C);out.put(payload.n,P)
    for mode,last,end,off in rows:out.put(mode,2);out.put(last,W);out.put(end,C);out.put(off,P)
    for i in range(payload.n):out.put((payload.data[i//8]>>(7-i%8))&1,1)
    return bytes(out.data)

class PartitionedEF:
    def __init__(self,data):
        r=Reader(data);self.r=r;self.n=r.gamma()-1;self.u=r.gamma();self.W,self.C,self.P,self.F=pef_widths(self.n,self.u)
        self.k=r.get(self.C);self.payload=r.get(self.P);self.directory_start=r.pos;self.payload_start=r.pos+self.k*self.F
        self.bit_length=self.payload_start+self.payload;r.finish(self.bit_length)
        need((self.n==0 and self.k==0) or 1<=self.k<=self.n,'block count')
        count=0;last=-1;offset=0
        for j in range(self.k):
            mode,mx,end,off=self.row(j);need(mode<3 and last<mx<self.u and count<end<=self.n and off==offset,'directory')
            cost,wanted=choice(end-count,mx-last);need(mode==wanted,'noncanonical local mode')
            vals=self.values(j);need(len(vals)==end-count and vals[-1]==mx,'last value')
            need(all(last<z<self.u for z in vals) and all(x<y for x,y in zip(vals,vals[1:])),'local sequence')
            count,last,offset=end,mx,offset+cost
        need(count==self.n and offset==self.payload,'totals')
    def row(self,j):
        integer(j,0,self.k);p=self.directory_start+j*self.F;r=self.r
        mode=r.at(p,2);p+=2;mx=r.at(p,self.W);p+=self.W;end=r.at(p,self.C);p+=self.C
        return mode,mx,end,r.at(p,self.P)
    def values(self,j):
        mode,mx,end,off=self.row(j);prev=self.row(j-1) if j else (0,-1,0,0);base=prev[1]+1;k=end-prev[2];u=mx-base+1
        integer(k,1);integer(u,k)
        if mode==0:return list(range(base,mx+1))
        at=self.payload_start+off
        if mode==1:return [base+i for i in range(u) if self.r.at(at+i,1)]
        ell,h,cost=ef_params(k,u);vals=[];high=at+k*ell
        for pos in range(k+h):
            if self.r.at(high+pos,1):
                i=len(vals);need(i<k,'too many high ones');low=self.r.at(at+i*ell,ell);vals.append(base+((pos-i)<<ell)+low)
        need(len(vals)==k,'high count');return vals
    def access(self,i):
        integer(i,0,self.n);lo,hi=0,self.k
        while lo<hi:
            m=(lo+hi)//2
            if self.row(m)[2]<=i:lo=m+1
            else:hi=m
        previous=self.row(lo-1)[2] if lo else 0
        return self.values(lo)[i-previous]
    def next_geq(self,x):
        integer(x,0,self.u+1);lo,hi=0,self.k
        while lo<hi:
            m=(lo+hi)//2
            if self.row(m)[1]<x:lo=m+1
            else:hi=m
        if lo==self.k:return None
        values=self.values(lo);j=bisect_left(values,x);previous=self.row(lo-1)[2] if lo else 0
        return previous+j,values[j]
    def decode(self):return [x for j in range(self.k) for x in self.values(j)]

def main():
    bits=[int(x) for x in '00000000111111111010101000010000'];packed=pack_rrr(bits);rrr=RRR(packed)
    need(rrr.decode()==bits,'RRR main');rows=[]
    for j in range(rrr.q):
        c,z=class_rank(rrr.block(j));rows.append({'block':j,'class':c,'offset':z,'width':width(comb(len(rrr.block(j)),c)),'directory':rrr.directory(j)})
    xs=[0,1,2,3,1000,1001,1002,1003];cuts,dp=optimal_partition(xs,1024);data=pack_pef(xs,1024);pef=PartitionedEF(data)
    need(pef.decode()==xs,'PEF main');whole=PartitionedEF(pack_pef(xs,1024,[len(xs)]));single=PartitionedEF(pack_pef(xs,1024,list(range(1,len(xs)+1))))
    need(cuts==[4,5,8] and dp[-1]==81,'main optimal partition')
    rrr_migrations=[]
    for label,bb in [('two runs',[0]*2048+[1]*2048),('alternating',[0,1]*2048)]:
        vv=RRR(pack_rrr(bb,64,16));need(vv.decode()==bb,'RRR density migration')
        rrr_migrations.append({'shape':label,'length':4096,'block_width':64,'blocks_per_superblock':16,'ones':vv.ones,'payload_bits':vv.payload,'total_bits':vv.bit_length,'bytes':len(vv.r.data)})
    migrations=[]
    for ys,u in [([],1),([0],1),([0,1,2,3,4,5,6,7],1024),([0,128,256,384,512,640,768,896],1024),([8,9,10,11,12,13,14,15,16,18,20,22,27],32)]:
        code=pack_pef(ys,u);p=PartitionedEF(code);cs,ds=optimal_partition(ys,u);need(p.decode()==ys,'migration decode')
        migrations.append({'values':ys,'universe':u,'cuts':cs,'cost_without_header':ds[-1],'bits':p.bit_length,'bytes':len(code),'queries':[p.next_geq(x) for x in (0,u//2,u)]})
    def rewrite(data,start,w,value):
        out=bytearray(data)
        for j in range(w):
            pos=start+j;mask=1<<(7-pos%8)
            out[pos//8]=(out[pos//8]&~mask)|(((value>>(w-1-j))&1)<<(7-pos%8))
        return bytes(out)
    rejected=[]
    def reject(name,fn):
        try:fn()
        except ValueError:rejected.append(name);return
        raise ValueError('invalid accepted: '+name)
    for name,fn in [('nonbit',lambda:pack_rrr([0,True])),('zero block',lambda:pack_rrr(bits,0)),('rank past end',lambda:rrr.rank(33)),('select zero',lambda:rrr.select(0)),('RRR extra byte',lambda:RRR(packed+b'\0')),('RRR truncation',lambda:RRR(packed[:-1])),('duplicates',lambda:pack_pef([1,1],4)),('outside universe',lambda:pack_pef([4],4)),('missing cut',lambda:pack_pef(xs,1024,[4])),('bool cut',lambda:pack_pef([0],1,[True])),('PEF extra byte',lambda:PartitionedEF(data+b'\0')),('PEF truncation',lambda:PartitionedEF(data[:-1])),('query past universe',lambda:pef.next_geq(1025))]:reject(name,fn)
    for name,fn in [
        ('unused combination code',lambda:RRR(rewrite(packed,rrr.offsets,7,127))),
        ('nonzero RRR padding',lambda:RRR(rewrite(packed,rrr.bit_length,1,1))),
        ('unused PEF mode',lambda:PartitionedEF(rewrite(data,pef.directory_start,2,3))),
        ('wrong payload pointer',lambda:PartitionedEF(rewrite(data,pef.directory_start+pef.F+2+pef.W+pef.C,pef.P,1))),
        ('nonzero PEF padding',lambda:PartitionedEF(rewrite(pack_pef([0],1),12,1,1)))]:reject(name,fn)
    print(json.dumps({'rrr':{'input':''.join(map(str,bits)),'rows':rows,'ones':rrr.ones,'payload_bits':rrr.payload,'total_bits':rrr.bit_length,'hex':packed.hex(),'rank_21':rrr.rank(21),'select_11':rrr.select(11),'select_14':rrr.select(14)},'partitioned_ef':{'values':xs,'universe':1024,'cuts':cuts,'dp':dp,'rows':[pef.row(j) for j in range(pef.k)],'payload_bits':pef.payload,'directory_bits':pef.k*pef.F,'total_bits':pef.bit_length,'hex':data.hex(),'one_block_bits':whole.bit_length,'singleton_bits':single.bit_length,'next_500':pef.next_geq(500),'next_1004':pef.next_geq(1004),'access_6':pef.access(6)},'rrr_migrations':rrr_migrations,'migrations':migrations,'rejections':rejected},ensure_ascii=False,indent=2))
if __name__=='__main__':main()
