"""Stable natural-run Powersort and shorter-side galloping merge.
Integer keys with original record identity. Standard library; stdout only.
Trace storage and cubic optimal-tree diagnostics are explicit optional extras.
"""
from math import log2
from itertools import product
import json

def check(ok,msg):
 if not ok:raise ValueError(msg)
def integer(x):return type(x) is int

def gallop(a,start,end,key,upper,cost):
 # Return first key >= target (upper=False) or > target (upper=True).
 # The array is assumed sorted. No reading of the virtual endpoint.
 def before(i):
  cost['comparisons']+=1
  return a[i][0]<=key if upper else a[i][0]<key
 lo=start;step=1;probe=start
 while probe<end:
  if not before(probe):break
  lo=probe+1;step*=2;probe=start+step-1
 hi=min(probe,end)
 while lo<hi:
  mid=(lo+hi)//2
  if before(mid):lo=mid+1
  else:hi=mid
 return lo

def merge_range(a,left,middle,right,method='gallop',search_trace=None):
 # a[left:middle] and a[middle:right] are sorted, adjacent original intervals.
 check(method in ('gallop','linear'),'merge method');out=[];cost={'comparisons':0};writes=0
 if method=='linear':
  i,j=left,middle
  while i<middle and j<right:
   cost['comparisons']+=1
   if a[i][0]<=a[j][0]:out.append(a[i]);i+=1
   else:out.append(a[j]);j+=1
  while i<middle:out.append(a[i]);i+=1
  while j<right:out.append(a[j]);j+=1
 else:
  short_left=middle-left<=right-middle
  i,iend,j,jend=(left,middle,middle,right) if short_left else (middle,right,left,middle)
  # Short left: emit strictly smaller long/right keys. Short right: emit <=.
  upper=not short_left
  while i<iend and j<jend:
   old=j;before=cost['comparisons'];j=gallop(a,j,jend,a[i][0],upper,cost)
   if search_trace is not None:search_trace.append(dict(pivot=a[i],start=old,end=jend,boundary=j,upper=upper,comparisons=cost['comparisons']-before))
   for k in range(old,j):out.append(a[k])
   out.append(a[i]);i+=1
  while i<iend:out.append(a[i]);i+=1
  while j<jend:out.append(a[j]);j+=1
 for k,record in enumerate(out,left):a[k]=record;writes+=1
 return dict(comparisons=cost['comparisons'],volume=right-left,buffer_writes=len(out),array_writes=writes)

def merge_checked(left,right,method='gallop'):
 # Public demonstration includes and COUNTS a separate sortedness pass.
 check(all(integer(k) for k in left+right),'integer keys')
 validation=0
 for run in (left,right):
  for i in range(1,len(run)):
   validation+=1;check(run[i-1]<=run[i],'sorted input runs')
 records=[(x,i) for i,x in enumerate(left+right)];trace=[]
 stats=merge_range(records,0,len(left),len(records),method,trace)
 return dict(output=records,stats=stats,sortedness_comparisons=validation,searches=trace)

def node_power(left,middle,right,n):
 check(all(integer(x) for x in (left,middle,right,n)) and 0<=left<middle<right<=n,'adjacent nonempty original runs')
 denominator=2*n;a=left+middle;b=middle+right;power=0
 while True:
  a*=2;b*=2;power+=1
  ah=a>=denominator;bh=b>=denominator
  if ah!=bh:return power
  if ah:a-=denominator;b-=denominator

