import * as Q from "./quantum.js";
import {shell,range,sel,fld,btn,table,matrix,stats,bars,amplitudes,complexPlot,sphere,plot,f,zstr,pct,esc,val,sharedShots,download} from "./labs.js";
function oracles(root){shell(root,sel("type","Algorithm",[["Deutsch","Deutsch"],["DJ-constant","Deutsch–Jozsa: constant"],["DJ-balanced","Deutsch–Jozsa: balanced parity family"],["BV","Bernstein–Vazirani"]],"BV")+sel("n","Input qubits",[1,2,3,4,5],3)+fld("secret","Secret / truth-table encoding",5,0,31,1)+sel("step","Evolution stage",[["0","Preparation"],["1","Hadamards"],["2","Phase oracle"],["3","Final Hadamards"]],3)+sharedShots(),{},s=>{
 const n=s.type==="Deutsch"?1:Number(s.n);if(s.type==="Deutsch"&&s.secret>3)throw Error("Deutsch truth table uses integers 0–3: f(0) is the low bit.");
 if(s.type!=="Deutsch"&&(s.secret>=2**n||(s.type==="DJ-balanced"&&s.secret===0)))throw Error("Use a secret fitting the register; balanced parity requires a nonzero mask.");
 const r=Q.oracleAlgorithm(n,s.type,s.secret),counts=Q.sample(r.probabilities,s.shots,s.seed);
 return '<p class="fine">Deutsch encodes f(0) in bit 0 and f(1) in bit 1 of the integer. DJ balanced instances use f(x)=s·x with nonzero s. The target |−⟩ is factored out using exact phase kickback.</p>'+table(["x","f(x)","Oracle phase"],r.table.map((v,x)=>[Q.bits(x,n),v,v?"−1":"+1"]))+'<h4>State at selected stage</h4>'+amplitudes(r.history[Number(s.step)],n)+'<h4>Final measurement</h4>'+bars(r.probabilities,Array.from({length:2**n},(_,x)=>Q.bits(x,n)),counts);
 },{resample:(s,b,r)=>val(r,"seed",s.seed+1)});}
function simonLab(root){shell(root,sel("n","Input bits",[2,3,4,5],3)+fld("mask","Hidden nonzero mask (integer)",5,1,31,1)+fld("shots","Oracle samples",8,1,100,1)+fld("seed","Random seed",9,0,1e9,1)+btn("resample","New sample"),{},s=>{
 const n=Number(s.n),r=Q.simon(n,s.mask,s.shots,s.seed);return stats([["Equation rank",r.rank+" / "+(n-1)],["Nonzero candidates",r.candidates.filter(Boolean).map(x=>Q.bits(x,n)).join(", ")||"None"]])+table(["Sample y","Constraint y·s = 0"],r.observations.map(y=>[Q.bits(y,n),Array.from({length:n},(_,j)=>((y>>(n-1-j))&1)?"s"+j:"").filter(Boolean).join(" ⊕ ")+" = 0"]))+
 '<h4>Reduced independent rows over GF(2)</h4><p class="mono">'+r.rref.map(x=>Q.bits(x,n)).join(" ; ")+'</p><h4>Exact marginal after Hadamards</h4>'+bars(r.probabilities,Array.from({length:2**n},(_,x)=>Q.bits(x,n)))+'<details><summary>Inspect the oracle pairs</summary>'+table(["x","f(x)=min(x,x⊕s)"],r.pairs.map(v=>[Q.bits(v.x,n),Q.bits(v.fx,n)]))+'</details>';
 },{resample:(s,b,r)=>val(r,"seed",s.seed+1)});}
