update

Author: mhjensen

doc/Programs/OOcodes/codes.do.txt

@@ -1812,3 +1812,255 @@ if __name__ == "__main__":
 
 !ec
 
+
+!split
+===== More codes =====
+
+!bc pycod
+import numpy as np
+
+# One-qubit gates
+def I():  return np.eye(2)
+def X():  return np.array([[0,1],[1,0]])
+def Y():  return np.array([[0,-1j],[1j,0]])
+def Z():  return np.array([[1,0],[0,-1]])
+def H():  return (1/np.sqrt(2))*np.array([[1,1],[1,-1]])
+def S():  return np.array([[1,0],[0,1j]])
+def T():  return np.array([[1,0],[0,np.exp(1j*np.pi/4)]])
+
+def Rx(theta):
+    return np.array([
+        [np.cos(theta/2), -1j*np.sin(theta/2)],
+        [-1j*np.sin(theta/2), np.cos(theta/2)]
+    ])
+
+def Ry(theta):
+    return np.array([
+        [np.cos(theta/2), -np.sin(theta/2)],
+        [np.sin(theta/2),  np.cos(theta/2)]
+    ])
+
+def Rz(theta):
+    return np.array([
+        [np.exp(-1j*theta/2), 0],
+        [0, np.exp(1j*theta/2)]
+    ])
+
+# Two-qubit gates
+def CNOT():
+    return np.array([
+        [1,0,0,0],
+        [0,1,0,0],
+        [0,0,0,1],
+        [0,0,1,0]
+    ])
+
+def CZ():
+    return np.array([
+        [1,0,0,0],
+        [0,1,0,0],
+        [0,0,1,0],
+        [0,0,0,-1]
+    ])
+
+def SWAP():
+    return np.array([
+        [1,0,0,0],
+        [0,0,1,0],
+        [0,1,0,0],
+        [0,0,0,1]
+    ])
+
+
+!ec
+
+
+!split
+===== Circuit =====
+!bc pycod
+import random
+from collections import Counter
+
+class Gate:
+    def __init__(self, matrix, targets):
+        self.matrix = matrix
+        self.targets = targets
+
+class Circuit:
+    def __init__(self, num_qubits):
+        self.n = num_qubits
+        self.reset()
+
+    def reset(self):
+        self.state = np.zeros(2**self.n, dtype=np.complex128)
+        self.state[0] = 1.0
+        self.gates = []
+
+    def add_gate(self, gate):
+        self.gates.append(gate)
+
+    def run(self):
+        for gate in self.gates:
+            self.apply_gate(gate)
+
+    def apply_gate(self, gate):
+        full_U = self.expand_gate(gate)
+        self.state = full_U @ self.state
+
+    def expand_gate(self, gate):
+        if len(gate.targets) == 1:
+            return self.tensor_product_gate(gate)
+        elif len(gate.targets) == 2:
+            return self.tensor_product_two_qubit_gate(gate)
+        else:
+            raise ValueError("Only 1 or 2 qubit gates supported.")
+
+    def tensor_product_gate(self, gate):
+        ops = [np.eye(2) for _ in range(self.n)]
+        ops[gate.targets[0]] = gate.matrix
+        U = ops[0]
+        for op in ops[1:]:
+            U = np.kron(U, op)
+        return U
+
+    def tensor_product_two_qubit_gate(self, gate):
+        ops = [np.eye(2) for _ in range(self.n)]
+        idx = sorted(gate.targets)
+        # naive expansion for generality:
+        full_U = np.eye(1)
+        for i in range(self.n):
+            if i in gate.targets:
+                continue
+            full_U = np.kron(full_U, np.eye(2))
+        full_U = np.kron(gate.matrix, full_U)
+        return full_U  # for small N this works; we can optimize later
+
+    def get_statevector(self):
+        return self.state
+
+    def get_probabilities(self):
+        return np.abs(self.state)**2
+
+    def measure(self, shots=1024):
+        probs = self.get_probabilities()
+        basis_states = [format(i, f'0{self.n}b') for i in range(2**self.n)]
+        samples = random.choices(basis_states, weights=probs, k=shots)
+        return dict(Counter(samples))
+
+
+!ec
+
+!split
+===== Noise model =====
+!bc pycod
+
+import numpy as np
+import random
+
+def bit_flip(state, p):
+    noisy_state = state.copy()
+    for i in range(len(state)):
+        if random.random() < p:
+            flipped = i ^ 1  # flips qubit 0
+            noisy_state[flipped] += state[i]
+            noisy_state[i] = 0
+    return noisy_state / np.linalg.norm(noisy_state)
+
+def depolarizing(state, p):
+    d = len(state)
+    noisy_state = (1-p) * state + p/d * np.ones(d)
+    return noisy_state / np.linalg.norm(noisy_state)
+
+def bit_flip(state, p):
+    noisy_state = state.copy()
+    for i in range(len(state)):
+        if random.random() < p:
+            flipped = i ^ 1  # flips qubit 0
+            noisy_state[flipped] += state[i]
+            noisy_state[i] = 0
+    return noisy_state / np.linalg.norm(noisy_state)
+
+def depolarizing(state, p):
+    d = len(state)
+    noisy_state = (1-p) * state + p/d * np.ones(d)
+    return noisy_state / np.linalg.norm(noisy_state)
+
+
+!ec
+
+!split
+===== Bell states =====
+!bc pycod
+from core.circuit import Circuit, Gate
+from core import gates
+
+def bell_state(label="Phi+"):
+    c = Circuit(2)
+    c.add_gate(Gate(gates.H(), [0]))
+    c.add_gate(Gate(gates.CNOT(), [0,1]))
+
+    if label == "Phi+":
+        pass
+    elif label == "Phi-":
+        c.add_gate(Gate(gates.Z(), [0]))
+    elif label == "Psi+":
+        c.add_gate(Gate(gates.X(), [1]))
+    elif label == "Psi-":
+        c.add_gate(Gate(gates.X(), [1]))
+        c.add_gate(Gate(gates.Z(), [0]))
+    else:
+        raise ValueError("Unknown Bell state label")
+
+    c.run()
+    # Assuming the Circuit class has a measure method for simulation
+    # If not, this line will need to be adjusted based on how measurement is handled
+    return c
+
+!ec
+
+!split
+===== Example =====
+!bc pycod
+import matplotlib.pyplot as plt
+import numpy as np # Add numpy import
+
+# Assuming Circuit and Gate are available in the environment from previous cell execution
+# If not, you might need to re-import them here or ensure the 'core' module is structured correctly.
+# For now, let's assume they are available after running the previous cell.
+
+labels = ["Phi+", "Phi-", "Psi+", "Psi-"]
+shots = 1000
+
+for label in labels:
+    print(f"\nBell state {label}")
+    c = bell_state(label)
+    probs = c.get_probabilities()
+    print("Statevector probabilities:", probs)
+    # The Circuit class defined earlier does NOT have a 'measure' method.
+    # You need to implement a simulation of measurement based on probabilities.
+    # Here's a possible implementation:
+    counts = {}
+    statevector = c.get_statevector()
+    num_qubits = int(np.log2(len(statevector))) # infer number of qubits
+    basis_states = [bin(i)[2:].zfill(num_qubits) for i in range(len(statevector))]
+
+    # Simulate shots based on probabilities
+    measured_indices = np.random.choice(len(statevector), size=shots, p=probs)
+    for idx in measured_indices:
+        basis_state_str = basis_states[idx]
+        counts[basis_state_str] = counts.get(basis_state_str, 0) + 1
+
+    print("Measurement counts:", counts)
+    plt.bar(counts.keys(), counts.values())
+    plt.title(f"Bell state {label}")
+    plt.show()
+
+!ec
+
+
+
+
+
+
+
+