def powersort(keys,method='gallop',record_trace=True):
 check(all(integer(x) for x in keys),'integer keys');check(method in ('gallop','linear'),'method');check(type(record_trace)is bool,'trace switch')
 a=[(key,i) for i,key in enumerate(keys)];n=len(a);runs=[];nodes=[];events=[];boundaries=[];stack=[];serial=0
 stats=dict(n=n,runs=0,run_comparisons=0,merge_comparisons=0,merge_volume=0,buffer_writes=0,array_merge_writes=0,reversal_writes=0,initial_record_writes=n,power_bit_steps=0,stack_peak=0)
 def node(data):
  nonlocal serial
  i=serial;serial+=1
  if record_trace:nodes.append(dict(id=i,**data))
  return i
 def first_run(left):
  if left+1==n:right=n;descending=False
  else:
   stats['run_comparisons']+=1;descending=a[left+1][0]<a[left][0];right=left+2
   while right<n:
    stats['run_comparisons']+=1
    if (a[right][0]<a[right-1][0]) if descending else (a[right][0]>=a[right-1][0]):right+=1
    else:break
   if descending:
    i,j=left,right-1
    while i<j:a[i],a[j]=a[j],a[i];stats['reversal_writes']+=2;i+=1;j-=1
  ordinal=stats['runs'];stats['runs']+=1
  if record_trace:runs.append(dict(index=ordinal,left=left,right=right,descending=descending))
  return (left,right,node(dict(kind='leaf',run=ordinal,left=left,right=right)))
 def merge(u,v):
  left,middle,li=u;mid,right,ri=v;check(middle==mid,'merge adjacency')
  c=merge_range(a,left,middle,right,method)
  stats['merge_comparisons']+=c['comparisons'];stats['merge_volume']+=c['volume'];stats['buffer_writes']+=c['buffer_writes'];stats['array_merge_writes']+=c['array_writes']
  parent=node(dict(kind='merge',left=left,right=right,children=(li,ri)))
  if record_trace:events.append(dict(action='merge',node=parent,left=left,middle=middle,right=right,**c))
  return (left,right,parent)
 if not n:return dict(output=a,root=None,runs=runs,nodes=nodes,events=events,boundaries=boundaries,stats=stats)
 current=first_run(0)
 while current[1]<n:
  following=first_run(current[1]);left,middle=current[:2];right=following[1]
  power=node_power(left,middle,right,n);stats['power_bit_steps']+=power
  if record_trace:boundaries.append(dict(left=left,middle=middle,right=right,power=power))
  while stack and stack[-1][1]>power:
   previous,old=stack.pop()
   if record_trace:events.append(dict(action='pop',node=previous[2],power=old))
   current=merge(previous,current)
  check(not stack or stack[-1][1]<power,'strictly increasing pending powers')
  stack.append((current,power));stats['stack_peak']=max(stats['stack_peak'],len(stack))
  if record_trace:events.append(dict(action='push',node=current[2],power=power))
  current=following
 while stack:
  previous,old=stack.pop()
  if record_trace:events.append(dict(action='pop',node=previous[2],power=old))
  current=merge(previous,current)
 check(stats['run_comparisons']==n-1,'single pass adjacent comparisons')
 return dict(output=a,root=current[2],runs=runs,nodes=nodes,events=events,boundaries=boundaries,stats=stats)

def tree_certificate(result):
 # Linear in stored nodes/runs. A diagnostic, not part of the quiet sorter.
 nodes=result['nodes'];n=result['stats']['n'];root=result['root']
 if not n:check(root is None and not nodes,'empty tree');return dict(depths=[],volume=0)
 check(root==len(nodes)-1,'last node root');parents=[0]*len(nodes);volume=0;leaves=[];runs=result['runs'];seen=set();end=0
 for j,run in enumerate(runs):
  check(run['index']==j and run['left']==end and end<run['right']<=n,'original run partition');end=run['right']
 check(end==n and len(nodes)==2*len(runs)-1,'full leaf tree size')
 for i,row in enumerate(nodes):
  check(row['id']==i and 0<=row['left']<row['right']<=n,'node interval')
  if row['kind']=='merge':
   u,v=row['children'];check(0<=u<i and 0<=v<i and u!=v,'earlier children');parents[u]+=1;parents[v]+=1
   check(nodes[u]['left']==row['left'] and nodes[u]['right']==nodes[v]['left'] and nodes[v]['right']==row['right'],'ordered adjacent children');volume+=row['right']-row['left']
  else:
   check(row['kind']=='leaf','leaf kind');j=row['run'];check(integer(j) and 0<=j<len(runs) and j not in seen,'unique leaf ordinal');seen.add(j)
   check((row['left'],row['right'])==(runs[j]['left'],runs[j]['right']),'leaf original interval');leaves.append(row)
 check(len(seen)==len(runs),'every original leaf');check(parents[root]==0 and all(parents[i]==1 for i in range(root)),'tree parent counts');check(nodes[root]['left']==0 and nodes[root]['right']==n,'whole interval')
 depths=[0]*len(leaves);pending=[(root,0)]
 while pending:
  i,d=pending.pop();row=nodes[i]
  if row['kind']=='leaf':depths[row['run']]=d
  else:pending.extend((j,d+1)for j in row['children'])
 lengths=[r['right']-r['left']for r in result['runs']]
 check(sum(w*d for w,d in zip(lengths,depths))==volume==result['stats']['merge_volume'],'weighted external length')
 return dict(depths=depths,volume=volume)