function groverLab(root){shell(root,sel("n","Qubits",[2,3,4,5,6],3)+'<label>Marked integers, separated by commas<input name="marked" value="5" type="text"></label>'+fld("iterations","Grover iterations",2,0,20,1)+sharedShots(),{},s=>{
 const n=Number(s.n),parts=s.marked.split(",").map(x=>x.trim());if(parts.some(x=>!/^\d+$/.test(x)))throw Error("Enter marked integers separated by commas.");const marked=[...new Set(parts.map(Number))],r=Q.grover(n,marked,s.iterations),counts=Q.sample(r.probabilities,s.shots,s.seed),curve=Array.from({length:21},(_,k)=>Math.sin((2*k+1)*Math.asin(Math.sqrt(marked.length/2**n)))**2);
 return stats([["Marked success",pct(r.success)],["Rotation prediction",pct(r.theory)],["Oracle calls",s.iterations]])+plot(curve,{title:"Success probability versus number of Grover iterations"})+bars(r.probabilities,Array.from({length:2**n},(_,x)=>Q.bits(x,n)),counts)+'<details><summary>Amplitudes at each iteration</summary>'+r.history.map((v,k)=>'<h4>Iteration '+k+'</h4>'+amplitudes(v,n)).join("")+'</details>';
 },{resample:(s,b,r)=>val(r,"seed",s.seed+1)});}
function qftLab(root){shell(root,sel("n","Qubits",[2,3,4,5],3)+sel("input","Input pattern",[["basis","Basis |x⟩"],["comb","Periodic comb"],["phase","Phase ramp"]],"basis")+fld("x","Basis index / comb offset",1,0,31,1)+fld("period","Comb spacing",2,1,32,1)+range("phi","Phase ramp φ",.25,0,1,.001)+sel("inverse","Transform sign",[["no","Forward QFT (+)"],["yes","Inverse QFT (−)"]],"no"),{},s=>{
 const n=Number(s.n),N=2**n;if(s.x>=N)throw Error("Input index must fit the register.");
 let v;if(s.input==="basis")v=Q.ket(n,s.x);else if(s.input==="comb")v=Q.normalize(Array.from({length:N},(_,x)=>Q.C(x>=s.x&&(x-s.x)%s.period===0?1:0)));else v=Array.from({length:N},(_,x)=>Q.scale(Q.phase(2*Math.PI*s.phi*x),1/Math.sqrt(N)));
 const out=Q.qft(v,s.inverse==="yes"),back=Q.qft(out,s.inverse!=="yes");
 return stats([["Input norm²",f(Q.norm(v)**2)],["Inverse recovery fidelity",f(Q.abs2(Q.inner(v,back)),9)]])+'<div class="visual-pair"><div><h4>Input amplitudes</h4>'+complexPlot(v)+'</div><div><h4>Output amplitudes</h4>'+complexPlot(out)+'</div></div>'+amplitudes(out,n)+'<p class="fine">Complex plots show the first 16 amplitudes. Computational probabilities alone can hide Fourier phase information.</p>';
 });}
function qpeLab(root){shell(root,range("phi","Eigenphase φ",.375,0,.999,.001)+sel("m","Control qubits",[2,3,4,5,6,7],3)+sharedShots(),{},s=>{
 const m=Number(s.m),r=Q.qpe(m,s.phi),counts=Q.sample(r.probabilities,s.shots,s.seed),mode=counts.indexOf(Math.max(...counts)),err=Math.min(Math.abs(mode/2**m-s.phi),1-Math.abs(mode/2**m-s.phi));
 return stats([["Most frequent y / 2ᵐ",mode+" / "+2**m],["Phase estimate",f(mode/2**m)],["Circular error",f(err)],["U calls if powers repeated",2**m-1]])+bars(r.probabilities,Array.from({length:2**m},(_,y)=>Q.bits(y,m)+" · "+f(y/2**m,3)),counts,32)+'<details><summary>Phase-kicked control amplitudes before inverse QFT</summary>'+amplitudes(r.before,m)+'</details>';
 },{resample:(s,b,r)=>val(r,"seed",s.seed+1)});}
