Update library to time dependent hamiltonians and damping operators

Signed-off-by: Budai <[email protected]>

Fermi_Hub_multi_qubit.ipynb

@@ -65,7 +65,7 @@
     "\n",
     "\n",
     "# Fermi-Hubbard Hamiltonian\n",
-    "def H_operator():\n",
+    "def H_operator(t = 0):\n",
     "    return Qobj(\n",
     "        qib.operator.FermiHubbardHamiltonian(\n",
     "            field_hamil, FermiHubbard_t, FermiHubbard_u, False\n",

qsimulations.py

@@ -1,5 +1,4 @@
 # Utility function, system and constants
-import matplotlib.pyplot as plt
 import scipy.linalg as la
 import numpy as np
 from qutip import *
@@ -154,7 +153,7 @@ def set_nr_of_ancillas(self, nr):
         self._nrAncilla = nr
         self._update_module_varibles()
 
-    def H_op():
+    def H_op(t = 0):
         """Time dependent periodic Hamiltonian for TFIM model
 
         Args:
@@ -165,7 +164,7 @@ def H_op():
         """
         return 0
 
-    def V_op(i):
+    def V_op(i, t = 0):
         """Damping operators, should be overwritten when system is designed
 
         Args:
@@ -174,10 +173,10 @@ def V_op(i):
         """
         return 0
 
-    def _prep_energy_states(self):
+    def _prep_energy_states(self, t = 0):
         """Prepare pare highest energy level state and ground state"""
-        eigenValues, eigenVectors = la.eig(self.H_op().full())
-        self._systemSize = (int)(np.round(np.sqrt(np.size(self.H_op().full()[:, 0]))))
+        eigenValues, eigenVectors = la.eig(self.H_op(t).full())
+        self._systemSize = (int)(np.round(np.sqrt(np.size(self.H_op(t).full()[:, 0]))))
         self._update_module_varibles()
 
         idx = eigenValues.argsort()
@@ -190,7 +189,7 @@ def _prep_energy_states(self):
         self.rho_ground = Qobj(self.psi_ground.conj().T @ self.psi_ground)
         self.rho_0 = Qobj(self.psi_0.conj().T @ self.psi_0)
 
-    def H_op_derivative(self):
+    def H_op_derivative(t = 0):
         """Time derivative of periodic Hamiltonian
 
         Args:
@@ -199,9 +198,9 @@ def H_op_derivative(self):
         Returns:
             np.array: The drived Hamiltonian n requested time stamp
         """
-        return Qobj(np.zeros((self._systemSizeDim, self._systemSizeDim)))
+        return 0
 
-    def V_op_derivative(self, i):
+    def V_op_derivative(i, t = 0):
         """Time derivative of jump operators
 
         Args:
@@ -211,7 +210,7 @@ def V_op_derivative(self, i):
         Returns:
             Qobj: Time derivative of the selected jump operator
         """
-        return Qobj(np.zeros((self._systemSizeDim, self._systemSizeDim)))
+        return 0
 
     def sum_of_V_dag_V(self):
         """Summation of jump operators
@@ -227,7 +226,7 @@ def sum_of_V_dag_V(self):
             sum = sum + self.V_op(j).full().conj().T @ self.V_op(j).full()
         return Qobj(sum)
 
-    def H_tilde_first_order(self, dt):
+    def H_tilde_first_order(self, dt, t = 0):
         """Form of H tilde:
         H =
         sqrt(dt)H   |  V_1^{\\dag}   |   V_2^{\\dag}  |       0
@@ -241,7 +240,7 @@ def H_tilde_first_order(self, dt):
         sum = Qobj(
             np.kron(
                 outer_prod(0, 0, self._nrAncillaDim).full(),
-                (np.sqrt(dt) * self.H_op()).full(),
+                (np.sqrt(dt) * self.H_op(t)).full(),
             )
         )
         # First Row
@@ -270,16 +269,16 @@ def H_tilde_first_order(self, dt):
             ],
         )
 
