524 lines
9.2 KiB
Markdown
524 lines
9.2 KiB
Markdown
# QuTiP Analysis and Measurement
|
||
|
||
## Expectation Values
|
||
|
||
### Basic Expectation Values
|
||
|
||
```python
|
||
from qutip import *
|
||
import numpy as np
|
||
|
||
# Single operator
|
||
psi = coherent(N, 2)
|
||
n_avg = expect(num(N), psi)
|
||
|
||
# Multiple operators
|
||
ops = [num(N), destroy(N), create(N)]
|
||
results = expect(ops, psi) # Returns list
|
||
```
|
||
|
||
### Expectation Values for Density Matrices
|
||
|
||
```python
|
||
# Works with both pure states and density matrices
|
||
rho = thermal_dm(N, 2)
|
||
n_avg = expect(num(N), rho)
|
||
```
|
||
|
||
### Variance
|
||
|
||
```python
|
||
# Calculate variance of observable
|
||
var_n = variance(num(N), psi)
|
||
|
||
# Manual calculation
|
||
var_n = expect(num(N)**2, psi) - expect(num(N), psi)**2
|
||
```
|
||
|
||
### Time-Dependent Expectation Values
|
||
|
||
```python
|
||
# During time evolution
|
||
result = mesolve(H, psi0, tlist, c_ops, e_ops=[num(N)])
|
||
n_t = result.expect[0] # Array of ⟨n⟩ at each time
|
||
```
|
||
|
||
## Entropy Measures
|
||
|
||
### Von Neumann Entropy
|
||
|
||
```python
|
||
from qutip import entropy_vn
|
||
|
||
# Density matrix entropy
|
||
rho = thermal_dm(N, 2)
|
||
S = entropy_vn(rho) # Returns S = -Tr(ρ log₂ ρ)
|
||
```
|
||
|
||
### Linear Entropy
|
||
|
||
```python
|
||
from qutip import entropy_linear
|
||
|
||
# Linear entropy S_L = 1 - Tr(ρ²)
|
||
S_L = entropy_linear(rho)
|
||
```
|
||
|
||
### Entanglement Entropy
|
||
|
||
```python
|
||
# For bipartite systems
|
||
psi = bell_state('00')
|
||
rho = psi.proj()
|
||
|
||
# Trace out subsystem B to get reduced density matrix
|
||
rho_A = ptrace(rho, 0)
|
||
|
||
# Entanglement entropy
|
||
S_ent = entropy_vn(rho_A)
|
||
```
|
||
|
||
### Mutual Information
|
||
|
||
```python
|
||
from qutip import entropy_mutual
|
||
|
||
# For bipartite state ρ_AB
|
||
I = entropy_mutual(rho, [0, 1]) # I(A:B) = S(A) + S(B) - S(AB)
|
||
```
|
||
|
||
### Conditional Entropy
|
||
|
||
```python
|
||
from qutip import entropy_conditional
|
||
|
||
# S(A|B) = S(AB) - S(B)
|
||
S_cond = entropy_conditional(rho, 0) # Entropy of subsystem 0 given subsystem 1
|
||
```
|
||
|
||
## Fidelity and Distance Measures
|
||
|
||
### State Fidelity
|
||
|
||
```python
|
||
from qutip import fidelity
|
||
|
||
# Fidelity between two states
|
||
psi1 = coherent(N, 2)
|
||
psi2 = coherent(N, 2.1)
|
||
|
||
F = fidelity(psi1, psi2) # Returns value in [0, 1]
|
||
```
|
||
|
||
### Process Fidelity
|
||
|
||
```python
|
||
from qutip import process_fidelity
|
||
|
||
# Fidelity between two processes (superoperators)
|
||
U_ideal = (-1j * H * t).expm()
|
||
U_actual = mesolve(H, basis(N, 0), [0, t], c_ops).states[-1]
|
||
|
||
F_proc = process_fidelity(U_ideal, U_actual)
|
||
```
|
||
|
||
### Trace Distance
|
||
|
||
```python
|
||
from qutip import tracedist
|
||
|
||
# Trace distance D = (1/2) Tr|ρ₁ - ρ₂|
|
||
rho1 = coherent_dm(N, 2)
|
||
rho2 = thermal_dm(N, 2)
|
||
|
||
D = tracedist(rho1, rho2) # Returns value in [0, 1]
|
||
```
|
||
|
||
### Hilbert-Schmidt Distance
|
||
|
||
```python
|
||
from qutip import hilbert_dist
|
||
|
||
# Hilbert-Schmidt distance
|
||
D_HS = hilbert_dist(rho1, rho2)
|
||
```
|
||
|
||
### Bures Distance
|
||
|
||
```python
|
||
from qutip import bures_dist
|
||
|
||
# Bures distance
|
||
D_B = bures_dist(rho1, rho2)
|
||
```
|
||
|
||
### Bures Angle
|
||
|
||
```python
|
||
from qutip import bures_angle
|
||
|
||
# Bures angle
|
||
angle = bures_angle(rho1, rho2)
|
||
```
|
||
|
||
## Entanglement Measures
|
||
|
||
### Concurrence
|
||
|
||
```python
|
||
from qutip import concurrence
|
||
|
||
# For two-qubit states
|
||
psi = bell_state('00')
|
||
rho = psi.proj()
|
||
|
||
C = concurrence(rho) # C = 1 for maximally entangled states
|
||
```
|
||
|
||
### Negativity
|
||
|
||
```python
|
||
from qutip import negativity
|
||
|
||
# Negativity (partial transpose criterion)
|
||
N_ent = negativity(rho, 0) # Partial transpose w.r.t. subsystem 0
|
||
|
||
# Logarithmic negativity
|
||
from qutip import logarithmic_negativity
|
||
E_N = logarithmic_negativity(rho, 0)
|
||
```
|
||
|
||
### Entangling Power
|
||
|
||
```python
|
||
from qutip import entangling_power
|
||
|
||
# For unitary gates
|
||
U = cnot()
|
||
E_pow = entangling_power(U)
|
||
```
|
||
|
||
## Purity Measures
|
||
|
||
### Purity
|
||
|
||
```python
|
||
# Purity P = Tr(ρ²)
|
||
P = (rho * rho).tr()
|
||
|
||
# For pure states: P = 1
|
||
# For maximally mixed: P = 1/d
|
||
```
|
||
|
||
### Checking State Properties
|
||
|
||
```python
|
||
# Is state pure?
|
||
is_pure = abs((rho * rho).tr() - 1.0) < 1e-10
|
||
|
||
# Is operator Hermitian?
|
||
H.isherm
|
||
|
||
# Is operator unitary?
|
||
U.check_isunitary()
|
||
```
|
||
|
||
## Measurement
|
||
|
||
### Projective Measurement
|
||
|
||
```python
|
||
from qutip import measurement
|
||
|
||
# Measure in computational basis
|
||
psi = (basis(2, 0) + basis(2, 1)).unit()
|
||
|
||
# Perform measurement
|
||
result, state_after = measurement.measure(psi, None) # Random outcome
|
||
|
||
# Specific measurement operator
|
||
M = basis(2, 0).proj()
|
||
prob = measurement.measure_povm(psi, [M, qeye(2) - M])
|
||
```
|
||
|
||
### Measurement Statistics
|
||
|
||
```python
|
||
from qutip import measurement_statistics
|
||
|
||
# Get all possible outcomes and probabilities
|
||
outcomes, probabilities = measurement_statistics(psi, [M0, M1])
|
||
```
|
||
|
||
### Observable Measurement
|
||
|
||
```python
|
||
from qutip import measure_observable
|
||
|
||
# Measure observable and get result + collapsed state
|
||
result, state_collapsed = measure_observable(psi, sigmaz())
|
||
```
|
||
|
||
### POVM Measurements
|
||
|
||
```python
|
||
from qutip import measure_povm
|
||
|
||
# Positive Operator-Valued Measure
|
||
E_0 = Qobj([[0.8, 0], [0, 0.2]])
|
||
E_1 = Qobj([[0.2, 0], [0, 0.8]])
|
||
|
||
result, state_after = measure_povm(psi, [E_0, E_1])
|
||
```
|
||
|
||
## Coherence Measures
|
||
|
||
### l1-norm Coherence
|
||
|
||
```python
|
||
from qutip import coherence_l1norm
|
||
|
||
# l1-norm of off-diagonal elements
|
||
C_l1 = coherence_l1norm(rho)
|
||
```
|
||
|
||
## Correlation Functions
|
||
|
||
### Two-Time Correlation
|
||
|
||
```python
|
||
from qutip import correlation_2op_1t, correlation_2op_2t
|
||
|
||
# Single-time correlation ⟨A(t+τ)B(t)⟩
|
||
A = destroy(N)
|
||
B = create(N)
|
||
taulist = np.linspace(0, 10, 200)
|
||
|
||
corr = correlation_2op_1t(H, rho0, taulist, c_ops, A, B)
|
||
|
||
# Two-time correlation ⟨A(t)B(τ)⟩
|
||
tlist = np.linspace(0, 10, 100)
|
||
corr_2t = correlation_2op_2t(H, rho0, tlist, taulist, c_ops, A, B)
|
||
```
|
||
|
||
### Three-Operator Correlation
|
||
|
||
```python
|
||
from qutip import correlation_3op_1t
|
||
|
||
# ⟨A(t)B(t+τ)C(t)⟩
|
||
C_op = num(N)
|
||
corr_3 = correlation_3op_1t(H, rho0, taulist, c_ops, A, B, C_op)
|
||
```
|
||
|
||
### Four-Operator Correlation
|
||
|
||
```python
|
||
from qutip import correlation_4op_1t
|
||
|
||
# ⟨A(0)B(τ)C(τ)D(0)⟩
|
||
D_op = create(N)
|
||
corr_4 = correlation_4op_1t(H, rho0, taulist, c_ops, A, B, C_op, D_op)
|
||
```
|
||
|
||
## Spectrum Analysis
|
||
|
||
### FFT Spectrum
|
||
|
||
```python
|
||
from qutip import spectrum_correlation_fft
|
||
|
||
# Power spectrum from correlation function
|
||
w, S = spectrum_correlation_fft(taulist, corr)
|
||
```
|
||
|
||
### Direct Spectrum Calculation
|
||
|
||
```python
|
||
from qutip import spectrum
|
||
|
||
# Emission/absorption spectrum
|
||
wlist = np.linspace(0, 2, 200)
|
||
spec = spectrum(H, wlist, c_ops, A, B)
|
||
```
|
||
|
||
### Pseudo-Modes
|
||
|
||
```python
|
||
from qutip import spectrum_pi
|
||
|
||
# Spectrum with pseudo-mode decomposition
|
||
spec_pi = spectrum_pi(H, rho0, wlist, c_ops, A, B)
|
||
```
|
||
|
||
## Steady State Analysis
|
||
|
||
### Finding Steady State
|
||
|
||
```python
|
||
from qutip import steadystate
|
||
|
||
# Find steady state ∂ρ/∂t = 0
|
||
rho_ss = steadystate(H, c_ops)
|
||
|
||
# Different methods
|
||
rho_ss = steadystate(H, c_ops, method='direct') # Default
|
||
rho_ss = steadystate(H, c_ops, method='eigen') # Eigenvalue
|
||
rho_ss = steadystate(H, c_ops, method='svd') # SVD
|
||
rho_ss = steadystate(H, c_ops, method='power') # Power method
|
||
```
|
||
|
||
### Steady State Properties
|
||
|
||
```python
|
||
# Verify it's steady
|
||
L = liouvillian(H, c_ops)
|
||
assert (L * operator_to_vector(rho_ss)).norm() < 1e-10
|
||
|
||
# Compute steady-state expectation values
|
||
n_ss = expect(num(N), rho_ss)
|
||
```
|
||
|
||
## Quantum Fisher Information
|
||
|
||
```python
|
||
from qutip import qfisher
|
||
|
||
# Quantum Fisher information
|
||
F_Q = qfisher(rho, num(N)) # w.r.t. generator num(N)
|
||
```
|
||
|
||
## Matrix Analysis
|
||
|
||
### Eigenanalysis
|
||
|
||
```python
|
||
# Eigenvalues and eigenvectors
|
||
evals, ekets = H.eigenstates()
|
||
|
||
# Just eigenvalues
|
||
evals = H.eigenenergies()
|
||
|
||
# Ground state
|
||
E0, psi0 = H.groundstate()
|
||
```
|
||
|
||
### Matrix Functions
|
||
|
||
```python
|
||
# Matrix exponential
|
||
U = (H * t).expm()
|
||
|
||
# Matrix logarithm
|
||
log_rho = rho.logm()
|
||
|
||
# Matrix square root
|
||
sqrt_rho = rho.sqrtm()
|
||
|
||
# Matrix power
|
||
rho_squared = rho ** 2
|
||
```
|
||
|
||
### Singular Value Decomposition
|
||
|
||
```python
|
||
# SVD of operator
|
||
U, S, Vh = H.svd()
|
||
```
|
||
|
||
### Permutations
|
||
|
||
```python
|
||
from qutip import permute
|
||
|
||
# Permute subsystems
|
||
rho_permuted = permute(rho, [1, 0]) # Swap subsystems
|
||
```
|
||
|
||
## Partial Operations
|
||
|
||
### Partial Trace
|
||
|
||
```python
|
||
# Reduce to subsystem
|
||
rho_A = ptrace(rho_AB, 0) # Keep subsystem 0
|
||
rho_B = ptrace(rho_AB, 1) # Keep subsystem 1
|
||
|
||
# Keep multiple subsystems
|
||
rho_AC = ptrace(rho_ABC, [0, 2]) # Keep 0 and 2, trace out 1
|
||
```
|
||
|
||
### Partial Transpose
|
||
|
||
```python
|
||
from qutip import partial_transpose
|
||
|
||
# Partial transpose (for entanglement detection)
|
||
rho_pt = partial_transpose(rho, [0, 1]) # Transpose subsystem 0
|
||
|
||
# Check if entangled (PPT criterion)
|
||
evals = rho_pt.eigenenergies()
|
||
is_entangled = any(evals < -1e-10)
|
||
```
|
||
|
||
## Quantum State Tomography
|
||
|
||
### State Reconstruction
|
||
|
||
```python
|
||
from qutip_qip.tomography import state_tomography
|
||
|
||
# Prepare measurement results
|
||
# measurements = ... (experimental data)
|
||
|
||
# Reconstruct density matrix
|
||
rho_reconstructed = state_tomography(measurements, basis='Pauli')
|
||
```
|
||
|
||
### Process Tomography
|
||
|
||
```python
|
||
from qutip_qip.tomography import qpt
|
||
|
||
# Characterize quantum process
|
||
chi = qpt(U_gate, method='lstsq') # Chi matrix representation
|
||
```
|
||
|
||
## Random Quantum Objects
|
||
|
||
Useful for testing and Monte Carlo simulations.
|
||
|
||
```python
|
||
# Random state vector
|
||
psi_rand = rand_ket(N)
|
||
|
||
# Random density matrix
|
||
rho_rand = rand_dm(N)
|
||
|
||
# Random Hermitian operator
|
||
H_rand = rand_herm(N)
|
||
|
||
# Random unitary
|
||
U_rand = rand_unitary(N)
|
||
|
||
# With specific properties
|
||
rho_rank2 = rand_dm(N, rank=2) # Rank-2 density matrix
|
||
H_sparse = rand_herm(N, density=0.1) # 10% non-zero elements
|
||
```
|
||
|
||
## Useful Checks
|
||
|
||
```python
|
||
# Check if operator is Hermitian
|
||
H.isherm
|
||
|
||
# Check if state is normalized
|
||
abs(psi.norm() - 1.0) < 1e-10
|
||
|
||
# Check if density matrix is physical
|
||
rho.tr() ≈ 1 and all(rho.eigenenergies() >= 0)
|
||
|
||
# Check if operators commute
|
||
commutator(A, B).norm() < 1e-10
|
||
```
|