use super::{FusionBranch, FusionOutcome};
use crate::{DenseSimulator, Error, Observable, Result};
use burn::tensor::Tensor;
#[derive(Clone, Debug)]
pub struct HeraldedReadout {
probability: Tensor<1>,
numerator: Tensor<1>,
}
impl HeraldedReadout {
pub fn probability(&self) -> Tensor<1> {
self.probability.clone()
}
pub fn unconditional(&self) -> Tensor<1> {
self.numerator.clone()
}
pub fn conditional(&self, minimum_probability: f64) -> Result<Tensor<1>> {
crate::error::tolerance(minimum_probability)?;
if crate::complex::read(self.probability())?
.iter()
.any(|p| !p.is_finite() || *p <= minimum_probability)
{
return Err(Error::InvalidValue(
"herald probability below conditioning threshold".into(),
));
}
Ok(self.unconditional() / self.probability())
}
}
pub fn heralded_readout(
branches: &[FusionBranch],
engine: &DenseSimulator,
observable: &Observable,
accept: impl Fn(&[FusionOutcome]) -> bool,
) -> Result<HeraldedReadout> {
let first = branches
.first()
.ok_or_else(|| Error::InvalidValue("empty instrument".into()))?;
let mut probability = first.probability() * 0.;
let mut numerator = probability.clone();
for branch in branches {
let p = branch.probability();
if p.dims() != probability.dims() {
return Err(Error::Shape("inconsistent branch batches".into()));
}
crate::complex::compatible(&p, &probability)?;
if accept(branch.outcomes()) {
let rho = branch.unnormalized_state().ok_or_else(|| {
Error::InvalidValue("accepted branch has no quantum output".into())
})?;
numerator = numerator + rho.expectation(engine, observable)?;
probability = probability + p;
}
}
Ok(HeraldedReadout {
probability,
numerator,
})
}