Theory provides the blueprint, but software implementation drives the actual simulation of quantum systems. The mathematical behavior of qubits is described using Dirac notation and matrix algebra. To build a working quantum simulator or algorithm, these mathematical objects must be translated into code.Python, with its NumPy library, is the standard tool for this task because it handles vector and matrix operations efficiently. In this section, you will implement the linear algebra concepts covered in this chapter, transforming abstract mathematical states into concrete Python objects.Setting Up the EnvironmentBefore manipulating quantum states, import NumPy. We generally treat quantum mechanics as linear algebra over complex numbers. Therefore, when defining arrays, it is important to specify the data type as complex to support imaginary components (amplitudes).import numpy as np # Define the standard basis states |0> and |1> # We use column vectors (2 rows, 1 column) ket_0 = np.array([[1], [0]], dtype=complex) ket_1 = np.array([[0], [1]], dtype=complex) print(f"|0> state:\n{ket_0}")Notice the shape of these arrays. In quantum texts, a ket $|\psi\rangle$ is always a column vector. In NumPy, a 1D array (e.g., [1, 0]) behaves differently than a 2D column array (e.g., [[1], [0]]). Maintaining the 2D structure prevents shape mismatch errors when you start multiplying these vectors by matrices.Visualizing Matrix MultiplicationWhen a quantum gate acts on a qubit, it is mathematically a matrix multiplying a vector. It is helpful to visualize how the dimensions align during this operation using a standard matrix multiplication flow.digraph G { rankdir=LR; node [shape=box, style=filled, fontname="Arial"]; matrix [label="Quantum Gate (Matrix)\nShape: (2, 2)", fillcolor="#d0bfff", color="#7950f2"]; vector [label="Qubit State (Vector)\nShape: (2, 1)", fillcolor="#a5d8ff", color="#1c7ed6"]; result [label="New State (Vector)\nShape: (2, 1)", fillcolor="#96f2d7", color="#0ca678"]; matrix -> result [label="Operates on"]; vector -> result [label="Transforms to"]; }The flow of data when applying a quantum gate. The (2,2) matrix acts on the (2,1) vector to produce a new (2,1) vector.Calculating Probabilities with the Inner ProductThe probability of finding a system in a specific state is determined by the squared magnitude of the amplitude. Mathematically, for a state vector $| \psi \rangle$, the probability of measuring $|0\rangle$ is $|\langle 0 | \psi \rangle|^2$.To calculate this in NumPy, you need the Hermitian conjugate (conjugate transpose). In Dirac notation, if $|\psi\rangle$ is a column vector, then $\langle\psi|$ is a row vector where every element is the complex conjugate of the original.NumPy provides the .conj().T attribute chain to perform this operation. Alternatively, you can use the .H attribute if you use the matrix subclass, but standard arrays are preferred in modern NumPy.# Create a superposition state: 1/sqrt(2) * (|0> + |1>) psi = 1/np.sqrt(2) * (ket_0 + ket_1) # Calculate the Bra (row vector, complex conjugated) bra_0 = ket_0.conj().T # Calculate Inner Product <0|psi> amplitude = bra_0 @ psi # Calculate Probability (Amplitude squared) probability = np.abs(amplitude)**2 print(f"Amplitude for |0>: {amplitude[0][0]}") print(f"Probability of measuring |0>: {probability[0][0]}")Note the use of the @ operator. This is the standard operator for matrix multiplication in Python. If you use *, NumPy will attempt an element-wise multiplication, which is mathematically incorrect for this context.Modeling Quantum Gates as MatricesQuantum gates are unitary matrices. To simulate a gate, you define the matrix and multiply it by your state vector. Let us implement the Pauli-X gate, often called the "Quantum NOT" gate, which flips $|0\rangle$ to $|1\rangle$.$$X = \begin{bmatrix} 0 & 1 \ 1 & 0 \end{bmatrix}$$# Define the Pauli-X matrix X_gate = np.array([ [0, 1], [1, 0] ], dtype=complex) # Apply the gate to |0> new_state = X_gate @ ket_0 print(f"Original state:\n{ket_0}") print(f"State after X gate:\n{new_state}")If you run this code, new_state will equal [[0.+0.j], [1.+0.j]], which corresponds to the $|1\rangle$ state.Verifying Unitary MatricesA fundamental property of quantum mechanics is that operations must be unitary to preserve probability. A matrix $U$ is unitary if:$$U^\dagger U = I$$where $I$ is the identity matrix. When defining custom gates or checking your math, it is important to verify this property.def is_unitary(matrix): # Calculate U dagger (conjugate transpose) dagger = matrix.conj().T # Calculate the product U dagger * U product = dagger @ matrix # Create an identity matrix of the same size identity = np.eye(len(matrix)) # Check if they are close (handling floating point errors) return np.allclose(product, identity) print(f"Is Pauli-X unitary? {is_unitary(X_gate)}")The function np.allclose is used instead of == because floating-point arithmetic can result in minute differences (e.g., 0.9999999 instead of 1.0).Finding Eigenvalues and EigenvectorsMeasurement in quantum mechanics relates to eigenvalues. When you measure a state, the system collapses to one of the eigenstates of the operator, and the value you read out is the eigenvalue.NumPy simplifies this with np.linalg.eig.# Define the Pauli-Z matrix Z_gate = np.array([ [1, 0], [0, -1] ], dtype=complex) # Calculate eigenvalues and eigenvectors eigenvalues, eigenvectors = np.linalg.eig(Z_gate) print("Eigenvalues:", eigenvalues) print("Eigenvectors:\n", eigenvectors)For the Z-gate, you should see eigenvalues of $1$ and $-1$. These correspond to the basis states $|0\rangle$ and $|1\rangle$. This confirms that if you measure a qubit in the Z-basis (the standard computational basis), you will get outcomes distinguishing the zero state from the one state.Calculating the Outer ProductWhile the inner product (bra-ket) produces a scalar (a number), the outer product (ket-bra) produces a matrix. This is used to define projection operators. The notation $|\psi\rangle\langle\psi|$ represents the outer product.In NumPy, you can calculate this using np.outer or standard matrix multiplication of a column vector by a row vector.# Calculate projection operator |0><0| # ket_0 is (2,1) and bra_0 is (1,2) projection_0 = ket_0 @ ket_0.conj().T print("Projection Operator |0><0|:\n", projection_0)The result is the matrix: $$ \begin{bmatrix} 1 & 0 \ 0 & 0 \end{bmatrix} $$If you apply this operator to a superposition state, it projects the state onto the $|0\rangle$ subspace, effectively simulating the mathematical part of a "filtering" measurement.By mastering these NumPy operations, vector definition, matrix multiplication, conjugate transposes, and eigendecomposition, you have built a miniature quantum simulator. In the next chapters, you will use these tools to simulate more complex systems involving superposition and entanglement.