{
"cells": [
{
"cell_type": "markdown",
"id": "fed75000",
"metadata": {},
"source": [
"# Exploring measurement\n",
"\n",
"(Srikumar K. S., 15 Apr 2023)\n",
"\n",
"The goal of this exploratory notebook is to examine measurement as a process of\n",
"interaction between a high dimensional system and a qubit being measured."
]
},
{
"cell_type": "code",
"execution_count": 23,
"id": "33e17d3a-418a-4c7a-a9a2-745a9e08ebb1",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"MersenneTwister(1234)"
]
},
"execution_count": 23,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Random, LinearAlgebra\n",
"rng = MersenneTwister(1234)"
]
},
{
"cell_type": "markdown",
"id": "fffc7288-54dd-4d5a-b09d-f4a6b5c33f4b",
"metadata": {},
"source": [
"\n",
"$$\n",
"\\newcommand{\\ket}[1]{\\left|{#1}\\right\\rangle}\n",
"\\newcommand{\\bra}[1]{\\left\\langle{#1}\\right|}\n",
"\\newcommand{\\braket}[1]{\\left\\langle{#1}\\right\\rangle}\n",
"$$\n"
]
},
{
"cell_type": "markdown",
"id": "b24e2c32",
"metadata": {},
"source": [
"## Approach\n",
"\n",
"The idea is to see how a qubit in a random superposition of its two possibilities evolves when \"measured\" by a random observer that can distinguish between the two qubit states. Towards this, we need the ability to construct a random state, a random observable and a random observer. In this exploration, we're not initially concerned with the nature of the distribution when creating random states/observables and observers. We'll get to that later on ... if necessary.\n",
"\n",
"## Random states\n",
"\n",
"We can make any N-dimensional random vector with complex components such that $\\sum_{i}|c_i|^2 = 1$ where $c_i$ are the amplitudes of the basis states that make up the random state vector. Note that the constraint implies $|c_i|^2 \\leq 1$ for all $i$."
]
},
{
"cell_type": "code",
"execution_count": 24,
"id": "4d495dfd-43e6-4ebd-9241-5958da11c9c2",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"4-element Vector{ComplexF64}:\n",
" 0.13894291414503196 + 0.4080544353672953im\n",
" 0.10130724069115052 - 0.06104772848254118im\n",
" 0.44969950485897653 + 0.5416517697082908im\n",
" -0.4579406055551145 - 0.30801068114907093im"
]
},
"execution_count": 24,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"function random_state(rng, N)\n",
" @assert N >= 2\n",
" z = zeros(Complex, N)\n",
" for i in 1:N\n",
" a = 2.0 * (rand(rng, Float64) - 0.5)\n",
" b = 2.0 * (rand(rng, Float64) - 0.5)\n",
" z[i] = a+b*im\n",
" end\n",
" return z / sqrt(sum(abs2, z))\n",
"end\n",
"\n",
"random_state(rng, 4)"
]
},
{
"cell_type": "markdown",
"id": "1f111eee",
"metadata": {},
"source": [
"## Random density matrix\n",
"\n",
"A random density matrix must meet the criterion that its eigenvalues must be real and they must sum to unity."
]
},
{
"cell_type": "code",
"execution_count": 25,
"id": "573165b7",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"4×4 Matrix{ComplexF64}:\n",
" 0.0981778-5.40725e-18im … -0.011408+0.0348368im\n",
" -0.0607316-0.00682211im 0.0594336-0.180245im\n",
" 0.00579272-0.0772894im -0.00395489+0.00300028im\n",
" -0.011408-0.0348368im 0.191688+1.65149e-17im"
]
},
"execution_count": 25,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"\n",
"function random_ρ(rng, N; noise_floor=1e-10/N, bias=ones(Float64,N))\n",
" m = 2.0 * (rand(rng, Float64, (N,N)) + im * rand(rng, Float64, (N,N)) .- (0.5 + im*0.5))\n",
" m = 0.5*(m+m')\n",
" e = eigen(m)\n",
" eva = abs2.(e.values)\n",
" evam = maximum(eva)\n",
" ev = (eva .+ noise_floor * evam) .* bias\n",
" ev /= sum(ev)\n",
" e.vectors * Diagonal(ev) * inv(e.vectors)\n",
"end\n",
"\n",
"random_ρ(rng, 4)"
]
},
{
"cell_type": "markdown",
"id": "987b9bef",
"metadata": {},
"source": [
"## Random observers\n",
"\n",
"What is a random observer? Let's focus on observers for single qubits. A single qubit can be either in the $\\ket{0}$ state or in the $\\ket{1}$ state or in general in a superposition of the form $\\ket{q} = q_0\\ket{0}+q_1\\ket{1}$ where $|q_0|^2 + |q_1|^2 = 1$. Let's consider an \"ideal observer\" represented by a unitary operator $\\Omega$ and with $n$ dimensional internal state $\\ket{\\omega}$ which behaves like this -\n",
"\n",
"$$\\Omega\\ket{\\omega}\\ket{0} = \\ket{\\omega_0}\\ket{0}$$\n",
"and\n",
"$$\\Omega\\ket{\\omega}\\ket{1} = \\ket{\\omega_1}\\ket{1}$$\n",
"\n",
"Therefore such an $\\Omega$ will transform an arbitrary $\\ket{q}$ according to --\n",
"\n",
"$$\\Omega\\ket{\\omega}\\ket{q} = q_0\\ket{\\omega_0}\\ket{0} + q_1\\ket{\\omega_1}\\ket{1}$$\n",
"\n",
"Such an $\\Omega$ can be thought of as a combination of $\\Omega_0$ and $\\Omega_1$ where --\n",
"\n",
"$$\\Omega_0\\ket{\\omega}\\ket{q} = q_0\\ket{\\omega_0}\\ket{0} + q_1\\ket{\\omega}\\ket{1}$$\n",
"and\n",
"$$\\Omega_1\\ket{\\omega}\\ket{q} = q_0\\ket{\\omega}\\ket{0} + q_1\\ket{\\omega_1}\\ket{1}$$\n",
"\n",
"and we can write -\n",
"\n",
"$$\\Omega = \\Omega_1\\Omega_0 = \\Omega_0\\Omega_1$$\n",
"\n",
"We can therefore pick $\\Omega_0$ and $\\Omega_1$ to be unitary transformations due to random commuting observables $H_0$ and $H_1$.\n",
"\n",
"The combined observable (a self adjoint matrix) can therefore be written as --\n",
"\n",
"$$H = \\begin{bmatrix}\n",
"H_0 & \\mathbb{0} \\\\\n",
"\\mathbb{0} & H_1\n",
"\\end{bmatrix}$$"
]
},
{
"cell_type": "code",
"execution_count": 26,
"id": "7e3630af",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Eigen{ComplexF64, ComplexF64, Matrix{ComplexF64}, Vector{ComplexF64}}\n",
"values:\n",
"4-element Vector{ComplexF64}:\n",
" -0.9999999999999992 + 2.865976966174497e-17im\n",
" -0.7017467680014327 + 2.2852080636888333e-17im\n",
" 0.03676648290639895 + 3.66067959963041e-17im\n",
" 1.0000000000000004 + 3.2140735975966996e-17im\n",
"vectors:\n",
"4×4 Matrix{ComplexF64}:\n",
" 0.645644+0.0im -0.140488+0.620761im … -0.144119-0.257587im\n",
" 0.134682+0.435929im 0.666776+0.0im -0.133027-0.275052im\n",
" 0.288964+0.0572762im 0.313974-0.137567im 0.830699+0.0im\n",
" 0.535022+0.0440478im -0.144123-0.109784im -0.0465203+0.356796im"
]
},
"execution_count": 26,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"\n",
"function random_observable(rng, N; minmax=(-1.0,1.0))\n",
" m = 2.0 * (rand(rng, Float64, (N,N)) + im * rand(rng, Float64, (N,N)) .- (0.5 + im*0.5))\n",
" m = 0.5*(m+m')\n",
" e = eigen(m)\n",
" emin = minimum(real.(e.values))\n",
" emax = maximum(real.(e.values))\n",
" ev = minmax[1] .+ (e.values .- emin) * ((minmax[2] - minmax[1]) / (emax - emin))\n",
" e.vectors * Diagonal(ev) * inv(e.vectors)\n",
"end\n",
"\n",
"eigen(random_observable(rng, 4))"
]
},
{
"cell_type": "code",
"execution_count": 27,
"id": "a4060de9",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"4×4 Matrix{Complex}:\n",
" 0.623319+0.0im -0.239081-0.744522im … 0+0im\n",
" -0.239081+0.744522im -0.623319+6.1382e-18im 0+0im\n",
" 0+0im 0+0im -0.479245+0.812471im\n",
" 0+0im 0+0im -0.331987-6.7288e-17im"
]
},
"execution_count": 27,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"\n",
"function random_observer(rng, n; minmax=(-1.0,1.0))\n",
" H0 = random_observable(rng, n; minmax)\n",
" H1 = random_observable(rng, n; minmax)\n",
" Z = zeros(Complex, (n,n))\n",
" H = vcat(hcat(H0, Z), hcat(Z, H1))\n",
"end\n",
"\n",
"random_observer(rng, 2)"
]
},