use super::{MeasurementRelations, ParityInference};
use crate::photonic::stabilizer::{
Pauli, PauliNoise, PauliProduct, StabilizerSimulator, execute_trajectory,
};
use crate::photonic::{
CompiledFusionNetwork, FusionOutcome, FusionProgram, NodeId, PhotonId, ResourceState,
runtime::shot_seed,
};
use crate::{Error, Result};
use rand::{Rng, SeedableRng};
use rand_chacha::ChaCha8Rng;
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
pub enum PatchBoundary {
Output,
Measure(Pauli),
}
#[derive(Clone, Debug)]
pub struct PhotonLossBurst {
photons: Vec<PhotonId>,
probability: f64,
}
impl PhotonLossBurst {
pub fn new(photons: Vec<PhotonId>, probability: f64) -> Result<Self> {
if photons.is_empty()
|| photons
.iter()
.enumerate()
.any(|(i, q)| photons[..i].contains(q))
{
return Err(Error::Wires("invalid photon loss burst".into()));
}
if !probability.is_finite() || !(0.0..=1.).contains(&probability) {
return Err(Error::InvalidValue("invalid loss burst probability".into()));
}
Ok(Self {
photons,
probability,
})
}
}
#[derive(Clone, Debug, PartialEq)]
pub struct PatchShot {
pub shot_id: u64,
pub fusion_outcomes: Vec<FusionOutcome>,
pub records: Vec<Option<bool>>,
pub output_present: Vec<bool>,
pub inferred: Vec<ParityInference>,
pub expectations: Vec<i8>,
}
#[derive(Clone, Debug)]
pub struct SixRingPatch {
dimensions: [usize; 3],
program: FusionProgram,
boundary: Vec<PauliProduct>,
outputs: Vec<PhotonId>,
}
impl SixRingPatch {
pub(crate) fn from_parts(
dimensions: [usize; 3],
program: FusionProgram,
boundary: Vec<PauliProduct>,
mut outputs: Vec<PhotonId>,
max_photons: usize,
) -> Result<Self> {
let n = program.num_photons();
let limit = max_photons.min(16384);
if n == 0 || dimensions.contains(&0) {
return Err(Error::InvalidDimensions(
"empty finite fusion geometry".into(),
));
}
if n > limit {
return Err(Error::ResourceLimit {
resource: "open six-ring photons",
required: n,
limit,
});
}
let mut used = vec![false; n];
for q in program
.fusions()
.iter()
.flat_map(|f| f.photons.iter().map(|q| q.0))
.chain(
boundary
.iter()
.map(|p| p.terms().first().map_or(n, |(q, _)| *q)),
)
.chain(outputs.iter().map(|q| q.0))
{
if q >= n || used[q] {
return Err(Error::Wires(
"finite geometry reuses or omits a photon identity".into(),
));
}
used[q] = true;
}
if used.iter().any(|v| !v) || boundary.iter().any(|p| p.terms().len() != 1) {
return Err(Error::Wires(
"finite geometry must consume or retain every photon exactly once".into(),
));
}
outputs.sort_by_key(|q| q.0);
Ok(Self {
dimensions,
program,
boundary,
outputs,
})
}
pub fn new(
dimensions: [usize; 3],
faces: [PatchBoundary; 6],
source: NodeId,
analyzer: NodeId,
max_photons: usize,
) -> Result<Self> {
if dimensions.contains(&0) {
return Err(Error::InvalidDimensions(
"empty six-ring patch dimension".into(),
));
}
let cells = dimensions
.iter()
.try_fold(1usize, |v, &n| v.checked_mul(n).ok_or(Error::Overflow))?;
let photons = cells.checked_mul(6).ok_or(Error::Overflow)?;
let limit = max_photons.min(4096);
if photons > limit {
return Err(Error::ResourceLimit {
resource: "open six-ring photons",
required: photons,
limit,
});
}
let mut program = FusionProgram::new();
let mut rings = Vec::with_capacity(cells);
for _ in 0..cells {
rings.push(program.add_resource(ResourceState::six_ring(), source)?);
}
let index = |c: [usize; 3]| (c[0] * dimensions[1] + c[1]) * dimensions[2] + c[2];
let positive = [1, 2, 0];
let mut boundary = Vec::new();
let mut outputs = Vec::new();
for x in 0..dimensions[0] {
for y in 0..dimensions[1] {
for z in 0..dimensions[2] {
let c = [x, y, z];
let cell = index(c);
for axis in 0..3 {
if c[axis] + 1 < dimensions[axis] {
let mut next = c;
next[axis] += 1;
program.fuse(
rings[cell][positive[axis]],
rings[index(next)][positive[axis] + 3],
analyzer,
)?;
}
for negative in [false, true] {
let face = axis + if negative { 3 } else { 0 };
if (negative && c[axis] == 0)
|| (!negative && c[axis] + 1 == dimensions[axis])
{
let photon =
rings[cell][positive[axis] + if negative { 3 } else { 0 }];
match faces[face] {
PatchBoundary::Output => outputs.push(photon),
PatchBoundary::Measure(pauli) => {
boundary.push(PauliProduct::new(vec![(photon.0, pauli)])?)
}
}
}
}
}
}
}
}
outputs.sort_by_key(|p| p.0);
Ok(Self {
dimensions,
program,
boundary,
outputs,
})
}
pub fn bell_strip(
length: usize,
source: NodeId,
analyzer: NodeId,
max_photons: usize,
) -> Result<Self> {
Self::new(
[1, 1, length],
[
PatchBoundary::Measure(Pauli::X),
PatchBoundary::Measure(Pauli::X),
PatchBoundary::Output,
PatchBoundary::Measure(Pauli::Z),
PatchBoundary::Measure(Pauli::Z),
PatchBoundary::Output,
],
source,
analyzer,
max_photons,
)
}
pub fn dimensions(&self) -> [usize; 3] {
self.dimensions
}
pub fn program(&self) -> &FusionProgram {
&self.program
}
pub fn boundary_measurements(&self) -> &[PauliProduct] {
&self.boundary
}
pub fn outputs(&self) -> &[PhotonId] {
&self.outputs
}
pub fn observables(&self) -> Result<Vec<PauliProduct>> {
let mut observables: Vec<_> = self
.program
.fusions()
.iter()
.flat_map(|f| {
[Pauli::X, Pauli::Z]
.map(|p| PauliProduct::new(f.photons.iter().map(|q| (q.0, p)).collect()))
})
.collect::<Result<_>>()?;
observables.extend(self.boundary.clone());
Ok(observables)
}
fn matching_plan(&self, plan: &CompiledFusionNetwork) -> Result<()> {
if plan.program() != &self.program {
return Err(Error::InvalidValue(
"compiled plan differs from finite patch".into(),
));
}
Ok(())
}
pub fn relations(
&self,
plan: &CompiledFusionNetwork,
max_entries: usize,
) -> Result<MeasurementRelations> {
self.matching_plan(plan)?;
let initial = StabilizerSimulator::default().prepare(plan)?;
MeasurementRelations::derive(&initial, &self.observables()?, max_entries)
}
#[allow(clippy::too_many_arguments)]
pub fn sample(
&self,
plan: &CompiledFusionNetwork,
start: u64,
shots: usize,
seed: u64,
pauli: PauliNoise,
bursts: &[PhotonLossBurst],
targets: &[PauliProduct],
max_entries: usize,
) -> Result<Vec<PatchShot>> {
self.matching_plan(plan)?;
start.checked_add(shots as u64).ok_or(Error::Overflow)?;
let fields = self
.program
.num_photons()
.checked_mul(3)
.and_then(|v| v.checked_add(targets.len() * 2))
.ok_or(Error::Overflow)?;
let required = shots.checked_mul(fields).ok_or(Error::Overflow)?;
if required > max_entries {
return Err(Error::ResourceLimit {
resource: "finite patch result entries",
required,
limit: max_entries,
});
}
if bursts
.iter()
.any(|b| b.photons.iter().any(|q| q.0 >= self.program.num_photons()))
{
return Err(Error::Wires("loss burst references absent photon".into()));
}
if targets.iter().any(|p| {
p.terms()
.iter()
.any(|&(q, _)| !self.outputs.iter().any(|p| p.0 == q))
}) {
return Err(Error::Wires(
"patch target is not on retained outputs".into(),
));
}
let initial = StabilizerSimulator::default().prepare(plan)?;
let relations = MeasurementRelations::derive(&initial, &self.observables()?, max_entries)?;
let domain = shot_seed(seed, 0, 0x4f50_454e_4642_0001);
let mut result = Vec::with_capacity(shots);
for index in 0..shots {
let shot_id = start + index as u64;
let mut burst_loss = vec![false; initial.qubits()];
for (i, burst) in bursts.iter().enumerate() {
let mut rng =
ChaCha8Rng::seed_from_u64(shot_seed(domain, shot_id, 0x1000_0000 + i as u64));
if rng.random::<f64>() < burst.probability {
for q in &burst.photons {
burst_loss[q.0] = true;
}
}
}
let (mut state, lost, fusion_outcomes) =
execute_trajectory(plan, &initial, shot_id, domain, pauli, Some(&burst_loss))?;
let mut records = MeasurementRelations::fusion_records(&fusion_outcomes);
for (i, p) in self.boundary.iter().enumerate() {
let mut rng =
ChaCha8Rng::seed_from_u64(shot_seed(domain, shot_id, 0x2000_0000 + i as u64));
let q = p.terms()[0].0;
if lost[q] {
records.push(None);
} else {
records.push(Some(state.measure(p, rng.random())?.negative));
}
}
let mut expectations = Vec::with_capacity(targets.len());
let mut inferred = Vec::with_capacity(targets.len());
for p in targets {
if p.terms().iter().any(|&(q, _)| lost[q]) {
expectations.push(0);
inferred.push(ParityInference::Erased);
} else {
expectations.push(state.expectation(p)?);
inferred.push(relations.infer(&records, p)?);
}
}
result.push(PatchShot {
shot_id,
fusion_outcomes,
records,
output_present: self.outputs.iter().map(|q| !lost[q.0]).collect(),
inferred,
expectations,
});
}
Ok(result)
}
}