GKP qudit

bosonic/codes/gkp.py

@@ -13,7 +13,7 @@
 
 class GKPQubit(BosonicQubit):
     """
-    GKP Qubit Class.
+    GKP Qudit Class.
     """
     name = "gkp"
 
@@ -22,7 +22,11 @@ def _params_validation(self):
 
         if "delta" not in self.params:
             self.params["delta"] = 0.25
-        self.params["l"] = 2.0 * jnp.sqrt(jnp.pi)
+
+        if 'd' not in self.params:
+            self.params['d'] = 2
+
+        self.params["l"] = jnp.sqrt(jnp.pi * self.params['d'])
         s_delta = jnp.sinh(self.params["delta"] ** 2)
         self.params["epsilon"] = s_delta * self.params["l"]
 
@@ -33,71 +37,39 @@ def _gen_common_gates(self) -> None:
         super()._gen_common_gates()
 
         # phase space
-        self.common_gates["x"] = (
-            self.common_gates["a_dag"] + self.common_gates["a"]
-        ) / jnp.sqrt(2.0)
-        self.common_gates["p"] = (
-            1.0j * (self.common_gates["a_dag"] - self.common_gates["a"]) / jnp.sqrt(2.0)
-        )
+        self.common_gates["x"] = (self.common_gates["a_dag"] + self.common_gates["a"]) / jnp.sqrt(self.params['d'])
+        self.common_gates["p"] = (1.0j * (self.common_gates["a_dag"] - self.common_gates["a"]) / jnp.sqrt(self.params['d']))
 
         # finite energy
-        self.common_gates["E"] = jqt.expm(
-            -self.params["delta"] ** 2
-            * self.common_gates["a_dag"]
-            @ self.common_gates["a"]
-        )
-        self.common_gates["E_inv"] = jqt.expm(
-            self.params["delta"] ** 2
-            * self.common_gates["a_dag"]
-            @ self.common_gates["a"]
-        )
+        self.common_gates["E"] = jqt.expm(-self.params["delta"] ** 2* self.common_gates["a_dag"]@ self.common_gates["a"])
+        self.common_gates["E_inv"] = jqt.expm(self.params["delta"] ** 2* self.common_gates["a_dag"]@ self.common_gates["a"])
 
         # axis
         x_axis, z_axis = self._get_axis()
         y_axis = x_axis + z_axis
 
         # gates
-        X_0 = jqt.expm(1.0j * self.params["l"] / 2.0 * z_axis)
-        Z_0 = jqt.expm(1.0j * self.params["l"] / 2.0 * x_axis)
+        X_0 = jqt.expm(1.0j * self.params["l"] * jnp.sqrt(2) / self.params['d'] * z_axis)
+        Z_0 = jqt.expm(1.0j * self.params["l"] * jnp.sqrt(2) / self.params['d']* x_axis)
         Y_0 = 1.0j * X_0 @ Z_0
+
         self.common_gates["X"] = self._make_op_finite_energy(X_0)
         self.common_gates["Z"] = self._make_op_finite_energy(Z_0)
         self.common_gates["Y"] = self._make_op_finite_energy(Y_0)
 
         # symmetric stabilizers and gates
