diff --git a/tests/test_chemistry.py b/tests/test_chemistry.py
index 86b2dc33..c4369442 100644
--- a/tests/test_chemistry.py
+++ b/tests/test_chemistry.py
@@ -8,6 +8,8 @@
 from tequila.objective import ExpectationValue
 from tequila.quantumchemistry.encodings import known_encodings
 from tequila.simulators.simulator_api import simulate
+import tequila.quantumchemistry.qc_base as qcb
+import tequila.tools.random_generators as rg
 
 HAS_PYSCF = "pyscf" in qc.INSTALLED_QCHEMISTRY_BACKENDS
 HAS_PSI4 = "psi4" in qc.INSTALLED_QCHEMISTRY_BACKENDS
@@ -725,3 +727,55 @@ def test_orbital_optimization_hcb(geometry):
     assert numpy.isclose(opt1.energy,opt2.energy,atol=1.e-5)
     assert time1 < time2
     assert (numpy.isclose(opt1.mo_coeff,opt2.mo_coeff, atol=1.e-5)).all()
+
+@pytest.mark.parametrize("transformation", ["JordanWigner", "ReorderedJordanWigner", "BravyiKitaev", "BravyiKitaevTree"])
+@pytest.mark.parametrize("size", [2, 8])
+def test_givens_on_molecule(size, transformation):
+    # dummy one-electron integrals
+    h = numpy.ones(shape=[size,size])
+    # dummy two-electron integrals
+    g = numpy.ones(shape=[size, size, size, size])
+
+    U = rg.generate_random_unitary(size)
+
+    # transformed integrals
+    th = (U.T.dot(h)).dot(U)
+    tg = numpy.einsum("ijkx, xl -> ijkl", g, U, optimize='greedy')
+    tg = numpy.einsum("ijxl, xk -> ijkl", tg, U, optimize='greedy')
+    tg = numpy.einsum("ixkl, xj -> ijkl", tg, U, optimize='greedy')
+    tg = numpy.einsum("xjkl, xi -> ijkl", tg, U, optimize='greedy')
+
+    # original molecule/H
+    mol = tq.Molecule(geometry="He 0.0 0.0 0.0", nuclear_repulsion=0.0, one_body_integrals=h, two_body_integrals=g, basis_set="dummy", transformation=transformation)
+    H = mol.make_hamiltonian()
+    # transformed molecule/H
+    tmol = tq.Molecule(geometry="He 0.0 0.0 0.0", nuclear_repulsion=0.0, one_body_integrals=th, two_body_integrals=tg,basis_set="dummy", transformation=transformation)
+    tH = tmol.make_hamiltonian()
+
+    # transformation in qubit space (this corresponds to the U above)
+    UR = mol.get_givens_circuit(U) # Works!
+
+    # test circuit
+    circuit = rg.make_random_circuit(size)
+
+    # create expectation values and see if they are the same
+    E1 = tq.ExpectationValue(U=circuit, H=tH)
+    E2 = tq.ExpectationValue(U=circuit + UR, H=H)
+
+    result1 = tq.simulate(E1)
+    result2 = tq.simulate(E2)
+    
+    assert numpy.isclose(result1, result2)
+
+@pytest.mark.parametrize("size", [2, 8])
+def test_givens_decomposition(size):
+    # generate random unitary
+    unitary = rg.generate_random_unitary(size)
+
+    # decompose givens
+    theta_list, phi_list = qcb.get_givens_decomposition(unitary)
+
+    # reconstruct original unitary from givens
+    reconstructed_matrix = qcb.reconstruct_matrix_from_givens(unitary.shape[0], theta_list, phi_list)
+    
+    assert numpy.allclose(unitary, reconstructed_matrix)
diff --git a/tests/test_givens.py b/tests/test_givens.py
deleted file mode 100644
index 399c4d6d..00000000
--- a/tests/test_givens.py
+++ /dev/null
@@ -1,60 +0,0 @@
-import numpy
-import tequila as tq
-import tequila.quantumchemistry.qc_base as qc
-import tequila.tools.random_generators as rg
-import random
-
-transformations = ["JordanWigner", "ReorderedJordanWigner", "BravyiKitaev", "BravyiKitaevTree"]
-def test_givens_on_molecule():
-    # random size and transformation
-    size = random.randint(2, 10)
-    transformation = random.choice(transformations)
-
-    # dummy one-electron integrals
-    h = numpy.ones(shape=[size,size])
-    # dummy two-electron integrals
-    g = numpy.ones(shape=[size, size, size, size])
-
-    U = rg.generate_random_unitary(size)
-
-    # transformed integrals
-    th = (U.T.dot(h)).dot(U)
-    tg = numpy.einsum("ijkx, xl -> ijkl", g, U, optimize='greedy')
-    tg = numpy.einsum("ijxl, xk -> ijkl", tg, U, optimize='greedy')
-    tg = numpy.einsum("ixkl, xj -> ijkl", tg, U, optimize='greedy')
-    tg = numpy.einsum("xjkl, xi -> ijkl", tg, U, optimize='greedy')
-
-    # original molecule/H
-    mol = tq.Molecule(geometry="He 0.0 0.0 0.0", nuclear_repulsion=0.0, one_body_integrals=h, two_body_integrals=g, basis_set="dummy", transformation=transformation)
-    H = mol.make_hamiltonian()
-    # transformed molecule/H
-    tmol = tq.Molecule(geometry="He 0.0 0.0 0.0", nuclear_repulsion=0.0, one_body_integrals=th, two_body_integrals=tg,basis_set="dummy", transformation=transformation)
-    tH = tmol.make_hamiltonian()
-
-    # transformation in qubit space (this corresponds to the U above)
-    UR = mol.get_givens_circuit(U) # Works!
-
-    # test circuit
-    circuit = rg.make_random_circuit(size)
-
-    # create expectation values and see if they are the same
-    E1 = tq.ExpectationValue(U=circuit, H=tH)
-    E2 = tq.ExpectationValue(U=circuit + UR, H=H)
-
-    result1 = tq.simulate(E1)
-    result2 = tq.simulate(E2)
-    
-    assert numpy.isclose(result1, result2)
-
-def test_givens_decomposition():
-    # random unitary of random size
-    size = random.randint(2, 10)
-    unitary = rg.generate_random_unitary(size)
-
-    # decompose givens
-    theta_list, phi_list = qc.get_givens_decomposition(unitary)
-
-    # reconstruct original unitary from givens
-    reconstructed_matrix = qc.reconstruct_matrix_from_givens(unitary.shape[0], theta_list, phi_list)
-    
-    assert numpy.allclose(unitary, reconstructed_matrix)