
{
"cell_type": "markdown",
"id": "9fdd5928",
"metadata": {},
"source": [
"## Utilities\n",
"\n",
"We need some utilities like ability to make pure state vectors, take scalar products, tensor products and such."
]
},
{
"cell_type": "code",
"execution_count": 28,
"id": "b51b4716",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"rot (generic function with 1 method)"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"\n",
"function pure_state(N, m)\n",
" z = zeros(Complex, N)\n",
" m = mod(m, N)\n",
" z[m+1] = 1.0 + 0im\n",
" return z\n",
"end\n",
"\n",
"function bra(state)\n",
" state'\n",
"end\n",
"\n",
"function braket(state1, state2)\n",
" bra(state1) * state2\n",
"end\n",
"\n",
"function ketbra(state1, state2)\n",
" state1 * bra(state2)\n",
"end\n",
"\n",
"function ρ(probs)\n",
" aprobs = abs.(probs)\n",
" nprobs = aprobs / sum(aprobs)\n",
" return Diagonal{Complex,Vector{Complex}}(nprobs)\n",
"end\n",
"\n",
"# The purity of the state represented by a sensity matrix is\n",
"# the trace of its square. The purity will be 1 if it represents\n",
"# a pure state and less than 1 if the state is a mixed state.\n",
"function purity(density)\n",
" tr(density * density)\n",
"end\n",
"\n",
"# We extend that notion of purity of a density matrix to that\n",
"# of a state, relative of the chosen set of basis states as well.\n",
"# Straightforward extension ... and I'm not even sure I'll be using it.\n",
"function spurity(state)\n",
" s = state .* conj.(state)\n",
" sum(abs2, s)\n",
"end\n",
"\n",
"# Tensor product of states is just the kron operator.\n",
"⊗ = kron\n",
"\n",
"# Exepectation of an operator given the density matrix.\n",
"function E(densityop, op)\n",
" tr(densityop * op)\n",
"end\n",
"\n",
"# Pauli matrices and identity matrices.\n",
"σx = Complex.([0 1; 1 0])\n",
"σy = Complex.([0 (-im); (im) 0])\n",
"σz = Complex.([1 0; 0 -1])\n",
"Ik(k) = Diagonal(ones(Complex, k))\n",
"I = Ik(2)\n",
"\n",
"# Pauli rotation. Computes e^{iPa}\n",
"# Assumes that pauli * pauli = I, i.e. pauli is its own inverse.\n",
"function rot(pauli, alpha)\n",
" id = Ik(size(pauli)[1])\n",
" cos(alpha) * id + im * sin(alpha) * pauli\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 29,
"id": "dc2216d0",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"2-element Vector{Complex}:\n",
" 1.0 + 0.0im\n",
" 0 + 0im"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"4×4 Diagonal{Complex, Vector{Complex}}:\n",
" 1+0im ⋅ ⋅ ⋅ \n",
" ⋅ 1+0im ⋅ ⋅ \n",
" ⋅ ⋅ 1+0im ⋅ \n",
" ⋅ ⋅ ⋅ 1+0im"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"4×4 Matrix{ComplexF64}:\n",
" 0.809017+0.0im 0.0+0.0im … 0.0+0.0im\n",
" 0.0+0.0im 0.809017+0.0im 0.0+0.587785im\n",
" 0.0+0.587785im 0.0+0.0im 0.0+0.0im\n",
" 0.0+0.0im 0.0+0.587785im 0.809017+0.0im"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"display(pure_state(2,0))\n",
"display(Ik(4))\n",
"display(rot(kron(σx,I), pi/5))"
]
},
{
"cell_type": "markdown",
"id": "221b0d6f",
"metadata": {},
"source": [
"## Evolution due to a random observable\n",
"\n",
"We take a random observable $H$ and evolve a random state $\\ket{r}$ over time\n",
"using the unitary operation $\\ket{r(t)} = e^{i2\\pi Ht}\\ket{r}$. Then we compare the result state to the\n",
"original random state as a function of time using the dot product of the two\n",
"states $\\braket{r(t)|r}$. We use $2\\pi$ as a convenience so the H makes sense as a spectrum."
]
},