-    def H_tilde_second_order(self, dt):
+    def H_tilde_second_order(self, dt, t = 0):
         """Second order approximation fo H tilde
         https://arxiv.org/pdf/2311.15533 page 27 equation B10
 
         Args:
                 t (float): time stamp
         """
 
-        sum_tmp = np.sqrt(dt) * (self.H_op().full()) + np.power(dt, 3 / 2) * (
-            -1 / 12 * anti_commute(self.H_op().full(), self.sum_of_V_dag_V().full())
+        sum_tmp = np.sqrt(dt) * (self.H_op(t).full()) + np.power(dt, 3 / 2) * (
+            -1 / 12 * anti_commute(self.H_op(t).full(), self.sum_of_V_dag_V().full())
         )
         sum = Qobj(np.kron(outer_prod(0, 0, self._nrAncillaDim).full(), sum_tmp))
 
@@ -288,7 +287,7 @@ def H_tilde_second_order(self, dt):
                 anti_commute(self.V_op(j).full(), self.V_op(0).full())
                 + self.V_op_derivative(j).full()
                 + 1 / 6 * self.V_op(j).full() @ self.sum_of_V_dag_V().full()
-                + 1j / 2 * self.V_op(j).full() @ self.H_op().full()
+                + 1j / 2 * self.V_op(j).full() @ self.H_op(t).full()
             )
             sum = (
                 sum

test_qsimulations.py

@@ -1,14 +1,7 @@
-from scipy.linalg import sqrtm, cosm, sinm, expm
-import matplotlib.pyplot as plt
-from scipy.special import erf
-import scipy.linalg as la
 import qsimulations
 import numpy as np
 from qutip import *
 import unittest
-import math
-import qib
-
 
 class Test(unittest.TestCase):
     """
@@ -18,7 +11,7 @@ class Test(unittest.TestCase):
     testObj = qsimulations.qsimulations(0, 0, 0)
 
     def test_0_set_Hamiltonian(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(np.array([[0.0, 0.0], [0.0, 1.0]]))
 
         self.testObj.H_op = H_test
@@ -39,7 +32,7 @@ def test_3_set_nr_of_ancillas(self):
         self.assertEqual(self.testObj._nrAncilla, 1)
 
     def test_4_qobj_hamiltonian_test(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(np.array([[0.0, 0.0], [0.0, 1.0]]))
 
         self.testObj.H_op = H_test
@@ -48,7 +41,7 @@ def H_test():
         )
 
     def test_5_qobj_hamiltonian_type_test(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(np.array([[0.0, 0.0], [0.0, 1.0]]))
 
         self.testObj.H_op = H_test
@@ -100,7 +93,7 @@ def test_7_Pauli_arrays(self):
         np.testing.assert_array_equal(test_IZ, correctValue)
 
     def test_8_test_damping_operator(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [
@@ -156,7 +149,7 @@ def test_11_anti_commute_test(self):
         np.testing.assert_array_equal(result, qsimulations.anti_commute(tmp1, tmp2))
 
     def test_12_prepare_energy_states(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [
@@ -181,7 +174,7 @@ def H_test():
         self.assertIsInstance(self.testObj.rho_0, Qobj)
 
     def test_13_sum_of_V(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [
@@ -222,7 +215,7 @@ def _test_damping_Fermi_Hubbard(i):
         np.testing.assert_array_equal(self.testObj.sum_of_V_dag_V().full(), testValue)
 
     def test_14_T_first_ord_type_check(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [
@@ -259,7 +252,7 @@ def _test_damping_Fermi_Hubbard(i):
         self.assertIsInstance(self.testObj.H_tilde_first_order(0.1), Qobj)
 
     def test_15_T_second_ord_type_check(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [
@@ -297,7 +290,7 @@ def _test_damping_Fermi_Hubbard(i):
         self.assertIsInstance(self.testObj.H_tilde_second_order(0.1), Qobj)
 
     def test_16_check_H_psi_ground_state(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [
@@ -316,7 +309,7 @@ def H_test():
         )
 
     def test_17_check_H_psi_0_state(self):
-        def H_test():
+        def H_test(t=0):
             return Qobj(
                 np.array(
                     [