Quantum Fisher Information

Author

Orjan Ameye

Last Update

2025-07-30

This notebook demonstrates the computation of Quantum Fisher Information (QFI) for a driven-dissipative Kerr Parametric Oscillator (KPO) using automatic differentiation. The QFI quantifies the ultimate precision limit for parameter estimation in quantum systems.

We import the necessary packages for quantum simulations and automatic differentiation:

using QuantumToolbox   # Quantum optics simulations
using ForwardDiff      # For automatic differentiation
using LinearAlgebra    # Linear algebra operations
using CairoMakie       # Plotting
CairoMakie.activate!(type = "svg")

System Parameters and Hamiltonian

The KPO system is governed by the Hamiltonian: \[H = -p_1 a^\dagger a + K (a^\dagger)^2 a^2 - G (a^\dagger a^\dagger + a a)\]

where:

\(p_1\) is the parameter we want to estimate (detuning)
\(K\) is the Kerr nonlinearity
\(G\) is the parametric drive strength
\(\gamma\) is the decay rate

const G, K, γ = 0.002, 0.001, 0.01

const N = 20 # cutoff of the Hilbert space dimension
const a = destroy(N) # annihilation operator

const c_ops = [sqrt(γ)*a]
const ψ0 = fock(N, 0) # initial state

const tlist = range(0, 2000, 100)

function final_state(p, t)
    H = - p[1] * a' * a + K * a' * a' * a * a - G * (a' * a' + a * a)

    sol = mesolve(H, ψ0, tlist, c_ops; params = p, progress_bar = Val(false), saveat = [t])
    ρ_vec = vec(sol.states[end].data)

    # Split the real and imaginary parts of the density matrix to a real vector since ForwardDiff doesn't handle complex numbers directly
    return vcat(real.(ρ_vec), imag.(ρ_vec))
end

final_state (generic function with 1 method)

Quantum Fisher Information Calculation

The QFI is computed using the symmetric logarithmic derivative (SLD). For a density matrix \(\rho(\theta)\) parametrized by \(\theta\):

\[F_Q = \text{Tr}[\partial_\theta \rho \cdot L]\]

where \(L\) is the SLD satisfying: \(\partial_\theta \rho = \frac{1}{2}(\rho L + L \rho)\)

function compute_fisher_information(ρ, dρ)
    # Add small regularization to avoid numerical issues with zero eigenvalues
    reg = 1e-12 * I
    ρ_reg = ρ + reg

    # Solve for the symmetric logarithmic derivative L
    # dρ = (1/2)(ρL + Lρ)
    # This is a Sylvester equation: ρL + Lρ = 2*dρ
    L = sylvester(ρ_reg, ρ_reg, -2*dρ)

    # Fisher information F = Tr(dρ * L)
    F = real(tr(dρ * L))

    return F
end

compute_fisher_information (generic function with 1 method)

Automatic Differentiation Setup

We use finite differences through DifferentiationInterface.jl to compute the derivative of the quantum state with respect to the parameter. This is a key step that enables efficient QFI computation without manual derivative calculations.

# Test the system
final_state([0], 100)

# Define state function for automatic differentiation
state(p) = final_state(p, 2000)

# Compute both the state and its derivative
ρ, dρ = state([0.0]), ForwardDiff.jacobian(state, [0.0])

# Reshape back to complex matrix form
ρ = reshape(ρ[1:end÷2] + 1im * ρ[end÷2+1:end], N, N)

# Reshape the derivative back to matrix form
dρ = reshape(dρ[1:end÷2] + 1im * dρ[end÷2+1:end], N, N)

# Compute QFI at final time
qfi_final = compute_fisher_information(ρ, dρ)
println("QFI at final time: ", qfi_final)

QFI at final time: 19795.69738041981

Time Evolution of Quantum Fisher Information

Now we compute how the QFI evolves over time to understand the optimal measurement time for parameter estimation:

ts = range(0, 2000, 100)

@time QFI_t = map(ts) do t
    state(p) = final_state(p, t)
    ρ, dρ = state([0.0]), ForwardDiff.jacobian(state, [0.0])

    ρ = reshape(ρ[1:end÷2] + 1im * ρ[end÷2+1:end], N, N)
    dρ = reshape(dρ[1:end÷2] + 1im * dρ[end÷2+1:end], N, N)

    compute_fisher_information(ρ, dρ)
end

println("QFI computed for ", length(ts), " time points")

 11.930142 seconds (15.09 M allocations: 947.459 MiB, 1.65% gc time, 47.33% compilation time)
QFI computed for 100 time points

Visualization

Plot the time evolution of the Quantum Fisher Information:

fig = Figure()
ax = Axis(
 fig[1,1],
 xlabel = "Time",
 ylabel = "Quantum Fisher Information"
)
lines!(ax, ts, QFI_t)
fig

Version Information

QuantumToolbox.versioninfo()


 QuantumToolbox.jl: Quantum Toolbox in Julia
≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡
Copyright © QuTiP team 2022 and later.
Current admin team:
    Alberto Mercurio and Yi-Te Huang

Package information:
====================================
Julia              Ver. 1.11.6
QuantumToolbox     Ver. 0.34.0
SciMLOperators     Ver. 1.4.0
LinearSolve        Ver. 3.24.0
OrdinaryDiffEqCore Ver. 1.26.2

System information:
====================================
OS       : Linux (x86_64-linux-gnu)
CPU      : 4 × AMD EPYC 7763 64-Core Processor
Memory   : 15.621 GB
WORD_SIZE: 64
LIBM     : libopenlibm
LLVM     : libLLVM-16.0.6 (ORCJIT, znver3)
BLAS     : libopenblas64_.so (ilp64)
Threads  : 4 (on 4 virtual cores)

using Pkg
Pkg.status()

Status `~/work/qutip-julia-tutorials/qutip-julia-tutorials/Project.toml`
  [13f3f980] CairoMakie v0.15.4
  [f6369f11] ForwardDiff v1.0.1
  [a62dbcb7] HierarchicalEOM v2.5.6
  [6c2fb7c5] QuantumToolbox v0.34.0
  [4c0109c6] QuartoNotebookRunner v0.17.3