
{
"cell_type": "code",
"execution_count": 30,
"id": "23104103",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"\"/Users/sriku/Downloads/state-evolution.svg\""
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"function random_evolution(rng, n, T, dt)\n",
" H = random_observable(rng, n)\n",
" r = random_state(rng, n)\n",
" r0 = r\n",
" evolve = exp(im * 2 * pi * dt * H)\n",
" result = []\n",
" for t in 0:dt:T\n",
" push!(result, braket(r, r0))\n",
" r = evolve * r\n",
" end\n",
" \n",
" return real.(result), imag.(result)\n",
"end\n",
"\n",
"using Plots\n",
"p = plot(random_evolution(rng, 4, 100, 0.01)[1])\n",
"savefig(p, \"state-evolution.svg\")"
]
},
{
"cell_type": "markdown",
"id": "0bb1f98e",
"metadata": {},
"source": [
"That is kind of what you'd expect - a superposition of oscillations at a number of different frequencies."
]
},
{
"cell_type": "markdown",
"id": "f491b616",
"metadata": {},
"source": [
"## Using an observer\n",
"\n",
"Now what happens to the states if we evolve a random state similarly using an\n",
"observer. We'll try to look at the observation process as a series of small unitary\n",
"transformations and keep track of how the two states of the original measured qubit\n",
"evolve when considered along with the n states of the observer."
]
},
{
"cell_type": "code",
"execution_count": 31,
"id": "2d71dea6",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"\"/Users/sriku/Downloads/evolve-observer.svg\""
]
},
"execution_count": 31,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"function evolve_observer(rng, n, T, dt)\n",
" H = random_observer(rng, n)\n",
" w = random_state(rng, n)\n",
" q = random_state(rng, 2)\n",
" r = kron(q,w)\n",
" r0 = r\n",
" result = []\n",
" evolve = exp(im * 2 * pi * H * dt)\n",
" for t in 0:dt:T\n",
" costate_0 = r[1:n]\n",
" costate_1 = r[n+1:end]\n",
" costate_0 /= sqrt(sum(abs2, costate_0))\n",
" costate_1 /= sqrt(sum(abs2, costate_1))\n",
" push!(result, braket(costate_0, costate_1))\n",
" r = evolve * r\n",
" end\n",
" return abs.(result)\n",
"end\n",
"\n",
"evols = evolve_observer.(rng, 2 .^ (4:8), 10, 0.1)\n",
"p = plot(evols)\n",
"savefig(p, \"evolve-observer.svg\")\n"
]
},
{
"cell_type": "markdown",
"id": "d5b7b64a",
"metadata": {},
"source": [
"That looks interesting to me. While initially there is some correlation between the two parts of the\n",
"state vector $r_{1\\ldots n}$ and $r_{(n+1)\\ldots 2n}$, it rapidly declines to near 0 and lingers around 0 forever, the larger $n$ gets. For small $n$, we essentially have a small coherent superposition, but the superpositions get complex for observers with larger n. Let's split the plots into a series and see it more clearly."
]
},
{
"cell_type": "code",
"execution_count": 32,
"id": "7f9770cb",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"\"/Users/sriku/Downloads/evolve-observer-stack.svg\""
]
},
"execution_count": 32,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"tsteps = 0:0.1:100\n",
"evols = evolve_observer.(rng, 2 .^ (4:10), tsteps[end], tsteps[2]-tsteps[1])\n",
"p = plot(tsteps, [evols...], layout=(length(evols),1), size=(500,1000))\n",
"savefig(p, \"evolve-observer-stack.svg\")"
]
},
{
"cell_type": "markdown",
"id": "fa69042c-405a-4afe-b2a8-8119af855187",
"metadata": {},
"source": [
"Let's now zoom in to the first few \"seconds\" to see what happens there."
]
},
{
"cell_type": "code",
"execution_count": 33,
"id": "e7a36d51-0716-4f73-8ed7-0a08dabb4e5b",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"\"/Users/sriku/Downloads/decoherence-process.svg\""
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"zoom=1:30\n",
"p = plot(tsteps[zoom], [e[zoom] for e in evols], layout=(length(evols),1), size=(500,1000))\n",
"savefig(p, \"decoherence-process.svg\")"
]
},
{
"cell_type": "markdown",
"id": "1183835c",
"metadata": {},
"source": [
"## Observations\n",
"\n",
"The dot products of the two co-states goes to zero after about a time of $1$\n",
"for larger $n$.\n",
"\n",
"The dot product of the co-states represents the internal states taken on by the\n",
"observer for each of the two states of the qubit - i.e. the\n",
"two possibilities starting with $\\ket{\\omega}\\ket{0}$ and\n",
"$\\ket{\\omega}\\ket{1}$.\n",
"\n",
"I interpret the dot product of the co-states going to 0 as $n$ becomes larger,\n",
"as the states evolving in their own \"worlds\" and not crossing each other.\n",
"If they did cross, it would mean that at some point in the future of the\n",
"observer the qubit and the observer would become separable."
]
},