Skip to content

Commit

Permalink
WIP
Browse files Browse the repository at this point in the history
  • Loading branch information
@chemix-lunacy
chemix-lunacy committed Aug 27, 2024
1 parent 77d3764 commit 0738595
Show file tree
Hide file tree
Showing 2 changed files with 208 additions and 41 deletions.
243 changes: 202 additions & 41 deletions src/rasqal/src/analysis.rs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@ use crate::{with_mutable, with_mutable_self};
use log::{log, Level};
use ndarray::{array, Array2};
use num::range;
use num_complex::Complex;
use num_complex::{Complex, Complex64, ComplexFloat};
use std::cmp::Ordering;
use std::collections::hash_map::Keys;
use std::collections::{HashMap, HashSet};
Expand Down Expand Up @@ -139,7 +139,7 @@ impl AnalysisQubit {
/// Returns 0...x depending upon what this qubit would be measured as.
pub fn measure(&self) -> MeasureAnalysis {
// TODO: Add entanglement information for clusters.
MeasureAnalysis::qubit(self.state.matrix.get((1, 1)).unwrap().re)
MeasureAnalysis::qubit(self.state.get((1, 1)).re)
}

/// Applies this gate to this qubit and all tangles.
Expand Down Expand Up @@ -177,7 +177,7 @@ impl AnalysisQubit {
self.index,
self.measure().result
));
for matrix_fragment in stringified_matrix_lines(&self.state.matrix) {
for matrix_fragment in stringified_matrix_lines(&self.state.matrix_fragment.matrix) {
result.push(format!("{}{}\n", indent, matrix_fragment));
}

Expand All @@ -188,7 +188,7 @@ impl AnalysisQubit {
"{}<{}~{}>:\n",
indent, state.upper_left, state.bottom_right
));
for matrix_fragment in stringified_matrix_lines(&state.state.matrix) {
for matrix_fragment in stringified_matrix_lines(&state.state.matrix_fragment.matrix) {
result.push(format!("{}{}\n", indent, matrix_fragment));
}
}
Expand Down Expand Up @@ -391,9 +391,32 @@ impl MatrixFragment {
}
}

pub fn get(&self, key: (usize, usize)) -> &Complex<f64> {
self.matrix.get(key).unwrap()
}

/// Multiplies current matrix by ID, expanding it to fit more qubits.
pub fn expand(&self) -> MatrixFragment { self * &MatrixFragment::id() }

pub fn transpose_conjugate(&self) -> MatrixFragment {
// TODO: Cache and re-use conjugate transposes.
if self.affected_qubits == 1 {
Self::new(array![
[self.get((0, 0)).conj(), self.get((1, 0)).conj()],
[self.get((0, 1)).conj(), self.get((1, 1)).conj()],
], self.affected_qubits)
} else if self.affected_qubits == 2 {
Self::new(array![
[self.get((0, 0)).conj(), self.get((1, 0)).conj(), self.get((2, 0)).conj(), self.get((3, 0)).conj()],
[self.get((0, 1)).conj(), self.get((1, 1)).conj(), self.get((2, 1)).conj(), self.get((3, 1)).conj()],
[self.get((0, 2)).conj(), self.get((1, 2)).conj(), self.get((2, 2)).conj(), self.get((3, 2)).conj()],
[self.get((0, 3)).conj(), self.get((1, 3)).conj(), self.get((2, 3)).conj(), self.get((3, 3)).conj()],
], self.affected_qubits)
} else {
panic!("Can't transpose a matrix covering {} qubits.", self.affected_qubits);
}
}

