quantrs2_device/process_tomography/reconstruction/
bayesian.rs1use scirs2_core::ndarray::{Array1, Array2, Array4};
4use scirs2_core::Complex64;
5
6use super::super::core::SciRS2ProcessTomographer;
7use super::super::results::{ExperimentalData, ReconstructionQuality};
8use super::utils::calculate_reconstruction_quality;
9use crate::DeviceResult;
10
11pub fn reconstruct_bayesian(
13 tomographer: &SciRS2ProcessTomographer,
14 experimental_data: &ExperimentalData,
15) -> DeviceResult<(Array4<Complex64>, ReconstructionQuality)> {
16 let num_qubits = (experimental_data.input_states[0].dim().0 as f64).log2() as usize;
17 let dim = 1 << num_qubits;
18
19 let mut process_matrix = Array4::zeros((dim, dim, dim, dim));
21 for i in 0..dim {
22 process_matrix[[i, i, i, i]] = Complex64::new(1.0, 0.0);
23 }
24
25 let log_likelihood = calculate_bayesian_log_likelihood(&process_matrix, experimental_data)?;
29 let reconstruction_quality =
30 calculate_reconstruction_quality(&process_matrix, experimental_data, log_likelihood);
31
32 Ok((process_matrix, reconstruction_quality))
33}
34
35fn calculate_bayesian_log_likelihood(
37 process_matrix: &Array4<Complex64>,
38 experimental_data: &ExperimentalData,
39) -> DeviceResult<f64> {
40 let mut log_likelihood = 0.0;
41
42 for (observed, &uncertainty) in experimental_data
44 .measurement_results
45 .iter()
46 .zip(experimental_data.measurement_uncertainties.iter())
47 {
48 let predicted = 0.5; let diff = observed - predicted;
50 let variance = uncertainty * uncertainty;
51 log_likelihood -= 0.5 * (diff * diff / variance);
52 }
53
54 let dim = process_matrix.dim().0;
56 for i in 0..dim {
57 for j in 0..dim {
58 for k in 0..dim {
59 for l in 0..dim {
60 if i == j && k == l && i == k {
61 let deviation = (process_matrix[[i, j, k, l]].re - 1.0).abs();
63 log_likelihood -= 0.5 * deviation * deviation;
64 } else {
65 let magnitude = process_matrix[[i, j, k, l]].norm();
67 log_likelihood -= 0.5 * magnitude * magnitude;
68 }
69 }
70 }
71 }
72 }
73
74 Ok(log_likelihood)
75}