
{
"cell_type": "markdown",
"id": "4983281d",
"metadata": {},
"source": [
"## Analysis\n",
"\n",
"Given the observer $H = \\begin{bmatrix} H_0 & 0 \\\\ 0 & H_1 \\end{bmatrix}$, and $H_0$ and $H_1$ are both self adjoint (i.e. $H_0 = H_0^\\dagger$ and $H_1 = H_1^\\dagger$), we have --\n",
"\n",
"$$ e^{iHt} = \\begin{bmatrix} e^{iH_0t} & 0 \\\\ 0 & e^{iH_1t} \\end{bmatrix} $$\n",
"\n",
"If we notate the combined state vector of the observer+qubit system as $\\ket{\\Psi} = \\ket{q} \\otimes \\ket{psi} = \\begin{bmatrix}q_0\\ket{\\psi} \\\\ q_1\\ket{\\psi}\\end{bmatrix}$ where $\\ket{q} = q_0\\ket{0} + q_1\\ket{1}$ (implying $|q_0|^2 + |q_1|^2 = 1$), The evolution of the joint state can be written --\n",
"\n",
"$$e^{iHt}\\ket{\\Psi} = \\begin{bmatrix}q_0e^{iH_0t}\\ket{\\psi} \\\\ q_1e^{iH_1t}\\ket{\\psi}\\end{bmatrix} = \\begin{bmatrix} q_0\\ket{\\psi_0(t)} \\\\ q_1\\ket{\\psi_1(t)}\\end{bmatrix}$$\n",
"\n",
"So if we're examining the inner product of the two co-state vectors, we get -\n",
"\n",
"$$ \\alpha = \\braket{\\psi\\left|e^{-iH_1t}e^{iH_0t}\\right|\\psi} = \\braket{\\psi\\left|e^{i(H_0-H_1)t}\\right|\\psi} = \\braket{\\psi\\left|e^{i\\Delta Ht}\\right|\\psi}$$ \n",
"\n",
"If $H_0$ and $H_1$ are very close to each other, then $\\Delta H$ ends up being small -- i.e. its eigenvalues end up small relative to the eigenvalues of $H_0$ and $H_1$ and the inner product becomes close to $1.0$. Therefore for $H$ to actually represent an observer, it must be able to differentiate between the two qubit states and therefore $\\Delta H$ itself must have random eigen values comparable in magnitude to both $H_0$ and $H_1$. \n",
"\n",
"If we then consider $\\ket{\\psi}$ in an eigen basis $\\ket{\\phi_k}$ of $\\Delta H$ with $\\braket{\\phi_i|\\phi_j} = \\delta_{ij}$ and $\\Delta H \\ket{\\phi_k} = 2\\pi f_k\\ket{\\phi_k}$, then taking $\\ket{\\psi} = \\sum_k{\\beta_k\\ket{\\phi_k}}$, we see that -\n",
"\n",
"$$ \\alpha = \\sum_k{|\\beta_k|^2e^{i(2\\pi f_kt)}} $$"
]
},
{
"cell_type": "markdown",
"id": "eb3b270f",
"metadata": {},
"source": [
"So if a) $\\ket{\\psi}$ is mixed enough in that basis and b) the \"frequencies\" $f_k$ are different enough, then for $t \\gtrapprox 0$, we see that $\\alpha \\approx 1$ and as $t$ increases, all sorts of frequencies will get mixed up and $\\alpha \\to 0$ with larger dimensions $n \\gg 1$. So the behaviour of the observer relative to the qubit states must also not have degenerate eigenvalues, or at least have sufficient number of different eigenvalues (\"frequencies\"). If we start out in a pure state of the observer system in that eigen basis though, we'll get $\\alpha = 1$. The interval during which $\\alpha \\sim 1$ (say $\\alpha > 0.5$) can then be called the \"coherence interval\" of the observer, beyond which you get \"decoherence\".\n",
"\n",
"So this is what we're seeing in the trend plot in the previous section using a random observer starting out at a random state.\n",
"\n",
"All these points correlate well with the mental model of \"observer\" that I started with, so I'm wondering whether all this should've been \"obvious\" from the get go :) I'm still irked by not being able to see why the Born rule has to hold based on this model of an observer, but that problem remains unsolved for now."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d4e4ac68-47ba-426c-b552-357f4bbb7f13",
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.8.5",
"language": "julia",
"name": "julia-1.8"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.8.5"
}
},
"nbformat": 4,
"nbformat_minor": 5
}