function shorLab(root){
 let cache=null,outcome=64;
 shell(root,sel("N","Composite N",[15,21],15)+fld("a","Base a",2,2,20,1)+sel("m","Control qubits",[4,5,6,7,8,9,10],8)+fld("seed","Sample seed",4,0,1e9,1)+btn("sample","Sample Fourier outcome")+fld("y","Inspect outcome y",64,0,1023,1)+btn("save","Export exact probabilities"),{},s=>{
 const N=Number(s.N),m=Number(s.m);if(s.a>=N)throw Error("Choose 1 < a < N.");const key=[N,s.a,m].join("-");
 if(!cache||cache.key!==key)cache={...Q.shor(N,s.a,m),key};
 if(cache.gcd)return '<p>A classical gcd already finds factor '+cache.gcd+' and '+N/cache.gcd+'. No quantum order finding is needed for this base.</p>';
 if(s.y>=cache.Q)throw Error("Outcome y must be smaller than Q=2ᵐ.");
 const d=Q.factorFromOutcome(N,s.a,m,s.y),top=cache.probabilities.map((p,y)=>({p,y})).sort((a,b)=>b.p-a.p).slice(0,16).sort((a,b)=>a.y-b.y);
 return stats([["Control range Q",cache.Q],["Inspected y / Q",s.y+" / "+cache.Q],["Factors from this sample",d.factors.length?d.factors.join(", "):"No verified factors"]])+
 (cache.Q<N*N?'<p class="notice">Q < N²: this coarse experiment does not meet the elementary continued-fraction precision guarantee.</p>':"")+
 table(["Convergent","Denominator"],d.cf.map(c=>[c.p+"/"+c.q,c.q]))+'<p>Tested period candidates (denominators and up to eight multiples): '+d.tried.join(", ")+'. Each accepted period passes aʳ mod N = 1 before gcd extraction.</p>'+
 '<h4>Largest exact Fourier probabilities</h4>'+bars(top.map(x=>x.p),top.map(x=>String(x.y)))+
 '<details><summary>Reveal the order for checking and inspect modular groups</summary><p>Reference order, found by classical enumeration for this small test: '+cache.referenceOrder+'. This reference is not used to generate the distribution or infer candidates.</p>'+table(["aˣ mod N","Compatible x (first 20)"],cache.groups.map(([y,xs])=>[y,xs.slice(0,20).join(", ")+(xs.length>20?" …":"")]))+'</details><p class="fine">Exact classical simulation of modular-exponentiation groups and the inverse-QFT marginal. It is not a hardware run or a cryptographic arithmetic gate decomposition.</p>';
 },{sample:(s,b,r)=>{if(cache?.probabilities){val(r,"y",Q.draw(cache.probabilities,Q.rng(s.seed)));val(r,"seed",s.seed+1);}},save:()=>{if(cache?.probabilities)download("shor-probabilities.csv","y,probability\n"+cache.probabilities.map((p,y)=>y+","+p).join("\n"),"text/csv");}});
}
function hamiltonianLab(root){shell(root,range("time","Evolution time t",1,0,5)+fld("steps","Product-formula steps",10,1,200,1),{},s=>{const r=Q.hamiltonian(s.time,s.steps);return stats([["State-vector norm error",f(r.error,8)],["Fidelity",f(r.fidelity,8)]])+table(["State","Amplitude 0","Amplitude 1"],[["Exact",...r.exact.map(zstr)],["Product formula",...r.approx.map(zstr)]])+plot(Array.from({length:30},(_,j)=>Q.hamiltonian(s.time,j+1).error),{max:Math.max(.01,...Array.from({length:30},(_,j)=>Q.hamiltonian(s.time,j+1).error)),title:"Product formula error versus steps"})+'<p class="fine">H=X+Z, ℏ=1. Each step applies Rz(2t/r), then Rx(2t/r). Plot covers r=1,…,30. Physical gate noise is excluded.</p>';});}
function vqeLab(root){shell(root,range("theta","θ",2,0,2*Math.PI)+range("phi","φ",Math.PI,-Math.PI,Math.PI)+range("rate","Gradient step size",.15,.01,.5)+btn("step","Take θ gradient step")+btn("solve","Set exact ground state"),{},s=>{
 const r=Q.vqe(s.theta,s.phi);return stats([["Energy ⟨Z+X/2⟩",f(r.energy,8)],["Exact ground energy",f(r.ground,8)],["θ gradient",f(r.gradient)]])+sphere(Q.bloch(r.state))+plot(Array.from({length:81},(_,j)=>Q.vqe(j/80*2*Math.PI,s.phi).energy),{min:-1.2,max:1.2,title:"Energy versus θ from 0 to 2π"})+'<p class="fine">Gradient steps optimize θ with φ fixed. An unsuitable fixed φ can restrict the attainable minimum. Exact expectations are used; hardware shots would add uncertainty.</p>';
 },{step:(s,b,r)=>val(r,"theta",(s.theta-s.rate*Q.vqe(s.theta,s.phi).gradient+2*Math.PI)%(2*Math.PI)),solve:(s,b,r)=>{val(r,"theta",Math.PI-Math.atan(.5));val(r,"phi",Math.PI);}});}
