Repository: Agony5757/QPanda-lite (language: python)

Issue:
[feat] Noisy Simulator random seed
The output of the noisy simulator under the same circuit will return the same sampling output.

It should return another similar output. @YunJ1e 


Discussion / hints:

Fix patch:
diff --git a/QPandaLiteCpp/Pybinder/QPandaLitePy.cpp b/QPandaLiteCpp/Pybinder/QPandaLitePy.cpp
index c9fd832..03f718d 100644
--- a/QPandaLiteCpp/Pybinder/QPandaLitePy.cpp
+++ b/QPandaLiteCpp/Pybinder/QPandaLitePy.cpp
@@ -77,6 +77,7 @@ PYBIND11_MODULE(QPandaLitePy, m)
 
 	py::enum_<qpandalite::SupportOperationType>(m, "SupportOperationType")
         .value("HADAMARD", qpandalite::SupportOperationType::HADAMARD)
+        .value("IDENTITY", qpandalite::SupportOperationType::HADAMARD)
         .value("U22", qpandalite::SupportOperationType::U22)
         .value("X", qpandalite::SupportOperationType::X)
         .value("Y", qpandalite::SupportOperationType::Y)
@@ -117,7 +118,7 @@ PYBIND11_MODULE(QPandaLitePy, m)
 		.def_readonly("measurement_error_matrices", &qpandalite::NoisySimulator::measurement_error_matrices)
 		.def("load_opcode", &qpandalite::NoisySimulator::load_opcode)
 		.def("insert_error", &qpandalite::NoisySimulator::insert_error)
-		.def("get_measure_no_readout_error", &qpandalite::NoisySimulator::get_measure_no_readout_error)
+		.def("id", &qpandalite::NoisySimulator::id, py::arg("qn"), py::arg("is_dagger") = false)		.def("get_measure_no_readout_error", &qpandalite::NoisySimulator::get_measure_no_readout_error)
 		.def("get_measure", &qpandalite::NoisySimulator::get_measure)
 
 		.def("hadamard", &qpandalite::NoisySimulator::hadamard, py::arg("qn"), py::arg("dagger") = false)
diff --git a/QPandaLiteCpp/src/noisy_simulator.cpp b/QPandaLiteCpp/src/noisy_simulator.cpp
index ecf57f2..d4bd24f 100644
--- a/QPandaLiteCpp/src/noisy_simulator.cpp
+++ b/QPandaLiteCpp/src/noisy_simulator.cpp
@@ -387,6 +387,11 @@ namespace qpandalite {
 	{
 		_insert_global_error(qubits);
 	}
+	
+	void NoisySimulator::id(size_t qn, bool is_dagger)
+	{
+		id_cont(qn, {}, is_dagger);
+	}
 
 	void NoisySimulator::hadamard(size_t qn, bool is_dagger)
 	{
@@ -425,7 +430,7 @@ namespace qpandalite {
 
 	void NoisySimulator::swap(size_t qn1, size_t qn2, bool is_dagger)
 	{
-		iswap_cont(qn1, qn2, {}, is_dagger);
+		swap_cont(qn1, qn2, {}, is_dagger);
 	}
 
 	void NoisySimulator::iswap(size_t qn1, size_t qn2, bool is_dagger)
@@ -471,8 +476,8 @@ namespace qpandalite {
 	void NoisySimulator::rphi(size_t qn, double phi, double theta, bool is_dagger)
 	{
 		rphi_cont(qn, phi, theta, {}, is_dagger);
-	}
-
+    }
+	
 	void NoisySimulator::toffoli(size_t qn1, size_t qn2, size_t target, bool is_dagger)
 	{
 		toffoli_cont(qn1, qn2, target, {}, is_dagger);
@@ -483,6 +488,19 @@ namespace qpandalite {
 		cswap_cont(qn1, qn2, target, {}, is_dagger);
 	}
 
+	void NoisySimulator::id_cont(size_t qn, const std::vector<size_t>& global_controller, bool is_dagger)
+	{
+		opcodes.emplace_back(
+			OpcodeType(
+			(uint32_t)SupportOperationType::IDENTITY,
+			{ qn },
+			{},
+			is_dagger,
+			global_controller)
+		);
+		insert_error({ qn }, SupportOperationType::IDENTITY);
+	}
+
 	void NoisySimulator::hadamard_cont(size_t qn, const std::vector<size_t>& global_controller, bool is_dagger)
 	{
 		opcodes.emplace_back(
@@ -810,10 +828,7 @@ namespace qpandalite {
 				simulator.cz_cont(opcode.qubits[0], opcode.qubits[1], opcode.global_controller, opcode.dagger);
 				break;
 			case (uint32_t)SupportOperationType::CNOT:
-				simulator.cnot_cont(opcode.qubits[0], opcode.qubits[1], opcode.global_controller, opcode.dagger);
-				break;
-			case (uint32_t)SupportOperationType::SWAP:
-				simulator.swap_cont(opcode.qubits[0], opcode.qubits[1], opcode.global_controller, opcode.dagger);
+				simulator.cnot(opcode.qubits[0], opcode.qubits[1]);
 				break;
 			case (uint32_t)SupportOperationType::ISWAP:
 				simulator.iswap_cont(opcode.qubits[0], opcode.qubits[1], opcode.global_controller, opcode.dagger);
diff --git a/QPandaLiteCpp/src/noisy_simulator.h b/QPandaLiteCpp/src/noisy_simulator.h
index a1b178e..0b98f9e 100644
--- a/QPandaLiteCpp/src/noisy_simulator.h
+++ b/QPandaLiteCpp/src/noisy_simulator.h
@@ -16,6 +16,7 @@ namespace qpandalite {
     {
         __SupportOperationTypeBegin = 1000,
         HADAMARD,
+        IDENTITY,
         U22,
         X,
         Y,
@@ -92,6 +93,7 @@ namespace qpandalite {
         static const std::map<std::string, SupportOperationType> op_type_map =
         {
            {"HADAMARD", SupportOperationType::HADAMARD},
+           {"IDENTITY", SupportOperationType::IDENTITY},
            {"U22", SupportOperationType::U22},
            {"X", SupportOperationType::X},
            {"Y", SupportOperationType::Y},
@@ -194,6 +196,7 @@ namespace qpandalite {
         // Perform bit flip based on measurement 
 
         void hadamard(size_t qn, bool is_dagger = false);
+        void id(size_t qn, bool is_dagger = false);
         void u22(size_t qn, const u22_t& unitary, bool is_dagger = false);
         void x(size_t qn, bool is_dagger = false);
         void z(size_t qn, bool is_dagger = false);
@@ -213,6 +216,7 @@ namespace qpandalite {
         void toffoli(size_t qn1, size_t qn2, size_t target, bool is_dagger = false);
         void cswap(size_t controller, size_t target1, size_t target2, bool is_dagger = false);
 
+        void id_cont(size_t qn, const std::vector<size_t>& global_controller, bool is_dagger = false);
         void hadamard_cont(size_t qn, const std::vector<size_t>& global_controller, bool is_dagger = false);
         void u22_cont(size_t qn, const u22_t& unitary, const std::vector<size_t>& global_controller, bool is_dagger = false);
         void x_cont(size_t qn, const std::vector<size_t>& global_controller, bool is_dagger = false);
diff --git a/QPandaLiteCpp/src/simulator.cpp b/QPandaLiteCpp/src/simulator.cpp
index d879583..e553608 100644
--- a/QPandaLiteCpp/src/simulator.cpp
+++ b/QPandaLiteCpp/src/simulator.cpp
@@ -273,6 +273,7 @@ namespace qpandalite{
             }
         }
     }
+    
     void Simulator::cnot(size_t controller, size_t target, bool is_dagger)
     {
         if (controller >= total_qubit)
@@ -890,8 +891,19 @@ namespace qpandalite{
         {
             if (!control_enable(i, global_controller))
                 continue;
-            if (((i >> controller) & 1) && ((i >> target) & 1))
+            bool control_active = true;
+            
+            // Check if all global controllers are active (in state 1)
+            for (size_t q : global_controller) {
+
+                if (!((i >> q) & 1)) {
+                    control_active = false;
+                    break;
+                }
+            }
+            if (control_active && ((i >> controller) & 1) && ((i >> target) & 1))
             {
+                // Printing i when all conditions are satisfied
                 std::swap(state[i], state[i - pow2(target)]);
             }
         }
@@ -1116,7 +1128,6 @@ namespace qpandalite{
         u22_cont(qn, unitary, global_controller);
     }
 
-
     void Simulator::toffoli_cont(size_t qn1, size_t qn2, size_t target, const std::vector<size_t>& global_controller, bool is_dagger)
     {
         if (global_controller.size() == 0)
@@ -1190,7 +1201,7 @@ namespace qpandalite{
             }
         }
     }
-
+    
     bool Simulator::control_enable(size_t idx, const std::vector<size_t>& global_controller)
     {
         for (auto qid : global_controller)
diff --git a/qpandalite/simulator/QPandaLitePy.pyi b/qpandalite/simulator/QPandaLitePy.pyi
index f8dd88a..ee3e3cf 100644
--- a/qpandalite/simulator/QPandaLitePy.pyi
+++ b/qpandalite/simulator/QPandaLitePy.pyi
@@ -11,6 +11,9 @@ from typing_extensions import override
 def seed(seed_: int) -> None: ...
 def rand() -> float: ...
 
+def seed(seed_: int) -> None: ...
+def rand() -> float: ...
+
 class Simulator:
     def __init__(self) -> None: 
         '''Create a simulator instance (implemented by C++).
@@ -113,6 +116,15 @@ class Simulator:
         '''
         ...
 
+    def swap(self, q1 : int, q2 : int, is_dagger : bool = False) -> None:
+        '''iSWAP gate.
+
+        Args:
+            q1 (int): The first qubit.
+            q2 (int): The second qubit.
+        '''
+        ...
+
     def xy(self, q1 : int, q2 : int, theta: float, is_dagger : bool = False) -> None:
         '''XX+YY gate.
 
@@ -438,6 +450,14 @@ class NoisySimulator:
         '''
         ...
     
+    def id(self, qn : int, is_dagger : bool = False) -> None:
+        '''Identity gate.
+
+        Args:
+            qn (int): The gate qubit.
+        '''
+        ...
+
     def hadamard(self, qn : int, is_dagger : bool = False) -> None:
         '''Hadamard gate.
 
@@ -515,6 +535,15 @@ class NoisySimulator:
         '''
         ...  
     
+    def swap(self, q1 : int, q2 : int, is_dagger : bool = False) -> None:
+        '''CNOT gate.
+
+        Args:
+            controller (int): The controller qubit.
+            target (int): The target qubit.
+        '''
+        ...  
+
     def rx(self, qn : int, angle : float, is_dagger : bool = False) -> None:
         '''Rx gate.
 
@@ -569,7 +598,37 @@ class NoisySimulator:
             theta (float): The rotation angle.
         '''
         ...
-        
+    
+    def toffoli(self, qn1 : int, qn2 : int, target : int, is_dagger : bool = False) -> None:
+        '''Toffoli gate.
+
+        Args:
+            qn1 (int): The controller qubit1.
+            qn2 (int): The controller qubit2.
+            target (int): The target qubit.
+        '''
+        ... 
+    
+    def cswap(self, controller : int, target1 : int, target2 : int, is_dagger : bool = False) -> None:
+        '''CSWAP gate.
+
+        Args:
+            controller (int): The controller qubit.
+            target1 (int): The target qubit1.
+            target2 (int): The target qubit2.
+        '''
+        ... 
+
+    def cnot_cont(self, controller : int, target : int, control : List[int], is_dagger : bool = False) -> None:
+        '''CNOT gate.
+
+        Args:
+            controller (int): The controller qubit.
+            target (int): The target qubit.
+        '''
+        ... 
+
+
     def measure_shots(self, shots: int) -> Dict[int, int]:
         """
         Simulate a quantum measurement multiple times and tally the results.
diff --git a/qpandalite/simulator/originir_simulator.py b/qpandalite/simulator/originir_simulator.py
index ae7d50e..7a814e1 100644
--- a/qpandalite/simulator/originir_simulator.py
+++ b/qpandalite/simulator/originir_simulator.py
@@ -83,8 +83,7 @@ def _simulate_controlled_gate(self, operation, qubit, cbit, parameter, is_dagger
         else:
             raise RuntimeError('Unknown OriginIR operation. '
                                 f'Operation: {operation}. '
-                                f'Full opcode: {(operation, qubit, cbit, 
-                                                    parameter, is_dagger, control_qubits_set)}')
+                                f'Full opcode: {(operation, qubit, cbit, parameter, is_dagger, control_qubits_set)}')
     
     def _simulate_common_gate(self, operation, qubit, cbit, parameter, is_dagger, ):
         if operation == 'RX':
@@ -150,8 +149,7 @@ def _simulate_common_gate(self, operation, qubit, cbit, parameter, is_dagger, ):
         else:
             raise RuntimeError('Unknown OriginIR operation. '
                                 f'Operation: {operation}.'
-                                f'Full opcode: {(operation, qubit, cbit, 
-                                                    parameter, is_dagger)}')
+                                f'Full opcode: {(operation, qubit, cbit, parameter, is_dagger)}')
 
     def simulate_gate(self, operation, qubit, cbit, parameter, is_dagger, control_qubits_set):
         if control_qubits_set:


Test patch:
diff --git a/test/simulator/gate_dependent_noise_101.py b/test/simulator/gate_dependent_noise_101.py
index 3267ac2..4a5dd8b 100644
--- a/test/simulator/gate_dependent_noise_101.py
+++ b/test/simulator/gate_dependent_noise_101.py
@@ -37,17 +37,7 @@
 # print(gate_noise_desc)
 
 
-simulator.hadamard(0)
-simulator.x(0)
-simulator.hadamard(1)
-simulator.z(1)
-simulator.iswap(0, 1)
-simulator.rx(0, np.pi/8)
-simulator.ry(0, np.pi/4)
-simulator.rz(0, np.pi/2)
-simulator.rphi90(1, np.pi/8)
-simulator.rphi180(1, np.pi/4)
-simulator.rphi(1, np.pi/2, np.pi/4)
+# simulator.cnot_cont()
 # Number of measurement shots
 shots = 1024
 
diff --git a/test/simulator/qram_simulator.py b/test/simulator/qram_simulator.py
new file mode 100644
index 0000000..4efc8cc
--- /dev/null
+++ b/test/simulator/qram_simulator.py
@@ -0,0 +1,57 @@
+import numpy as np
+import math
+
+from qpandalite.simulator import Simulator, NoisySimulator
+from qpandalite.qasm_origin import OpenQASM2_Parser
+from qpandalite.circuit_builder import Circuit
+from qpandalite.simulator import seed
+
+import time
+current_time_seed = int(time.time())
+
+# Use the current time as a seed
+seed(current_time_seed)
+
+# Define the noise descriptions
+noise_description = {
+    # "depolarizing": 0.01
+}
+
+# Define the gate noise descriptions
+gate_noise_description = {
+    "X": {"depolarizing": 0.03},
+    "HADAMARD": {"depolarizing": 0.02},
+    "ISWAP": {"depolarizing": 0.02}
+}
+
+# Define the measurement errors
+measurement_error = [(0.01, 0.01), (0.02, 0.02)]
+
+# Create an instance of the NoisySimulator
+simulator = NoisySimulator(2,noise_description, gate_noise_description)
+
+# noise_desc = simulator.noise
+# gate_noise_desc = simulator.gate_dependent_noise
+
+# print(noise_desc)
+# print(gate_noise_desc)
+
+
+simulator.hadamard(0)
+simulator.x(0)
+simulator.hadamard(1)
+simulator.z(1)
+simulator.iswap(0, 1)
+simulator.rx(0, np.pi/8)
+simulator.ry(0, np.pi/4)
+simulator.rz(0, np.pi/2)
+simulator.rphi90(1, np.pi/8)
+simulator.rphi180(1, np.pi/4)
+simulator.rphi(1, np.pi/2, np.pi/4)
+# Number of measurement shots
+shots = 1024
+
+# # Measure the state multiple times
+measurement_results = simulator.measure_shots(shots)
+
+print(measurement_results)
\ No newline at end of file
diff --git a/test/simulator/test_bug.py b/test/simulator/test_bug.py
new file mode 100644
index 0000000..598133a
--- /dev/null
+++ b/test/simulator/test_bug.py
@@ -0,0 +1,11 @@
+import qpandalite.simulator as qsim
+
+sim = qsim.Simulator()
+sim.init_n_qubit(1)
+# sim.x(0)
+sim.rx(0, -1.57)
+sim.rz(0, 1.57)
+sim.hadamard(0)
+
+print(sim.state)
+print(sim.pmeasure(0))
\ No newline at end of file
diff --git a/test/simulator/test_clifford.py b/test/simulator/test_clifford.py
new file mode 100644
index 0000000..6148916
--- /dev/null
+++ b/test/simulator/test_clifford.py
@@ -0,0 +1,70 @@
+from qpandalite.simulator import Simulator, NoisySimulator
+from qpandalite.simulator import seed
+import numpy as np
+import time
+current_time_seed = int(time.time())
+
+# Use the current time as a seed
+seed(current_time_seed)
+
+single_qubit_cliffords = [
+ '',
+ 'H', 'S',
+ 'HS', 'SH', 'SS',
+ 'HSH', 'HSS', 'SHS', 'SSH', 'SSS',
+ 'HSHS', 'HSSH', 'HSSS', 'SHSS', 'SSHS',
+ 'HSHSS', 'HSSHS', 'SHSSH', 'SHSSS', 'SSHSS',
+ 'HSHSSH', 'HSHSSS', 'HSSHSS'
+]
+
+def apply_single_clifford(clifford_string, qsim, target_qbit, inverse=False):
+	if inverse:
+		clifford_string_actual_order = clifford_string[::-1]
+	else:
+		clifford_string_actual_order = clifford_string
+
+	for gate in clifford_string_actual_order:
+		if gate == 'H':
+			qsim.hadamard(target_qbit)
+		else:
+			sign = -1 if inverse else 1
+			qsim.rz(target_qbit, sign * np.pi/2)
+
+
+# Define the noise descriptions
+noise_description = {
+    # "depolarizing": 0.01
+}
+
+# Define the gate noise descriptions
+gate_noise_description = {
+    "X": {"depolarizing": 0.00},
+    "HADAMARD": {"depolarizing": 0.00},
+    "CNOT": {"depolarizing": 0.00},
+    "Y": {"depolarizing": 0.1}
+}
+
+# Define the measurement errors
+measurement_error = [(0.01, 0.01), (0.02, 0.02)]
+culmu = 0
+shots = 2**16
+for C in single_qubit_cliffords:
+	# print(C)
+	simulator = NoisySimulator(1, noise_description, gate_noise_description)
+	# simulator = Simulator()
+
+	# simulator.init_n_qubit(1)
+	apply_single_clifford(C, simulator, 0, inverse=False)
+	simulator.y(0)
+	simulator.y(0)
+	apply_single_clifford(C, simulator, 0, inverse=True)
+	
+	# Measure the state multiple times
+	measurement_results = simulator.measure_shots(shots)
+	# print(simulator.state)
+	# culmu += measurement_results[0]
+	# print(measurement_results[0] / shots)
+
+	# print( / shots)
+
+print(culmu / shots / len(single_qubit_cliffords))
diff --git a/test/simulator/test_plot.py b/test/simulator/test_plot.py
new file mode 100644
index 0000000..6fbdafa
--- /dev/null
+++ b/test/simulator/test_plot.py
@@ -0,0 +1,92 @@
+import json
+import matplotlib.pyplot as plt
+import numpy as np
+
+# # Load the dictionary from the JSON file
+# with open(r'C:\Users\RGZN90601\Desktop\QPanda-lite\data_sup.json', 'r') as file:
+#     data_sup = json.load(file)
+# with open(r'C:\Users\RGZN90601\Desktop\QPanda-lite\data_normal.json', 'r') as file:
+#     data_normal = json.load(file)
+
+# # Assuming 'data' is your dictionary and 'k_values' are the keys as floats
+# k_values_sup = [float(k) for k in data_sup.keys()]
+# predicted_values_sup = [1 - (1 - 2 * k/3) for k in k_values_sup]
+# # predicted_values_sup = [k/2 for k in k_values_sup]
+# k_values_normal = [float(k) for k in data_normal.keys()]
+# predicted_values_normal = [1 - (1 - 2 * k/3) for k in k_values_normal]
+
+# # Invert the data by calculating 1 - value for each value in the datasets
+# inverted_data_sup = {key: [1 - val for val in values] for key, values in data_sup.items()}
+# inverted_data_normal = {key: [1 - val for val in values] for key, values in data_normal.items()}
+
+# plt.figure(figsize=(12, 8))
+
+# # Create the first subplot with shared y-axis
+# ax1 = plt.subplot(1, 2, 1)
+# plt.boxplot(inverted_data_sup.values(), labels=inverted_data_sup.keys())
+# plt.plot(range(1, len(predicted_values_sup) + 1), predicted_values_sup, 'g-o', label='Predicted Suppressed')
+# plt.legend()
+# plt.title('(1 - F) Infidelity Error Suppression')
+# plt.xlabel('Depolarization Strength')
+# plt.ylabel('Value')
+# # plt.yscale('log')
+
+# # Create the second subplot, sharing the y-axis with the first
+# ax2 = plt.subplot(1, 2, 2, sharey=ax1)
+# plt.boxplot(inverted_data_normal.values(), labels=inverted_data_normal.keys())
+# plt.plot(range(1, len(predicted_values_normal) + 1), predicted_values_normal, 'g-o', label='Predicted Normal')
+# plt.title('(1 - F) Infidelity Raw Circuit')
+# plt.xlabel('Depolarization Strength')
+
+# # Hide y-axis labels on the second plot to avoid label duplication
+# plt.setp(ax2.get_yticklabels(), visible=False)
+# plt.legend()
+# plt.tight_layout()
+# plt.show()
+
+# Data
+# Given data
+x = [1e-1, 5e-2, 1e-2, 5e-3, 1e-3, 5e-4, 1e-4]
+p = np.array(x)
+predicted_raw = 2 * p / 3
+predicted_suppression = 1 - (1 - p + 2 * p**2 / 9) / (1 + 2 * p * (2 * p - 3) / 9)
+
+suppression = [0.038031953888349034, 0.017755342411343078, 0.003380587203426303, 
+               0.0016816904917188412, 0.0003341679612717067, 0.00016901651807362436, 
+               3.560617082634665e-05]
+raw = [0.06665436426798499, 0.03324747085571289, 0.006716569264729855, 
+       0.0033316612243652344, 0.0006745656331380578, 0.0003256797790527344, 
+       6.516774495446409e-05]
+
+# Calculate the ratio of suppression to raw
+ratio = np.array(suppression) / np.array(raw)
+predicted_ratio = predicted_suppression / predicted_raw
+
+# Create subplots
+fig, axs = plt.subplots(1, 2, figsize=(15, 6))
+
+# First subplot for Suppression, Raw, and their predicted values
+axs[0].scatter(x, suppression, label='Suppression', marker='o')
+axs[0].scatter(x, raw, label='Raw', marker='o')
+axs[0].plot(x, predicted_suppression, label='Predicted Suppression', marker='x', linestyle='--')
+axs[0].plot(x, predicted_raw, label='Predicted Raw', marker='x', linestyle='--')
+axs[0].set_xscale('log')
+axs[0].set_yscale('log')
+axs[0].set_title('Suppression and Raw with Predictions')
+axs[0].set_xlabel('Depolarization Strength')
+axs[0].set_ylabel('Value')
+axs[0].legend()
+axs[0].grid(True)
+
+# Second subplot for Suppression/Raw ratio
+axs[1].scatter(x, ratio, label='Suppression/Raw Ratio', marker='o', color='red')
+axs[1].plot(x, predicted_ratio, label='Predicted Ratio', marker='x', linestyle='--')
+axs[1].set_xscale('log')
+axs[1].set_title('Suppression/Raw Ratio')
+axs[1].set_xlabel('Depolarization Strength')
+axs[1].set_ylabel('Ratio')
+axs[1].legend()
+axs[1].grid(True)
+
+plt.tight_layout()
+plt.show()
\ No newline at end of file
diff --git a/test/simulator/test_toffoli_2.py b/test/simulator/test_toffoli_2.py
new file mode 100644
index 0000000..54193f9
--- /dev/null
+++ b/test/simulator/test_toffoli_2.py
@@ -0,0 +1,66 @@
+from qpandalite.simulator import Simulator, NoisySimulator
+from qpandalite.simulator import seed
+import numpy as np
+import json
+
+import time
+current_time_seed = int(time.time())
+
+
+noise_strength = 0
+shots = 2**18
+
+# Define the noise descriptions
+noise_description = {
+    # "depolarizing": 0.01
+}
+
+# Define the gate noise descriptions
+gate_noise_description = {
+    "X": {"depolarizing": 0.00},
+    "HADAMARD": {"depolarizing": 0.00},
+    "CNOT": {"depolarizing": 0.00},
+    "IDENTITY": {"depolarizing": noise_strength},
+    "Y": {"depolarizing": 0.00}
+}
+
+# Define the measurement errors
+measurement_error = [(0.01, 0.01), (0.02, 0.02)]
+
+
+# Create an instance of the NoisySimulator
+simulator = NoisySimulator(3, noise_description, gate_noise_description)
+# simulator = Simulator()
+# simulator.init_n_qubit(3)
+simulator.hadamard(0)
+# simulator.x(1)
+# simulator.cnot(0, 1)
+simulator.x(0)
+simulator.cnot(2, 1)
+simulator.cnot_cont(1, 2, [0])
+simulator.cnot(2, 1)
+simulator.x(0)
+
+simulator.x(2)
+
+simulator.x(0)
+simulator.cnot(2, 1)
+simulator.cnot_cont(1, 2, [0])
+# simulator.cnot(2, 1)
+# simulator.x(0)
+
+# simulator.cnot(2, 1)
+# simulator.cnot_cont(1, 2, [0])
+# simulator.cnot(2, 1)
+
+# simulator.cnot(2, 1)
+# simulator.cnot_cont(1, 2, [0])
+# simulator.cnot(2, 1)
+
+# simulator.hadamard(0)
+
+
+measurement_results = simulator.measure_shots(shots)
+
+# single_exps_sup[noise_strength].append(measurement_results[0] / shots)
+print(measurement_results)
diff --git a/test/simulator/test_toffoli_fredkin.py b/test/simulator/test_toffoli_fredkin.py
new file mode 100644
index 0000000..9dd082a
--- /dev/null
+++ b/test/simulator/test_toffoli_fredkin.py
@@ -0,0 +1,148 @@
+from qpandalite.simulator import Simulator, NoisySimulator
+from qpandalite.simulator import seed
+import numpy as np
+import json
+
+import time
+current_time_seed = int(time.time())
+
+# Use the current time as a seed
+seed(current_time_seed)
+
+single_qubit_cliffords = [
+ '',
+ 'H', 'S',
+ 'HS', 'SH', 'SS',
+ 'HSH', 'HSS', 'SHS', 'SSH', 'SSS',
+ 'HSHS', 'HSSH', 'HSSS', 'SHSS', 'SSHS',
+ 'HSHSS', 'HSSHS', 'SHSSH', 'SHSSS', 'SSHSS',
+ 'HSHSSH', 'HSHSSS', 'HSSHSS'
+]
+
+def apply_single_clifford(clifford_string, qsim, target_qbit, inverse=False):
+	if inverse:
+		clifford_string_actual_order = clifford_string[::-1]
+	else:
+		clifford_string_actual_order = clifford_string
+
+	for gate in clifford_string_actual_order:
+		if gate == 'H':
+			qsim.hadamard(target_qbit)
+		else:
+			sign = -1 if inverse else 1
+			qsim.rz(target_qbit, sign * np.pi/2)
+
+single_exps_sup = {}
+single_exps_normal = {}
+
+for noise_strength in [1e-1, 5e-2, 1e-2, 5e-3, 1e-3, 5e-4, 1e-4]:
+	print(noise_strength)
+	shots = 2**18
+
+	# Define the noise descriptions
+	noise_description = {
+	    # "depolarizing": 0.01
+	}
+
+	# Define the gate noise descriptions
+	gate_noise_description = {
+	    "X": {"depolarizing": 0.00},
+	    "HADAMARD": {"depolarizing": 0.00},
+	    "CNOT": {"depolarizing": 0.00},
+	    "IDENTITY": {"depolarizing": noise_strength},
+	    "Y": {"depolarizing": noise_strength},
+	    "Z": {"depolarizing": noise_strength}
+	}
+
+	# Define the measurement errors
+	measurement_error = [(0.01, 0.01), (0.02, 0.02)]
+
+	culmu = 0
+	condition_prob = []
+
+	single_exps_sup[noise_strength] = []
+	single_exps_normal[noise_strength] = []
+
+	for C in single_qubit_cliffords:
+		# Create an instance of the NoisySimulator
+		simulator = NoisySimulator(3, noise_description, gate_noise_description)
+
+		simulator.hadamard(0)
+		# apply_single_clifford(C, simulator, 1, inverse=False)
+
+		simulator.x(0)
+		simulator.cnot(2, 1)
+		simulator.cnot_cont(1, 2, [0])
+		simulator.cnot(2, 1)
+		simulator.x(0)
+
+		simulator.id(2)
+
+		simulator.x(0)
+		simulator.cnot(2, 1)
+		simulator.cnot_cont(1, 2, [0])
+		simulator.cnot(2, 1)
+		simulator.x(0)
+
+		simulator.cnot(2, 1)
+		simulator.cnot_cont(1, 2, [0])
+		simulator.cnot(2, 1)
+
+		simulator.id(2)
+
+		simulator.cnot(2, 1)
+		simulator.cnot_cont(1, 2, [0])
+		simulator.cnot(2, 1)
+
+		simulator.hadamard(0)
+
+		# apply_single_clifford(C, simulator, 1, inverse=True)
+		
+		measurement_results = simulator.measure_shots(shots)
+		aa = 0
+		culmu = 0
+		for key in list(measurement_results.keys()):  # Use list to avoid RuntimeError for changing dict size during iteration
+		    if int(key) % 2 == 0:  # Check if the key is even
+		        culmu += measurement_results[key]
+		    if int(key) in [0,4]:
+		        aa += measurement_results[key]
+		
+		condition_prob.append(1 - aa / culmu)
+	# single_exps_sup[noise_strength].append(measurement_results[0] / shots)
+		# print(measurement_results)
+
+	print('Suppression')
+	print(np.array(condition_prob).mean())
+	# single_exps_sup[noise_strength].append(culmu / shots / len(single_qubit_cliffords))
+	# print((2 * noise_strength **2 - 9 * noise_strength + 9) / (4 * noise_strength **2 - 6 * noise_strength + 9))
+	
+
+	culmu = 0
+
+	for C in single_qubit_cliffords:
+		# print(C)
+		simulator = NoisySimulator(1, noise_description, gate_noise_description)
+
+		apply_single_clifford(C, simulator, 0, inverse=False)
+		simulator.id(0)	
+		apply_single_clifford(C, simulator, 0, inverse=True)
+		
+		# Measure the state multiple times
+		measurement_results = simulator.measure_shots(shots)
+		culmu += measurement_results[0]
+		single_exps_normal[noise_strength].append(measurement_results[0] / shots)
+		# print(measurement_results)
+
+	print('Raw')
+	print(1 - culmu / shots / len(single_qubit_cliffords))
+	# single_exps_normal[noise_strength].append(measurement_results[0] / shots)
+	# print(1 - 2 * noise_strength /3)
+	print()
+
+# print(single_exps)
+# Save the dictionary to a JSON file
+# with open('data_sup.json', 'w') as file:
+#     json.dump(single_exps_sup_state0, file)
+
+# with open('data_normal.json', 'w') as file:
+#     json.dump(single_exps_normal, file)
\ No newline at end of file