-        self.common_gates["Z_s_0"] = self._symmetrized_expm(
-            1.0j * self.params["l"] / 2.0 * x_axis
-        )
-        self.common_gates["S_x_0"] = self._symmetrized_expm(
-            1.0j * self.params["l"] * z_axis
-        )
-        self.common_gates["S_z_0"] = self._symmetrized_expm(
-            1.0j * self.params["l"] * x_axis
-        )
-        self.common_gates["S_y_0"] = self._symmetrized_expm(
-            1.0j * self.params["l"] * y_axis
-        )
+        self.common_gates["Z_s_0"] = self._symmetrized_expm(1.0j * self.params["l"] / self.params['d'] * x_axis)
+        self.common_gates["S_x_0"] = self._symmetrized_expm(1.0j * self.params["l"] * z_axis * jnp.sqrt(2))
+        self.common_gates["S_z_0"] = self._symmetrized_expm(1.0j * self.params["l"] * x_axis * jnp.sqrt(2))
+        self.common_gates["S_y_0"] = self._symmetrized_expm(1.0j * self.params["l"] * y_axis * jnp.sqrt(2))
 
     def _get_basis_z(self) -> Tuple[jqt.Qarray, jqt.Qarray]:
         """
-        Construct basis states |+-x>, |+-y>, |+-z>.
-        step 1: use ideal GKP stabilizers to find ideal GKP |+z> state
-        step 2: make ideal eigenvector finite energy
-            We want the groundstate of H = E H_0 E⁻¹.
-            So, we can begin by find the groundstate of H_0 -> |λ₀⟩
-            Then, we know that E|λ₀⟩ = |λ⟩ is the groundstate of H.
-            pf. H|λ⟩ = (E H_0 E⁻¹)(E|λ₀⟩) = E H_0 |λ₀⟩ = λ₀ (E|λ₀⟩) = λ₀|λ⟩
-
-        TODO (if necessary):
-            Alternatively, we could construct a hamiltonian using
-            finite energy stabilizers S_x, S_y, S_z, Z_s. However,
-            this would make H = - S_x - S_y - S_z - Z_s non-hermitian.
-            Currently, JAX does not support derivatives of jnp.linalg.eig,
-            while it does support derivatives of jnp.linalg.eigh.
-            Discussion: https://github.com/google/jax/issues/2748
-        """
+        Construct basis states |+z> and |-z> for the GKP qubit.
 