function qaoaLab(root){shell(root,range("gamma","Cost angle γ",.6,-Math.PI,Math.PI)+range("beta","Mixer angle β",-.3,-Math.PI/2,Math.PI/2),{},s=>{const r=Q.qaoa(s.gamma,s.beta);return '<p>Triangle graph: edges (0,1), (1,2), (0,2). Cost = number of unequal endpoint bits.</p>'+stats([["Expected cut",f(r.expectation)],["Best possible cut",r.optimum]])+bars(r.probabilities,Array.from({length:8},(_,x)=>Q.bits(x,3)+" · cut "+r.costs[x]))+'<p class="fine">One exact QAOA layer: exp(−iγC), then exp(−iβΣX). The mixer uses Rx(2β) on each qubit.</p>';});}
function bytesText(log2){const units=["B","KiB","MiB","GiB","TiB","PiB","EiB","ZiB","YiB"];const i=Math.min(units.length-1,Math.floor(log2/10));return log2>100?"2^"+f(log2)+" bytes":f(2**(log2-10*i),3)+" "+units[i];}
function costLab(root){shell(root,range("n","Qubits",30,1,70,1)+sel("bytes","Bytes per complex amplitude",[[8,"8 — complex single"],[16,"16 — complex double"]],16),{},s=>{const b=Math.log2(Number(s.bytes));return stats([["State vector",bytesText(b+s.n)],["Density matrix",bytesText(b+2*s.n)],["Amplitude count","2^"+s.n]])+table(["Representation","Raw storage formula","Important qualification"],[["State vector","bytes × 2ⁿ","Pure state, no work arrays"],["Density matrix","bytes × 4ⁿ","Dense mixed-state storage"],["Stabilizer tableau","O(n²) bits","Clifford and compatible measurement/noise"],["Matrix product state","O(nχ²) parameters","Bond dimension depends on entanglement"]])+'<p class="fine">This computes storage estimates, not a large state. Raw array costs exclude metadata and temporary buffers.</p>';});}
function repetitionLab(root){shell(root,sel("d","Odd repetition length",[3,5,7,9,11],3)+fld("error","Injected X-error bit mask",2,0,2047,1)+range("p","Independent physical flip probability",.1,0,.5),{},s=>{
 const d=Number(s.d);if(s.error>=2**d)throw Error("Error mask must fit the code length.");const r=Q.repetitionSyndrome(s.error,d);
 return table(["Object","Bits"],[["Injected error",r.error.join("")],["Neighbor parity syndrome",r.syndrome.join("")],["Minimum-weight recovery",r.correction.join("")],["Residual",r.residual.join("")]])+stats([["Logical X failure?",r.logicalFailure?"Yes":"No"],["Exact independent failure",pct(Q.repetition(s.p,d))],["Unencoded failure",pct(s.p)]])+table(["Length d","Exact majority failure"],[3,5,7,9,11].map(d=>[d,pct(Q.repetition(s.p,d))]))+'<p class="fine">Perfect syndrome measurement and recovery; independent X errors only. This is not protection against arbitrary quantum phase errors.</p>';
 });}