def optimal_volume(lengths):
 # Cubic interval DP for a small-instance oracle, not a Powersort component.
 r=len(lengths);check(all(integer(x) and x>0 for x in lengths),'positive leaf weights')
 if r==0:return dict(volume=0,splits=[])
 prefix=[0]
 for x in lengths:prefix.append(prefix[-1]+x)
 D=[[0]*(r+1) for _ in range(r)];split={};candidates=0
 for width in range(2,r+1):
  for i in range(r-width+1):
   j=i+width;options=[]
   for k in range(i+1,j):options.append((D[i][k]+D[k][j],k));candidates+=1
   value,k=min(options);D[i][j]=prefix[j]-prefix[i]+value;split[i,j]=k
 return dict(volume=D[0][r],splits=[dict(left=i,right=j,split=k)for(i,j),k in sorted(split.items())],candidates=candidates)

def cost_bounds(result):
 lengths=[x['right']-x['left']for x in result['runs']];n=sum(lengths);M=result['stats']['merge_volume']
 if not n:return dict(entropy=0,lower_bound=0,upper_bound=0)
 H=sum(w/n*log2(n/w)for w in lengths)
 # Exact integer comparisons, including non-dyadic lengths. Logs display only.
 product_weight=1
 for w in lengths:product_weight*=w**w
 check((1<<M)*product_weight>=n**n,'entropy lower bound')
 check((1<<M)*product_weight<(1<<(2*n))*n**n,'midpoint tree entropy upper bound')
 return dict(entropy=H,lower_bound=n*H,upper_bound=n*(H+2))

def main():
 blocks=[[8,14,27],[1,3,9,11,17,19,29],[6,25],[0,2,4,7,10,13,20,24,31],[5,12,18,30],[15,16,21,22,23,26,28]]
 keys=[x for b in blocks for x in b];main=powersort(keys);linear=powersort(keys,'linear');check(main['output']==sorted((x,i)for i,x in enumerate(keys))==linear['output'],'stable sorted identities');cert=tree_certificate(main);bounds=cost_bounds(main);optimal=optimal_volume(list(map(len,blocks)))
 duplicates=[5,4,4,3,2,2,1];dup=powersort(duplicates);check(dup['output']==sorted((x,i)for i,x in enumerate(duplicates)),'strict reversal stability');tree_certificate(dup);cost_bounds(dup)
 pair_examples=[]
 for left,right in [([31,63,95],list(range(128))),(list(range(128)),[31,63,95]),([2],[0,1,3,4,5])]:
  ga=merge_checked(left,right);li=merge_checked(left,right,'linear');check(ga['output']==li['output']==sorted((x,i)for i,x in enumerate(left+right)),'galloping stable ties');pair_examples.append(dict(left=left,right=right,galloping=ga,linear_comparisons=li['stats']['comparisons']))
 counterexample=None
 for weights in product(range(2,8),repeat=4):
  profile=[];end=sum(weights)
  for w in weights:profile.extend(range(end-w,end));end-=w
  result=powersort(profile);opt=optimal_volume(weights)
  if result['stats']['merge_volume']>opt['volume']:
   counterexample=dict(lengths=weights,powersort=result['stats']['merge_volume'],optimal=opt['volume'],bounds=cost_bounds(result));break
 check(counterexample is not None,'near optimal not exact example')
 # Same lengths, disjoint descending key bands: identical schedule volume.
 disjoint=[];end=len(keys)
 for block in blocks:w=len(block);disjoint.extend(range(end-w,end));end-=w
 other=powersort(disjoint);check(other['stats']['merge_volume']==main['stats']['merge_volume'],'length-only schedule');tree_certificate(other)
 heavy_lengths=[1024]+[2]*128;heavy=[];end=sum(heavy_lengths)
 for w in heavy_lengths:heavy.extend(range(end-w,end));end-=w
 heavy_result=powersort(heavy);tree_certificate(heavy_result);cost_bounds(heavy_result);prefix=heavy_lengths[0];leftfold=0
 for w in heavy_lengths[1:]:prefix+=w;leftfold+=prefix
 boundaries=[]
 for xs in [[],[4],[2,2,2,2],[9,8,7,6]]:
  r=powersort(xs);check(r['output']==sorted((x,i)for i,x in enumerate(xs)),'boundary');tree_certificate(r);cost_bounds(r);boundaries.append(r)
 quiet=powersort(keys,record_trace=False);check(quiet['output']==main['output'] and quiet['stats']==main['stats'] and not quiet['events'] and not quiet['nodes'],'quiet storage interface')
 print(json.dumps(dict(status='PASS',main=dict(input=keys,result=main,certificate=cert,bounds=bounds,optimal=optimal,linear_stats=linear['stats']),duplicates=dup,pairs=pair_examples,nonoptimal=counterexample,same_profile=dict(input=disjoint,stats=other['stats']),heavy_profile=dict(lengths=heavy_lengths,stats=heavy_result['stats'],leftfold_volume=leftfold),boundaries=boundaries),ensure_ascii=False,indent=2))
if __name__=='__main__':main()