#[rustfmt::skip]
pub fn X(radians: &f64) -> MatrixFragment {
MatrixFragment {
Expand Down Expand Up @@ -426,6 +449,19 @@ impl MatrixFragment {
affected_qubits: 1
}
}

#[rustfmt::skip]
pub fn Had() -> MatrixFragment {
let one_sq2 = C!(1. / 2.0f64.sqrt(), 0.);
MatrixFragment {
matrix: array![
[one_sq2 * C!(1.0, 0.), one_sq2 * C!(1.0, 0.)],
[one_sq2 * C!(1.0, 0.), one_sq2 * C!(-1.0, 0.)]
],
affected_qubits: 1
}
}

#[rustfmt::skip]
pub fn CX(radians: &f64) -> MatrixFragment {
MatrixFragment {
Expand Down Expand Up @@ -536,25 +572,71 @@ impl Display for MatrixFragment {
fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result { stringified_matrix(&self.matrix, f) }
}

/// Multiply both matrix fragments together. Dimensions are expected to be equal.
fn multiply(left: &MatrixFragment, right: &MatrixFragment) -> MatrixFragment {
let mult = |left_index, right_index| -> Complex64 {
left.get(left_index) * right.get(right_index)
};

// TODO: Find specialized MM library or algorithm, current one not doing it. This is just
// a brute-force implementation for testing.
if left.affected_qubits == 1 {
MatrixFragment::new(array![
[
mult((0, 0), (0, 0)) + mult((0, 1), (1, 0)),
mult((0, 0), (1, 0)) + mult((0, 1), (1, 1))
],
[
mult((1, 0), (0, 0)) + mult((1, 1), (1, 0)),
mult((1, 0), (1, 0)) + mult((1, 1), (1, 1))
]
], left.affected_qubits)
} else if left.affected_qubits == 2 {
MatrixFragment::new(array![
[
mult((0, 0), (0, 0)) + mult((0, 1), (1, 0)) + mult((0, 2), (2, 0)) + mult((0, 3), (3, 0)),
mult((0, 0), (0, 1)) + mult((0, 1), (1, 1)) + mult((0, 2), (2, 1)) + mult((0, 3), (3, 1)),
mult((0, 0), (0, 2)) + mult((0, 1), (1, 2)) + mult((0, 2), (2, 2)) + mult((0, 3), (3, 2)),
mult((0, 0), (0, 3)) + mult((0, 1), (1, 3)) + mult((0, 2), (2, 3)) + mult((0, 3), (3, 3))
],
[
mult((1, 0), (0, 0)) + mult((1, 1), (1, 0)) + mult((1, 2), (2, 0)) + mult((1, 3), (3, 0)),
mult((1, 0), (0, 1)) + mult((1, 1), (1, 1)) + mult((1, 2), (2, 1)) + mult((1, 3), (3, 1)),
mult((1, 0), (0, 2)) + mult((1, 1), (1, 2)) + mult((1, 2), (2, 2)) + mult((1, 3), (3, 2)),
mult((1, 0), (0, 3)) + mult((1, 1), (1, 3)) + mult((1, 2), (2, 3)) + mult((1, 3), (3, 3))
],
[
mult((2, 0), (0, 0)) + mult((2, 1), (1, 0)) + mult((2, 2), (2, 0)) + mult((2, 3), (3, 0)),
mult((2, 0), (0, 1)) + mult((2, 1), (1, 1)) + mult((2, 2), (2, 1)) + mult((2, 3), (3, 1)),
mult((2, 0), (0, 2)) + mult((2, 1), (1, 2)) + mult((2, 2), (2, 2)) + mult((2, 3), (3, 2)),
mult((2, 0), (0, 3)) + mult((2, 1), (1, 3)) + mult((2, 2), (2, 3)) + mult((2, 3), (3, 3))
],
[
mult((3, 0), (0, 0)) + mult((3, 1), (1, 0)) + mult((3, 2), (2, 0)) + mult((3, 3), (3, 0)),
mult((3, 0), (0, 1)) + mult((3, 1), (1, 1)) + mult((3, 2), (2, 1)) + mult((3, 3), (3, 1)),
mult((3, 0), (0, 2)) + mult((3, 1), (1, 2)) + mult((3, 2), (2, 2)) + mult((3, 3), (3, 2)),
mult((3, 0), (0, 3)) + mult((3, 1), (1, 3)) + mult((3, 2), (2, 3)) + mult((3, 3), (3, 3))
],

], left.affected_qubits)
} else {
panic!("Attempted multiplication on unknown size of matrix.")
}
}

impl Mul for MatrixFragment {
type Output = MatrixFragment;

fn mul(self, rhs: Self) -> Self::Output {
MatrixFragment::new(
self.matrix * rhs.matrix,
self.affected_qubits + rhs.affected_qubits
)
multiply(&self, &rhs)
}
}

impl Mul for &MatrixFragment {
type Output = MatrixFragment;

fn mul(self, rhs: Self) -> Self::Output {
MatrixFragment::new(
&self.matrix * &rhs.matrix,
self.affected_qubits + rhs.affected_qubits
)
multiply(&self, &rhs)
}
}

Expand All @@ -578,21 +660,22 @@ type EntangledFragment = StateFragment;
/// smaller ones.
#[derive(Clone)]
pub struct StateFragment {
matrix: Array2<Complex<f64>>,

/// This can also be interpreted as qubits-affected when the fragment is considered a gate.
represented_qubits: i32
matrix_fragment: MatrixFragment,
}

impl StateFragment {
pub fn new(matrix_fragment: MatrixFragment) -> StateFragment {
StateFragment { matrix_fragment }
}

#[rustfmt::skip]
pub fn DefaultQubit() -> QubitFragment {
StateFragment {
matrix: array![
[C!(1.0, 0.0), C!(0.0, 0.0)],
[C!(0.0, 0.0), C!(0.0, 0.0)]
],
represented_qubits: 1
matrix_fragment:
MatrixFragment::new(array![
[C!(1.0, 0.0), C!(0.0, 0.0)],
[C!(0.0, 0.0), C!(0.0, 0.0)]
], 1)
}
}

Expand All @@ -601,61 +684,78 @@ impl StateFragment {
#[rustfmt::skip]
pub fn DefaultEntangled() -> EntangledFragment {
StateFragment {
matrix: array![
matrix_fragment:
MatrixFragment::new(array![
[C!(1.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0)],
[C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0)],
[C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0)],
[C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0), C!(0.0, 0.0)]
],
represented_qubits: 2
], 2)
}
}

pub fn get(&self, key: (usize, usize)) -> &Complex<f64> {
self.matrix_fragment.get(key)
}

pub fn represented_qubits(&self) -> i32 {
self.matrix_fragment.affected_qubits
}

pub fn EntangledFromExisting(
top_left: &Ptr<QubitFragment>, bottom_right: &Ptr<QubitFragment>
) -> EntangledFragment {
if top_left.represented_qubits != 1 || bottom_right.represented_qubits != 1 {
if top_left.represented_qubits() != 1 || bottom_right.represented_qubits() != 1 {
panic!("To create an entangled state both arguments must be qubits.")
}

let zero_zero = *top_left.matrix.get((0, 0)).unwrap();
let zero_one = *top_left.matrix.get((0, 1)).unwrap();
let one_zero = *top_left.matrix.get((1, 0)).unwrap();
let one_one = *top_left.matrix.get((1, 1)).unwrap();
let zero_zero = *top_left.get((0, 0));
let zero_one = *top_left.get((0, 1));
let one_zero = *top_left.get((1, 0));
let one_one = *top_left.get((1, 1));

let two_two = *bottom_right.matrix.get((0, 0)).unwrap();
let two_three = *bottom_right.matrix.get((0, 1)).unwrap();
let three_two = *bottom_right.matrix.get((1, 0)).unwrap();
let three_three = *bottom_right.matrix.get((1, 1)).unwrap();
let two_two = *bottom_right.get((0, 0));
let two_three = *bottom_right.get((0, 1));
let three_two = *bottom_right.get((1, 0));
let three_three = *bottom_right.get((1, 1));
StateFragment {
matrix: array![
matrix_fragment:
MatrixFragment::new(array![
[zero_zero, zero_one, C!(0.0, 0.0), C!(0.0, 0.0)],
[one_zero, one_one, C!(0.0, 0.0), C!(0.0, 0.0)],
[C!(0.0, 0.0), C!(0.0, 0.0), two_two, two_three],
[C!(0.0, 0.0), C!(0.0, 0.0), three_two, three_three]
],
represented_qubits: 2
], 2)
}
}

pub fn apply(&mut self, gate: &MatrixFragment) -> Option<String> {
// If a fragment is larger than us we just ignore it, otherwise try and expand the gate to
// fit the state size.
if gate.affected_qubits < self.represented_qubits {
if gate.affected_qubits < self.represented_qubits() {
let gate = &gate.expand();
}

if self.represented_qubits != gate.affected_qubits {
if self.represented_qubits() != gate.affected_qubits {
return Some(String::from("Can't expand fragment to size of the state."));
}

self.matrix.mul_assign(&gate.matrix);
let mut result = gate * &self.matrix_fragment;
self.matrix_fragment = result * gate.transpose_conjugate();
None
}
}

impl Display for StateFragment {
fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result { stringified_matrix(&self.matrix, f) }
fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result { stringified_matrix(&self.matrix_fragment.matrix, f) }
}

impl Mul for StateFragment {
type Output = StateFragment;

fn mul(self, rhs: Self) -> Self::Output {
StateFragment::new(self.matrix_fragment * rhs.matrix_fragment)
}
}

#[derive(Clone)]
Expand Down Expand Up @@ -1257,7 +1357,6 @@ impl QuantumProjection {
for inst in self.instructions.iter() {
log!(Level::Info, "{}", inst.to_string());
}
log!(Level::Info, "Projection results:");

// Order results so you can easily compare two side-by-side.
let mut result_values = self
Expand All @@ -1268,6 +1367,12 @@ impl QuantumProjection {
.iter()
.collect::<Vec<_>>();
result_values.sort_by(|(left_key, _), (right_key, _)| left_key.cmp(right_key));
if result_values.is_empty() {
log!(Level::Info, "No results from this execution.");
} else {
log!(Level::Info, "Projection results:");
}

for (key, value) in result_values.iter() {
log!(Level::Info, " \"{}\": {}", key, value);
}
Expand Down Expand Up @@ -1385,3 +1490,59 @@ impl Display for AnalysisResult {
f.debug_map().entries(self.distribution.iter()).finish()
}
}

#[cfg(test)]
mod tests {
use std::borrow::Borrow;
use std::fmt::{Display, Formatter};
use crate::analysis::{GateFragment, QubitFragment};

#[test]
fn X() {
let mut qubit = QubitFragment::DefaultQubit();
let result = qubit.apply(&GateFragment::X(&90.0));
assert!(result.is_none());

assert_eq!(qubit.get((0, 0)).re, 0.);
assert_eq!(qubit.get((0, 1)).re, 0.);
assert_eq!(qubit.get((1, 0)).re, 0.);
assert_eq!(qubit.get((1, 1)).re, 1.);
}

#[test]
fn Z() {
let mut qubit = QubitFragment::DefaultQubit();
let result = qubit.apply(&GateFragment::Z(&90.0));
assert!(result.is_none());

assert_eq!(qubit.get((0, 0)).re, 1.);
assert_eq!(qubit.get((0, 1)).re, 0.);
assert_eq!(qubit.get((1, 0)).re, 0.);
assert_eq!(qubit.get((1, 1)).re, 0.);
}

#[test]
fn Y() {
let mut qubit = QubitFragment::DefaultQubit();
let result = qubit.apply(&GateFragment::Y(&90.0));
assert!(result.is_none());

assert_eq!(qubit.get((0, 0)).re, 0.);
assert_eq!(qubit.get((0, 1)).re, 0.);
assert_eq!(qubit.get((1, 0)).re, 0.);
assert_eq!(qubit.get((1, 1)).re, -1.);
}

#[test]
fn Had() {
let mut qubit = QubitFragment::DefaultQubit();
let result = qubit.apply(&GateFragment::Had());
qubit.apply(&GateFragment::Had());
assert!(result.is_none());

assert_eq!(qubit.get((0, 0)).re, 0.);
assert_eq!(qubit.get((0, 1)).re, 0.);
assert_eq!(qubit.get((1, 0)).re, 0.);
assert_eq!(qubit.get((1, 1)).re, 1.);
}
}
6 changes: 6 additions & 0 deletions src/rasqal/src/execution.rs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -289,4 +289,10 @@ mod tests {

#[test]
fn execute_bell_theta_minus() { run(&"../tests/files/qir/bell_theta_minus.ll"); }

#[test]
fn execute_basic_cudaq() {
let config = RasqalConfig::default().with_trace_projections();
run_with_config(&"../tests/files/qir/basic_cudaq.ll", config);
}
}

0 comments on commit 0738595

Please sign in to comment.