let decodeTable=null;
function decoder(){
 if(decodeTable)return decodeTable;decodeTable=new Map();const chars=["X","Y","Z"];
 const add=(err)=>{const s=Q.syndromeFor(err).join("");if(!decodeTable.has(s))decodeTable.set(s,err);};
 add("IIIIIIIII");
 function gen(pos,left,e){if(!left){add(e.join(""));return;}for(let q=pos;q<=9-left;q++)for(const p of chars){e[q]=p;gen(q+1,left-1,e);e[q]="I";}}
 for(let w=1;w<=9&&decodeTable.size<256;w++)gen(0,w,Array(9).fill("I"));
 return decodeTable;
}
const pauliProduct=(a,b)=>{const map={I:0,X:1,Z:2,Y:3},back=["I","X","Z","Y"];return a.split("").map((p,j)=>back[map[p]^map[b[j]]]).join("");};
function surfaceLab(root){
 let error=Array(9).fill("I");
 shell(root,btn("clear","Clear errors")+btn("logicalX","Set logical X string")+btn("logicalZ","Set logical Z string")+btn("center","Single Y at center"),{},s=>{
 const e=error.join(""),checks=Q.rotatedSurface3(),synd=Q.syndromeFor(e),recovery=decoder().get(synd.join("")),residual=pauliProduct(e,recovery),lx=!Q.pauliCommutes(residual,"ZZZIIIIII"),lz=!Q.pauliCommutes(residual,"XIIXIIXII");
 return '<div class="surface-grid" aria-label="Nine data qubits">'+error.map((p,j)=>'<button type="button" data-action="toggle" data-q="'+j+'" class="pauli-'+p+'" aria-label="Qubit '+j+' has '+p+'; click to cycle error"><small>q'+j+'</small><strong>'+p+'</strong></button>').join("")+'</div>'+
 table(["Check","Pauli support q0…q8","Syndrome bit"],checks.map((c,j)=>[c.name,c.p,synd[j]]))+stats([["Minimum-weight recovery",recovery],["Residual",residual],["Residual logical action",lx?(lz?"Y":"X"):(lz?"Z":"Identity")]])+
 '<p class="fine">Click cycles I→X→Y→Z. Nine data qubits, eight ideal commuting checks. Decoder minimizes Pauli weight with equal X,Y,Z costs; ties use deterministic enumeration. This code-capacity model excludes syndrome-circuit faults and repeated noisy checks. Recovery differs from the actual error when a stabilizer-equivalent pattern is sufficient.</p>';
 },{toggle:(s,b)=>{const j=Number(b.dataset.q);error[j]=["I","X","Y","Z"][(["I","X","Y","Z"].indexOf(error[j])+1)%4];},clear:()=>error=Array(9).fill("I"),logicalX:()=>error="XIIXIIXII".split(""),logicalZ:()=>error="ZZZIIIIII".split(""),center:()=>{error=Array(9).fill("I");error[4]="Y";}});
}
function resourcesLab(root){shell(root,fld("p","Physical error p",.001,.000001,.1)+fld("threshold","Assumed threshold",.01,.00001,.5)+fld("A","Model prefactor A",.1,.00001,1)+sel("d","Odd code distance",[3,5,7,9,11,13,15,17,19,21,25,31],17)+fld("logical","Data logical patches",100,1,1e6,1)+fld("factories","Additional equal-size patches",0,0,1e6,1)+fld("locations","Failure locations L",1e8,1,1e15,1)+fld("depth","Sequential logical depth",1e6,1,1e12,1)+fld("cycle","Code cycle (seconds)",.000001,.000000001,1),{},s=>{const r=Q.resource({...s,d:Number(s.d)});return (!r.below?'<p class="notice">p ≥ threshold. Favorable distance suppression is unavailable; this extrapolation is not a calibrated above-threshold model.</p>':"")+stats([["Modeled pL",r.logicalError.toExponential(3)],["Union failure bound",pct(r.unionBound)],["Bare physical qubits",r.physical.toLocaleString()],["Modeled time",f(r.time,3)+" s"]])+table(["Assumption","Value"],[["Physical qubits per patch","2d²−1 = "+r.physicalPerPatch],["Logical error equation","A(p/pth)^((d+1)/2)"],["Schedule","depth × d × cycle"],["Not included","Distinct factory internals, routing, retries, communication, calibration"]])+'<p class="fine">An illustrative resource model, not a prediction for a named machine. Additional patches are assumed equal-sized; this is not a magic-state factory throughput calculation.</p>';});}
function mitigationLab(root){shell(root,range("lambda","Base noise λ",.2,0,1)+fld("shots","Shots per noise scale",1000,1,10000,1)+fld("seed","Random seed",8,0,1e9,1)+btn("resample","New sample"),{},s=>{const r=Q.mitigation(s.lambda,s.shots,s.seed);return stats([["Ideal value",1],["Raw estimate",f(r.estimates[0])],["Linear extrapolation",f(r.linear)],["Quadratic extrapolation",f(r.quadratic)]])+table(["Noise scale","Exact noisy mean","Sample estimate"],r.scales.map((v,j)=>[v,f(r.means[j]),f(r.estimates[j])]))+stats([["Exact linear bias",f(r.linearExact-1)],["Exact quadratic bias",f(r.quadraticExact-1)],["Linear standard deviation",f(Math.sqrt(r.varianceLinear))]])+'<p class="fine">Toy model E(sλ)=exp(−sλ). Linear weights (2,−1); quadratic weights (3,−3,1). Estimates are not clipped, so fluctuations outside [−1,1] remain visible.</p>';},{resample:(s,b,r)=>val(r,"seed",s.seed+1)});}
function bb84Lab(root){shell(root,fld("shots","Transmitted signals",512,8,10000,1)+range("eve","Intercept–resend fraction",0,0,1)+range("noise","Independent postmeasurement bit flip",0,0,.5)+fld("seed","Random seed",10,0,1e9,1)+btn("resample","New transmission")+btn("export","Export full transcript"),{},s=>{
 const r=Q.bb84(s.shots,s.eve,s.noise,s.seed);s.transcript=r.rows;
 return stats([["Sifted signals",r.sifted],["Observed sifted QBER",r.qber===null?"No matching bases":pct(r.qber)],["Model expected QBER",pct(r.expectedQber)],["Asymptotic ideal secret fraction",pct(r.asymptoticFraction)]])+table(["#","Alice bit","Alice basis","Eve basis","Bob basis","Bob bit","Sifted?"],r.rows.slice(0,32).map((x,j)=>[j+1,x.bit,x.ab?"X":"Z",x.eb===null?"—":x.eb?"X":"Z",x.bb?"X":"Z",x.out,x.keep?(x.error?"Yes · error":"Yes"):"No"]))+'<p class="fine">First 32 signals shown; export contains all. Ideal single photons and ideal basis measurements, with the selected teaching noise model. QBER uses the entire simulated sifted key for diagnosis; real protocols sacrifice a sample. The asymptotic 1−2h₂(Q) fraction is not a finite-key guarantee or a produced secure key. Authentication, reconciliation, privacy amplification, loss and detector attacks are not implemented as a production protocol.</p>';
 },{resample:(s,b,r)=>val(r,"seed",s.seed+1),export:s=>download("bb84-transcript.json",JSON.stringify(s.transcript,null,2),"application/json")});}
export const advancedLabs={oracles,simon:simonLab,grover:groverLab,qft:qftLab,qpe:qpeLab,shor:shorLab,hamiltonian:hamiltonianLab,vqe:vqeLab,qaoa:qaoaLab,cost:costLab,repetition:repetitionLab,surface:surfaceLab,resources:resourcesLab,mitigation:mitigationLab,bb84:bb84Lab};


export {decoder,pauliProduct};
