fusion_blossom/
mwpm_solver.rs

1//! Minimum-Weight Perfect Matching Solver
2//!
3//! This module includes some common usage of primal and dual modules to solve MWPM problem.
4//! Note that you can call different primal and dual modules, even interchangeably, by following the examples in this file
5//!
6
7use std::collections::{BTreeMap, BTreeSet};
8use std::fs::File;
9use std::io::prelude::*;
10use std::io::BufWriter;
11
12use nonzero::nonzero as nz;
13#[cfg(feature = "python_binding")]
14use pyo3::prelude::*;
15
16use crate::blossom_v;
17use crate::complete_graph::*;
18use crate::derivative::Derivative;
19use crate::dual_module::*;
20
21use super::dual_module::{DualModuleImpl, DualModuleInterfacePtr};
22use super::dual_module_parallel::*;
23use super::dual_module_serial::DualModuleSerial;
24use super::pointers::*;
25use super::primal_module::{PerfectMatching, PrimalModuleImpl, SubGraphBuilder, VisualizeSubgraph};
26use super::primal_module_parallel::*;
27use super::primal_module_serial::PrimalModuleSerialPtr;
28use super::util::*;
29use super::visualize::*;
30
31/// a serial solver
32#[derive(Derivative)]
33#[derivative(Debug)]
34#[cfg_attr(feature = "python_binding", cfg_eval)]
35#[cfg_attr(feature = "python_binding", pyclass)]
36pub struct LegacySolverSerial {
37    #[cfg_attr(feature = "python_binding", pyo3(get, set))]
38    initializer: SolverInitializer,
39    /// a serial implementation of the primal module
40    #[derivative(Debug = "ignore")]
41    primal_module: PrimalModuleSerialPtr,
42    /// a serial implementation of the dual module
43    #[derivative(Debug = "ignore")]
44    dual_module: DualModuleSerial,
45    /// the interface between the primal and dual module
46    interface_ptr: DualModuleInterfacePtr,
47    /// subgraph builder for easier integration with decoder
48    subgraph_builder: SubGraphBuilder,
49}
50
51impl Clone for LegacySolverSerial {
52    fn clone(&self) -> Self {
53        Self::new(&self.initializer) // create independent instances of the solver
54    }
55}
56
57impl LegacySolverSerial {
58    /// create a new decoder
59    pub fn new(initializer: &SolverInitializer) -> Self {
60        let dual_module = DualModuleSerial::new_empty(initializer);
61        let primal_module = PrimalModuleSerialPtr::new_empty(initializer);
62        let interface_ptr = DualModuleInterfacePtr::new_empty();
63        let subgraph_builder = SubGraphBuilder::new(initializer);
64        Self {
65            initializer: initializer.clone(),
66            primal_module,
67            dual_module,
68            interface_ptr,
69            subgraph_builder,
70        }
71    }
72
73    pub fn solve_perfect_matching(
74        &mut self,
75        syndrome_pattern: &SyndromePattern,
76        visualizer: Option<&mut Visualizer>,
77    ) -> PerfectMatching {
78        self.primal_module.clear();
79        self.dual_module.clear();
80        self.interface_ptr.clear();
81        self.primal_module
82            .solve_visualizer(&self.interface_ptr, syndrome_pattern, &mut self.dual_module, visualizer);
83        self.primal_module
84            .perfect_matching(&self.interface_ptr, &mut self.dual_module)
85    }
86
87    /// solve subgraph directly
88    pub fn solve_subgraph(&mut self, syndrome_pattern: &SyndromePattern) -> Vec<EdgeIndex> {
89        self.solve_subgraph_visualizer(syndrome_pattern, None)
90    }
91
92    pub fn solve_subgraph_visualizer(
93        &mut self,
94        syndrome_pattern: &SyndromePattern,
95        visualizer: Option<&mut Visualizer>,
96    ) -> Vec<EdgeIndex> {
97        let perfect_matching = self.solve_perfect_matching(syndrome_pattern, visualizer);
98        self.subgraph_builder.clear();
99        self.subgraph_builder.load_perfect_matching(&perfect_matching);
100        self.subgraph_builder.get_subgraph()
101    }
102
103    /// solve the minimum weight perfect matching (legacy API, the same output as the blossom V library)
104    pub fn solve_legacy(&mut self, syndrome_pattern: &SyndromePattern) -> Vec<VertexIndex> {
105        self.solve_legacy_visualizer(syndrome_pattern, None)
106    }
107
108    pub fn solve_legacy_visualizer(
109        &mut self,
110        syndrome_pattern: &SyndromePattern,
111        visualizer: Option<&mut Visualizer>,
112    ) -> Vec<VertexIndex> {
113        self.primal_module.clear();
114        self.dual_module.clear();
115        self.interface_ptr.clear();
116        self.primal_module
117            .solve_visualizer(&self.interface_ptr, syndrome_pattern, &mut self.dual_module, visualizer);
118        let perfect_matching = self
119            .primal_module
120            .perfect_matching(&self.interface_ptr, &mut self.dual_module);
121        perfect_matching.legacy_get_mwpm_result(syndrome_pattern.defect_vertices.clone())
122    }
123}
124
125// static functions, not recommended because it doesn't reuse the data structure of dual module
126impl LegacySolverSerial {
127    pub fn mwpm_solve(initializer: &SolverInitializer, syndrome_pattern: &SyndromePattern) -> Vec<VertexIndex> {
128        Self::mwpm_solve_visualizer(initializer, syndrome_pattern, None)
129    }
130
131    pub fn mwpm_solve_visualizer(
132        initializer: &SolverInitializer,
133        syndrome_pattern: &SyndromePattern,
134        visualizer: Option<&mut Visualizer>,
135    ) -> Vec<VertexIndex> {
136        let mut solver = Self::new(initializer);
137        solver.solve_legacy_visualizer(syndrome_pattern, visualizer)
138    }
139}
140
141pub trait PrimalDualSolver {
142    fn clear(&mut self);
143    fn reset_profiler(&mut self) {} // only if profiler records some information that needs to be cleared, e.g. vec![]
144    fn solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>);
145    fn solve(&mut self, syndrome_pattern: &SyndromePattern) {
146        self.solve_visualizer(syndrome_pattern, None)
147    }
148    fn perfect_matching_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching;
149    fn perfect_matching(&mut self) -> PerfectMatching {
150        self.perfect_matching_visualizer(None)
151    }
152    fn subgraph_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex>;
153    fn subgraph(&mut self) -> Vec<EdgeIndex> {
154        self.subgraph_visualizer(None)
155    }
156    fn sum_dual_variables(&self) -> Weight;
157    fn generate_profiler_report(&self) -> serde_json::Value;
158    #[allow(clippy::unnecessary_cast)]
159    fn stim_integration_predict_bit_packed_data(
160        &mut self,
161        in_file: String,
162        out_file: String,
163        edge_masks: &[usize],
164        num_shots: usize,
165        num_dets: usize,
166        num_obs: usize,
167    ) {
168        let mut in_reader = std::io::BufReader::new(File::open(in_file).expect("in_file not found"));
169        let mut out_writer = std::io::BufWriter::new(File::create(out_file).expect("out_file not found"));
170        let num_det_bytes = (num_dets + 7) / 8; // ceil
171        let mut dets_bit_packed = vec![0; num_det_bytes];
172        assert!(num_obs <= 64, "too many observables");
173        let prediction_bytes = (num_obs + 7) / 8; // ceil
174        for _ in 0..num_shots {
175            in_reader.read_exact(&mut dets_bit_packed).expect("read success");
176            let mut defect_vertices = vec![];
177            for (i, &byte) in dets_bit_packed.iter().enumerate() {
178                if byte == 0 {
179                    continue;
180                }
181                for j in 0..8 {
182                    if byte & (1 << j) != 0 {
183                        // little endian
184                        defect_vertices.push((i * 8 + j) as VertexIndex);
185                    }
186                }
187            }
188            let syndrome_pattern = SyndromePattern::new_vertices(defect_vertices);
189            self.solve(&syndrome_pattern);
190            let subgraph = self.subgraph();
191            let mut prediction = 0;
192            for edge_index in subgraph {
193                prediction ^= edge_masks[edge_index as usize];
194            }
195            for j in 0..prediction_bytes {
196                let byte = ((prediction >> (j * 8)) & 0x0FF) as u8;
197                out_writer.write_all(&[byte]).unwrap();
198            }
199            self.clear();
200        }
201    }
202}
203
204#[cfg(feature = "python_binding")]
205macro_rules! bind_trait_primal_dual_solver {
206    ($struct_name:ident) => {
207        #[pymethods]
208        impl $struct_name {
209            #[pyo3(name = "clear")]
210            fn trait_clear(&mut self) {
211                self.clear()
212            }
213            #[pyo3(name = "solve_visualizer")]
214            fn trait_solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>) {
215                self.solve_visualizer(syndrome_pattern, visualizer)
216            }
217            #[pyo3(name = "solve")] // in Python, `solve` and `solve_visualizer` is the same because it can take optional parameter
218            fn trait_solve(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>) {
219                self.solve_visualizer(syndrome_pattern, visualizer)
220            }
221            #[pyo3(name = "perfect_matching_visualizer")]
222            fn trait_perfect_matching_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching {
223                self.perfect_matching_visualizer(visualizer)
224            }
225            #[pyo3(name = "perfect_matching")] // in Python, `perfect_matching` and `perfect_matching_visualizer` is the same because it can take optional parameter
226            fn trait_perfect_matching(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching {
227                self.perfect_matching_visualizer(visualizer)
228            }
229            #[pyo3(name = "subgraph_visualizer")]
230            fn trait_subgraph_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
231                self.subgraph_visualizer(visualizer)
232            }
233            #[pyo3(name = "subgraph")] // in Python, `subgraph` and `subgraph_visualizer` is the same because it can take optional parameter
234            fn trait_subgraph(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
235                self.subgraph_visualizer(visualizer)
236            }
237            #[pyo3(name = "sum_dual_variables")]
238            fn trait_sum_dual_variables(&self) -> Weight {
239                self.sum_dual_variables()
240            }
241            #[pyo3(name = "generate_profiler_report")]
242            fn trait_generate_profiler_report(&self) -> PyObject {
243                json_to_pyobject(self.generate_profiler_report())
244            }
245            #[pyo3(name = "stim_integration_predict_bit_packed_data")]
246            fn trait_stim_integration_predict_bit_packed_data(
247                &mut self,
248                in_file: String,
249                out_file: String,
250                edge_masks: Vec<usize>,
251                num_shots: usize,
252                num_dets: usize,
253                num_obs: usize,
254            ) {
255                self.stim_integration_predict_bit_packed_data(in_file, out_file, &edge_masks, num_shots, num_dets, num_obs)
256            }
257        }
258    };
259}
260
261#[cfg_attr(feature = "python_binding", cfg_eval)]
262#[cfg_attr(feature = "python_binding", pyclass)]
263pub struct SolverSerial {
264    pub dual_module: DualModuleSerial,
265    pub primal_module: PrimalModuleSerialPtr,
266    pub interface_ptr: DualModuleInterfacePtr,
267    pub subgraph_builder: SubGraphBuilder,
268}
269
270bind_trait_fusion_visualizer!(SolverSerial);
271impl FusionVisualizer for SolverSerial {
272    fn snapshot(&self, abbrev: bool) -> serde_json::Value {
273        let mut value = self.primal_module.snapshot(abbrev);
274        snapshot_combine_values(&mut value, self.dual_module.snapshot(abbrev), abbrev);
275        snapshot_combine_values(&mut value, self.interface_ptr.snapshot(abbrev), abbrev);
276        value
277    }
278}
279
280#[cfg(feature = "python_binding")]
281bind_trait_primal_dual_solver! {SolverSerial}
282
283#[cfg(feature = "python_binding")]
284#[pymethods]
285impl SolverSerial {
286    #[new]
287    #[pyo3(signature = (initializer, *, max_tree_size = None))]
288    pub fn new_python(initializer: &SolverInitializer, max_tree_size: Option<usize>) -> Self {
289        let mut solver = Self::new(initializer);
290        if let Some(max_tree_size) = max_tree_size {
291            solver.primal_module.write().max_tree_size = max_tree_size;
292        }
293        solver
294    }
295}
296
297impl SolverSerial {
298    pub fn new(initializer: &SolverInitializer) -> Self {
299        Self {
300            dual_module: DualModuleSerial::new_empty(initializer),
301            primal_module: PrimalModuleSerialPtr::new_empty(initializer),
302            interface_ptr: DualModuleInterfacePtr::new_empty(),
303            subgraph_builder: SubGraphBuilder::new(initializer),
304        }
305    }
306}
307
308impl PrimalDualSolver for SolverSerial {
309    fn clear(&mut self) {
310        self.primal_module.clear();
311        self.dual_module.clear();
312        self.interface_ptr.clear();
313        self.subgraph_builder.clear();
314    }
315    fn solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>) {
316        if !syndrome_pattern.erasures.is_empty() {
317            assert!(
318                syndrome_pattern.dynamic_weights.is_empty(),
319                "erasures and dynamic_weights cannot be provided at the same time"
320            );
321            self.subgraph_builder.load_erasures(&syndrome_pattern.erasures);
322        }
323        if !syndrome_pattern.dynamic_weights.is_empty() {
324            self.subgraph_builder.load_dynamic_weights(&syndrome_pattern.dynamic_weights);
325        }
326        self.primal_module
327            .solve_visualizer(&self.interface_ptr, syndrome_pattern, &mut self.dual_module, visualizer);
328    }
329    fn perfect_matching_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching {
330        let perfect_matching = self
331            .primal_module
332            .perfect_matching(&self.interface_ptr, &mut self.dual_module);
333        if let Some(visualizer) = visualizer {
334            visualizer
335                .snapshot_combined(
336                    "perfect matching".to_string(),
337                    vec![&self.interface_ptr, &self.dual_module, &perfect_matching],
338                )
339                .unwrap();
340        }
341        perfect_matching
342    }
343    fn subgraph_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
344        let perfect_matching = self.perfect_matching();
345        self.subgraph_builder.load_perfect_matching(&perfect_matching);
346        let subgraph = self.subgraph_builder.get_subgraph();
347        if let Some(visualizer) = visualizer {
348            visualizer
349                .snapshot_combined(
350                    "perfect matching and subgraph".to_string(),
351                    vec![
352                        &self.interface_ptr,
353                        &self.dual_module,
354                        &perfect_matching,
355                        &VisualizeSubgraph::new(&subgraph),
356                    ],
357                )
358                .unwrap();
359        }
360        subgraph
361    }
362    fn sum_dual_variables(&self) -> Weight {
363        self.interface_ptr.read_recursive().sum_dual_variables
364    }
365    fn generate_profiler_report(&self) -> serde_json::Value {
366        json!({
367            "dual": self.dual_module.generate_profiler_report(),
368            "primal": self.primal_module.generate_profiler_report(),
369        })
370    }
371}
372
373#[cfg_attr(feature = "python_binding", cfg_eval)]
374#[cfg_attr(feature = "python_binding", pyclass)]
375pub struct SolverDualParallel {
376    pub dual_module: DualModuleParallel<DualModuleSerial>,
377    pub primal_module: PrimalModuleSerialPtr,
378    pub interface_ptr: DualModuleInterfacePtr,
379    pub subgraph_builder: SubGraphBuilder,
380}
381
382bind_trait_fusion_visualizer!(SolverDualParallel);
383impl FusionVisualizer for SolverDualParallel {
384    fn snapshot(&self, abbrev: bool) -> serde_json::Value {
385        let mut value = self.primal_module.snapshot(abbrev);
386        snapshot_combine_values(&mut value, self.dual_module.snapshot(abbrev), abbrev);
387        snapshot_combine_values(&mut value, self.interface_ptr.snapshot(abbrev), abbrev);
388        value
389    }
390}
391
392#[cfg(feature = "python_binding")]
393bind_trait_primal_dual_solver! {SolverDualParallel}
394
395#[cfg(feature = "python_binding")]
396#[pymethods]
397impl SolverDualParallel {
398    #[new]
399    pub fn new_python(
400        initializer: &SolverInitializer,
401        partition_info: &PartitionInfo,
402        primal_dual_config: PyObject,
403    ) -> Self {
404        let primal_dual_config = pyobject_to_json(primal_dual_config);
405        Self::new(initializer, partition_info, primal_dual_config)
406    }
407}
408
409impl SolverDualParallel {
410    pub fn new(
411        initializer: &SolverInitializer,
412        partition_info: &PartitionInfo,
413        primal_dual_config: serde_json::Value,
414    ) -> Self {
415        let config: DualModuleParallelConfig = serde_json::from_value(primal_dual_config).unwrap();
416        Self {
417            dual_module: DualModuleParallel::new_config(initializer, partition_info, config),
418            primal_module: PrimalModuleSerialPtr::new_empty(initializer),
419            interface_ptr: DualModuleInterfacePtr::new_empty(),
420            subgraph_builder: SubGraphBuilder::new(initializer),
421        }
422    }
423}
424
425impl PrimalDualSolver for SolverDualParallel {
426    fn clear(&mut self) {
427        self.dual_module.clear();
428        self.primal_module.clear();
429        self.interface_ptr.clear();
430        self.subgraph_builder.clear();
431    }
432    fn solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>) {
433        if !syndrome_pattern.erasures.is_empty() {
434            assert!(
435                syndrome_pattern.dynamic_weights.is_empty(),
436                "erasures and dynamic_weights cannot be provided at the same time"
437            );
438            self.subgraph_builder.load_erasures(&syndrome_pattern.erasures);
439        }
440        if !syndrome_pattern.dynamic_weights.is_empty() {
441            self.subgraph_builder.load_dynamic_weights(&syndrome_pattern.dynamic_weights);
442        }
443        self.dual_module.static_fuse_all();
444        self.primal_module
445            .solve_visualizer(&self.interface_ptr, syndrome_pattern, &mut self.dual_module, visualizer);
446    }
447    fn perfect_matching_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching {
448        let perfect_matching = self
449            .primal_module
450            .perfect_matching(&self.interface_ptr, &mut self.dual_module);
451        if let Some(visualizer) = visualizer {
452            visualizer
453                .snapshot_combined(
454                    "perfect matching".to_string(),
455                    vec![&self.interface_ptr, &self.dual_module, &perfect_matching],
456                )
457                .unwrap();
458        }
459        perfect_matching
460    }
461    fn subgraph_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
462        let perfect_matching = self.perfect_matching();
463        self.subgraph_builder.load_perfect_matching(&perfect_matching);
464        let subgraph = self.subgraph_builder.get_subgraph();
465        if let Some(visualizer) = visualizer {
466            visualizer
467                .snapshot_combined(
468                    "perfect matching and subgraph".to_string(),
469                    vec![
470                        &self.interface_ptr,
471                        &self.dual_module,
472                        &perfect_matching,
473                        &VisualizeSubgraph::new(&subgraph),
474                    ],
475                )
476                .unwrap();
477        }
478        subgraph
479    }
480    fn sum_dual_variables(&self) -> Weight {
481        self.interface_ptr.read_recursive().sum_dual_variables
482    }
483    fn generate_profiler_report(&self) -> serde_json::Value {
484        json!({
485            "dual": self.dual_module.generate_profiler_report(),
486            "primal": self.primal_module.generate_profiler_report(),
487        })
488    }
489}
490
491#[cfg_attr(feature = "python_binding", cfg_eval)]
492#[cfg_attr(feature = "python_binding", pyclass)]
493pub struct SolverParallel {
494    pub dual_module: DualModuleParallel<DualModuleSerial>,
495    pub primal_module: PrimalModuleParallel,
496    pub subgraph_builder: SubGraphBuilder,
497}
498
499bind_trait_fusion_visualizer!(SolverParallel);
500impl FusionVisualizer for SolverParallel {
501    fn snapshot(&self, abbrev: bool) -> serde_json::Value {
502        let mut value = self.primal_module.snapshot(abbrev);
503        snapshot_combine_values(&mut value, self.dual_module.snapshot(abbrev), abbrev);
504        value
505    }
506}
507
508#[cfg(feature = "python_binding")]
509bind_trait_primal_dual_solver! {SolverParallel}
510
511#[cfg(feature = "python_binding")]
512#[pymethods]
513impl SolverParallel {
514    #[new]
515    pub fn new_python(
516        initializer: &SolverInitializer,
517        partition_info: &PartitionInfo,
518        primal_dual_config: PyObject,
519    ) -> Self {
520        let primal_dual_config = pyobject_to_json(primal_dual_config);
521        Self::new(initializer, partition_info, primal_dual_config)
522    }
523
524    #[pyo3(name = "defect_perfect_matching")]
525    pub fn defect_perfect_matching(&mut self) -> Vec<(VertexIndex, VertexIndex)> {
526        let perfect_matching = self.perfect_matching_visualizer(None);
527        let mut defect_matching = vec![];
528        // iterate over peer matching
529        for (a, b) in perfect_matching.peer_matchings.iter() {
530            let node_a = a.read_recursive();
531            let vertex_a = if let DualNodeClass::DefectVertex { defect_index } = &node_a.class {
532                *defect_index
533            } else {
534                unreachable!("can only be syndrome")
535            };
536            let node_b = b.read_recursive();
537            let vertex_b = if let DualNodeClass::DefectVertex { defect_index } = &node_b.class {
538                *defect_index
539            } else {
540                unreachable!("can only be syndrome")
541            };
542            defect_matching.push((vertex_a, vertex_b));
543        }
544        // iterate over virtual matching
545        for (a, virtual_vertex) in perfect_matching.virtual_matchings.iter() {
546            let node_a = a.read_recursive();
547            let vertex_a = if let DualNodeClass::DefectVertex { defect_index } = &node_a.class {
548                *defect_index
549            } else {
550                unreachable!("can only be syndrome")
551            };
552            defect_matching.push((vertex_a, *virtual_vertex));
553        }
554        defect_matching
555    }
556}
557
558impl SolverParallel {
559    pub fn new(
560        initializer: &SolverInitializer,
561        partition_info: &PartitionInfo,
562        mut primal_dual_config: serde_json::Value,
563    ) -> Self {
564        let primal_dual_config = primal_dual_config.as_object_mut().expect("config must be JSON object");
565        let mut dual_config = DualModuleParallelConfig::default();
566        let mut primal_config = PrimalModuleParallelConfig::default();
567        if let Some(value) = primal_dual_config.remove("dual") {
568            dual_config = serde_json::from_value(value).unwrap();
569        }
570        if let Some(value) = primal_dual_config.remove("primal") {
571            primal_config = serde_json::from_value(value).unwrap();
572        }
573        if !primal_dual_config.is_empty() {
574            panic!(
575                "unknown primal_dual_config keys: {:?}",
576                primal_dual_config.keys().collect::<Vec<&String>>()
577            );
578        }
579        Self {
580            dual_module: DualModuleParallel::new_config(initializer, partition_info, dual_config),
581            primal_module: PrimalModuleParallel::new_config(initializer, partition_info, primal_config),
582            subgraph_builder: SubGraphBuilder::new(initializer),
583        }
584    }
585}
586
587impl PrimalDualSolver for SolverParallel {
588    fn clear(&mut self) {
589        self.dual_module.clear();
590        self.primal_module.clear();
591        self.subgraph_builder.clear();
592    }
593    fn solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>) {
594        if !syndrome_pattern.erasures.is_empty() {
595            self.subgraph_builder.load_erasures(&syndrome_pattern.erasures);
596        }
597        self.primal_module
598            .parallel_solve_visualizer(syndrome_pattern, &self.dual_module, visualizer);
599    }
600    fn perfect_matching_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching {
601        let useless_interface_ptr = DualModuleInterfacePtr::new_empty(); // don't actually use it
602        let perfect_matching = self
603            .primal_module
604            .perfect_matching(&useless_interface_ptr, &mut self.dual_module);
605        if let Some(visualizer) = visualizer {
606            let last_interface_ptr = &self.primal_module.units.last().unwrap().read_recursive().interface_ptr;
607            visualizer
608                .snapshot_combined(
609                    "perfect matching".to_string(),
610                    vec![last_interface_ptr, &self.dual_module, &perfect_matching],
611                )
612                .unwrap();
613        }
614        perfect_matching
615    }
616    fn subgraph_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
617        let perfect_matching = self.perfect_matching();
618        self.subgraph_builder.load_perfect_matching(&perfect_matching);
619        let subgraph = self.subgraph_builder.get_subgraph();
620        if let Some(visualizer) = visualizer {
621            let last_interface_ptr = &self.primal_module.units.last().unwrap().read_recursive().interface_ptr;
622            visualizer
623                .snapshot_combined(
624                    "perfect matching and subgraph".to_string(),
625                    vec![
626                        last_interface_ptr,
627                        &self.dual_module,
628                        &perfect_matching,
629                        &VisualizeSubgraph::new(&subgraph),
630                    ],
631                )
632                .unwrap();
633        }
634        subgraph
635    }
636    fn sum_dual_variables(&self) -> Weight {
637        let last_unit = self.primal_module.units.last().unwrap().write(); // use the interface in the last unit
638        let sum_dual_variables = last_unit.interface_ptr.read_recursive().sum_dual_variables;
639        sum_dual_variables
640    }
641    fn generate_profiler_report(&self) -> serde_json::Value {
642        json!({
643            "dual": self.dual_module.generate_profiler_report(),
644            "primal": self.primal_module.generate_profiler_report(),
645        })
646    }
647}
648
649#[cfg_attr(feature = "python_binding", cfg_eval)]
650#[cfg_attr(feature = "python_binding", pyclass)]
651pub struct SolverErrorPatternLogger {
652    pub file: BufWriter<File>,
653}
654
655#[cfg(feature = "python_binding")]
656bind_trait_primal_dual_solver! {SolverErrorPatternLogger}
657
658impl SolverErrorPatternLogger {
659    pub fn new(initializer: &SolverInitializer, positions: &Vec<VisualizePosition>, mut config: serde_json::Value) -> Self {
660        let mut filename = "tmp/syndrome_patterns.txt".to_string();
661        let config = config.as_object_mut().expect("config must be JSON object");
662        if let Some(value) = config.remove("filename") {
663            filename = value.as_str().expect("filename string").to_string();
664        }
665        if !config.is_empty() {
666            panic!("unknown config keys: {:?}", config.keys().collect::<Vec<&String>>());
667        }
668        let file = File::create(filename).unwrap();
669        let mut file = BufWriter::new(file);
670        file.write_all(b"Syndrome Pattern v1.0   <initializer> <positions> <syndrome_pattern>*\n")
671            .unwrap();
672        serde_json::to_writer(&mut file, &initializer).unwrap(); // large object write to file directly
673        file.write_all(b"\n").unwrap();
674        serde_json::to_writer(&mut file, &positions).unwrap();
675        file.write_all(b"\n").unwrap();
676        Self { file }
677    }
678}
679
680impl PrimalDualSolver for SolverErrorPatternLogger {
681    fn clear(&mut self) {}
682    fn solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, _visualizer: Option<&mut Visualizer>) {
683        self.file
684            .write_all(
685                serde_json::to_string(&serde_json::json!(syndrome_pattern))
686                    .unwrap()
687                    .as_bytes(),
688            )
689            .unwrap();
690        self.file.write_all(b"\n").unwrap();
691    }
692    fn perfect_matching_visualizer(&mut self, _visualizer: Option<&mut Visualizer>) -> PerfectMatching {
693        panic!("error pattern logger do not actually solve the problem, please use Verifier::None by `--verifier none`")
694    }
695    fn subgraph_visualizer(&mut self, _visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
696        // panic!("error pattern logger do not actually solve the problem, please use Verifier::None by `--verifier none`")
697        vec![]
698    }
699    fn sum_dual_variables(&self) -> Weight {
700        panic!("error pattern logger do not actually solve the problem")
701    }
702    fn generate_profiler_report(&self) -> serde_json::Value {
703        json!({})
704    }
705}
706
707/// an exact solver calling blossom V library for benchmarking comparison
708#[derive(Clone)]
709pub struct SolverBlossomV {
710    pub initializer: SolverInitializer,
711    pub prebuilt_complete_graph: PrebuiltCompleteGraph,
712    pub subgraph_builder: SubGraphBuilder,
713    pub matched_pairs: Vec<(VertexIndex, VertexIndex)>,
714}
715
716impl SolverBlossomV {
717    pub fn new(initializer: &SolverInitializer) -> Self {
718        Self {
719            initializer: initializer.clone(),
720            prebuilt_complete_graph: PrebuiltCompleteGraph::new_threaded(initializer, 0),
721            subgraph_builder: SubGraphBuilder::new(initializer),
722            matched_pairs: vec![],
723        }
724    }
725}
726
727impl PrimalDualSolver for SolverBlossomV {
728    fn clear(&mut self) {
729        self.matched_pairs.clear();
730        self.subgraph_builder.clear();
731    }
732    fn solve_visualizer(&mut self, syndrome_pattern: &SyndromePattern, visualizer: Option<&mut Visualizer>) {
733        assert!(visualizer.is_none(), "not supported");
734        assert!(syndrome_pattern.erasures.is_empty(), "doesn't support erasure for now");
735        let defect_vertices = &syndrome_pattern.defect_vertices;
736        if defect_vertices.is_empty() {
737            return;
738        }
739        let mut mapping_to_defect_vertices: BTreeMap<VertexIndex, usize> = BTreeMap::new();
740        for (i, &defect_vertex) in defect_vertices.iter().enumerate() {
741            mapping_to_defect_vertices.insert(defect_vertex, i);
742        }
743        // for each real vertex, add a corresponding virtual vertex to be matched
744        let defect_num = defect_vertices.len();
745        let legacy_vertex_num = defect_num * 2;
746        let mut legacy_weighted_edges = Vec::<(usize, usize, u32)>::new();
747        for i in 0..defect_num - 1 {
748            for j in i + 1..defect_num {
749                if let Some(weight) = self
750                    .prebuilt_complete_graph
751                    .get_edge_weight(defect_vertices[i], defect_vertices[j])
752                {
753                    legacy_weighted_edges.push((i, j, weight as u32));
754                }
755            }
756        }
757        for (i, &defect_vertex) in defect_vertices.iter().enumerate() {
758            if let Some((_, weight)) = self.prebuilt_complete_graph.get_boundary_weight(defect_vertex) {
759                // connect this real vertex to it's corresponding virtual vertex
760                legacy_weighted_edges.push((i, i + defect_num, weight as u32));
761            }
762        }
763        for i in 0..defect_num - 1 {
764            for j in i + 1..defect_num {
765                // virtual boundaries are always fully connected with weight 0
766                legacy_weighted_edges.push((i + defect_num, j + defect_num, 0));
767            }
768        }
769        // run blossom V to get matchings
770        // println!("[debug] legacy_vertex_num: {:?}", legacy_vertex_num);
771        // println!("[debug] legacy_weighted_edges: {:?}", legacy_weighted_edges);
772        let matchings = blossom_v::safe_minimum_weight_perfect_matching(legacy_vertex_num, &legacy_weighted_edges);
773        let mut matched_pairs = Vec::new();
774        for i in 0..defect_num {
775            let j = matchings[i];
776            if j < defect_num {
777                // match to a real vertex
778                if i < j {
779                    // avoid duplicate matched pair
780                    matched_pairs.push((defect_vertices[i], defect_vertices[j]));
781                }
782            } else {
783                assert_eq!(
784                    j,
785                    i + defect_num,
786                    "if not matched to another real vertex, it must match to it's corresponding virtual vertex"
787                );
788                matched_pairs.push((
789                    defect_vertices[i],
790                    self.prebuilt_complete_graph
791                        .get_boundary_weight(defect_vertices[i])
792                        .expect("boundary must exist if match to virtual vertex")
793                        .0,
794                ));
795            }
796        }
797        self.matched_pairs = matched_pairs;
798    }
799    fn perfect_matching_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> PerfectMatching {
800        assert!(visualizer.is_none(), "not supported");
801        let virtual_vertices: BTreeSet<VertexIndex> = self.initializer.virtual_vertices.iter().cloned().collect();
802        let mut perfect_matching = PerfectMatching::new();
803        let mut counter = 0;
804        let interface_ptr = DualModuleInterfacePtr::new_empty();
805        let mut create_dual_node = |vertex_index: VertexIndex| {
806            counter += 1;
807            DualNodePtr::new_value(DualNode {
808                index: counter,
809                class: DualNodeClass::DefectVertex {
810                    defect_index: vertex_index,
811                },
812                grow_state: DualNodeGrowState::Grow,
813                parent_blossom: None,
814                dual_variable_cache: (0, 0),
815                belonging: interface_ptr.downgrade(),
816                defect_size: nz!(1usize),
817            })
818        };
819        for &(vertex_1, vertex_2) in self.matched_pairs.iter() {
820            assert!(!virtual_vertices.contains(&vertex_1)); // 1 is not virtual
821            if virtual_vertices.contains(&vertex_2) {
822                perfect_matching
823                    .virtual_matchings
824                    .push((create_dual_node(vertex_1), vertex_2));
825            } else {
826                perfect_matching
827                    .peer_matchings
828                    .push((create_dual_node(vertex_1), create_dual_node(vertex_2)));
829            }
830            self.subgraph_builder.add_matching(vertex_1, vertex_2);
831        }
832        perfect_matching
833    }
834    fn subgraph_visualizer(&mut self, visualizer: Option<&mut Visualizer>) -> Vec<EdgeIndex> {
835        assert!(visualizer.is_none(), "not supported");
836        self.subgraph_builder.clear();
837        for &(vertex_1, vertex_2) in self.matched_pairs.iter() {
838            self.subgraph_builder.add_matching(vertex_1, vertex_2);
839        }
840        self.subgraph_builder.subgraph.iter().copied().collect()
841    }
842    #[allow(clippy::unnecessary_cast)]
843    fn sum_dual_variables(&self) -> Weight {
844        let mut subgraph_builder = self.subgraph_builder.clone();
845        subgraph_builder.clear();
846        for &(vertex_1, vertex_2) in self.matched_pairs.iter() {
847            subgraph_builder.add_matching(vertex_1, vertex_2);
848        }
849        let subgraph: Vec<EdgeIndex> = subgraph_builder.subgraph.iter().copied().collect();
850        let mut weight = 0;
851        for &edge_index in subgraph.iter() {
852            weight += self.initializer.weighted_edges[edge_index as usize].2;
853        }
854        weight
855    }
856    fn generate_profiler_report(&self) -> serde_json::Value {
857        json!({})
858    }
859}
860
861#[cfg(feature = "python_binding")]
862#[pyfunction]
863pub(crate) fn register(_py: Python<'_>, m: &PyModule) -> PyResult<()> {
864    m.add_class::<LegacySolverSerial>()?;
865    m.add_class::<SolverSerial>()?;
866    m.add_class::<SolverDualParallel>()?;
867    m.add_class::<SolverParallel>()?;
868    m.add_class::<SolverErrorPatternLogger>()?;
869    Ok(())
870}
fusion_blossom/mwpm_solver.rs

fusion_blossom/
mwpm_solver.rs