-        # step 1: use ideal GKP stabilizers to find ideal GKP |+z> state
+        Returns:
+            Tuple[jqt.Qarray, jqt.Qarray]: The |+z> and |-z> basis states.
+        """
         H_0 = (
             -self.common_gates["S_x_0"]
             - self.common_gates["S_y_0"]
@@ -111,7 +83,6 @@ def _get_basis_z(self) -> Tuple[jqt.Qarray, jqt.Qarray]:
         # step 2: make ideal eigenvector finite energy
         gstate = self.common_gates["E"] @ gstate_ideal
 
-        N = self.params["N"]
         plus_z = jqt.unit(gstate)
         minus_z = self.common_gates["X"] @ plus_z
         return plus_z, minus_z

bosonic/codes/qubit.py

@@ -20,7 +20,6 @@ class Qubit(BosonicQubit):
     """
     FockQubit
     """
-
     def _params_validation(self):
         super()._params_validation()
         self.params["N"] = 2
@@ -34,15 +33,12 @@ def _get_basis_z(self) -> Tuple[jqt.Qarray, jqt.Qarray]:
         minus_z = jqt.basis(N, 1)
         return plus_z, minus_z
 
-    @property
     def x_U(self) -> jqt.Qarray:
         return jqt.sigmax()
 
-    @property
     def y_U(self) -> jqt.Qarray:
         return jqt.sigmay()
 
-    @property
     def z_U(self) -> jqt.Qarray:
         return jqt.sigmaz()
 

docs/gkpqubit.md

@@ -0,0 +1,64 @@
+# GKP Qubit Implementation
+
+This document provides an overview of the GKP (Gottesman-Kitaev-Preskill) qubit implementation in the `bosonic` library. The GKP qubit is a type of bosonic code used for quantum error correction, and this implementation supports rectangular, square, and hexagonal GKP qubits.
+
+---
+
+## Overview
+The GKP qubit encodes a logical qubit into the continuous-variable space of a harmonic oscillator. This implementation provides:
+- A base `GKPQubit` class for general GKP qubits.
+- Derived classes for specific GKP geometries: `RectangularGKPQubit`, `SquareGKPQubit`, and `HexagonalGKPQubit`.
+- Common gates (`X`, `Y`, `Z`) and methods for generating basis states.
+
+---
+
+## Installation
+To use the GKP qubit implementation, ensure you have the `bosonic` library installed. You can install it using pip:
+```bash
+pip install bosonic
+```
+
+## Usage
+### Creating a GKP Qubit
+To create a GKP qubit, initialize the GKPQubit class with the desired parameters:
+
+```python
+from bosonic.codes.gkp import GKPQubit
+
+# Create a GKP qubit with default parameters
+qubit = GKPQubit(params={"delta": 0.25, "d": 2, "N": 10})
+```
+
+#### Parameters:
+- **delta**: The finite-energy parameter (default: 0.25).
+- **d**: The dimension of the qudit (default: 2 for a qubit).
+- **N**: The truncation level for the harmonic oscillator Hilbert space.
+
+### Generating Basis States
+The `_get_basis_z` method generates the `|+z>` and `|-z>` basis states:
+
+```python
+plus_z, minus_z = qubit._get_basis_z()
+print("|+z> state:", plus_z)
+print("|-z> state:", minus_z)
+```
+
+### Applying Gates
+You can apply common gates like X, Y, and Z to the GKP qubit:
+
+```python
+# Apply the X gate to the |+z> state
+X_gate = qubit.x_U
+result = X_gate @ plus_z
+print("X gate applied to |+z>:", result)
+```
+
+### Custom GKP Geometries
+The library supports rectangular, square, and hexagonal GKP qubits. For example, to create a rectangular GKP qubit:
+
+```python
+from bosonic.codes.gkp import RectangularGKPQubit
+
+rectangular_qubit = RectangularGKPQubit(params={"a": 0.8, "delta": 0.25, "d": 2, "N": 10})
+```
+

tests/gkp.py

@@ -0,0 +1,59 @@
+import unittest
+import jax.numpy as jnp
+import jaxquantum as jqt
+from bosonic.codes.gkp import GKPQubit
+class TestGKPQubit(unittest.TestCase):
+
+    def setUp(self):
+        # Initialize GKPQubit with default parameters
+        self.qubit = GKPQubit(params={"N": 10})
+
+    def test_params_validation(self):
+        # Test if default parameters are set correctly
+        self.assertEqual(self.qubit.params["delta"], 0.25)
+        self.assertEqual(self.qubit.params["d"], 2)
+        self.assertAlmostEqual(self.qubit.params["l"], jnp.sqrt(jnp.pi * 2))
+        self.assertAlmostEqual(self.qubit.params["epsilon"], jnp.sinh(0.25**2) * jnp.sqrt(jnp.pi * 2))
+
+    def test_common_gates(self):
+        # Test if common gates are generated correctly
+        self.assertIn("x", self.qubit.common_gates)
+        self.assertIn("p", self.qubit.common_gates)
+        self.assertIn("E", self.qubit.common_gates)
+        self.assertIn("E_inv", self.qubit.common_gates)
+        self.assertIn("X", self.qubit.common_gates)
+        self.assertIn("Z", self.qubit.common_gates)
+        self.assertIn("Y", self.qubit.common_gates)
+        self.assertIn("Z_s_0", self.qubit.common_gates)
+        self.assertIn("S_x_0", self.qubit.common_gates)
+        self.assertIn("S_z_0", self.qubit.common_gates)
+        self.assertIn("S_y_0", self.qubit.common_gates)
+
+    def test_basis_z(self):
+        plus_z, minus_z = self.qubit._get_basis_z()
+        self.assertIsInstance(plus_z, jqt.Qarray)
+        self.assertIsInstance(minus_z, jqt.Qarray)
+        # Check the shape of the underlying data
+        self.assertEqual(plus_z.data.shape, (self.qubit.params["N"], 1))
+        self.assertEqual(minus_z.data.shape, (self.qubit.params["N"], 1))
+
+    def test_gates(self):
+        # Test if gate properties return correct types
+        self.assertIsInstance(self.qubit.x_U, jqt.Qarray)
+        self.assertIsInstance(self.qubit.y_U, jqt.Qarray)
+        self.assertIsInstance(self.qubit.z_U, jqt.Qarray)
+
+    def test_make_op_finite_energy(self):
+        # Test if finite energy operation is applied correctly
+        op = self.qubit.common_gates["X"]
+        finite_op = self.qubit._make_op_finite_energy(op)
+        self.assertIsInstance(finite_op, jqt.Qarray)
+
+    def test_symmetrized_expm(self):
+        # Test if symmetrized exponential is applied correctly
+        op = self.qubit.common_gates["x"]
+        symm_op = self.qubit._symmetrized_expm(op)
+        self.assertIsInstance(symm_op, jqt.Qarray)
+
+if __name__ == "__main__":
+    unittest.main()