Noise, Decoherence, $T_1$ and $T_2$

Open systems and quantum channels#

Every circuit so far was a sequence of unitaries acting on an isolated register. Physical qubits are not isolated. A superconducting transmon couples to the electromagnetic modes of its readout resonator and to two-level defects in its substrate, and a trapped ion scatters photons from stray light. Each such coupling correlates the qubit with degrees of freedom that we neither control nor measure. Once the qubit is entangled with them, its own state is the reduced density matrix obtained by tracing them out, and it is in general mixed, even if the joint evolution is unitary. This loss of coherence to an unobserved environment is decoherence . On noisy intermediate-scale (NISQ) hardware it limits the depth of every circuit, including the variational circuits of the previous modules. A model whose expressive layers need more time than the coherence time of the device produces outputs close to the maximally mixed state, whatever its parameters. The mathematical object that replaces the unitary is the quantum channel. We first characterise channels abstractly, then derive their Kraus form from a unitary model of system plus environment, and finally build the three channels that describe most hardware noise.

Definition

Quantum channel (CPTP map) and Kraus representation

\mathcal E

Theorem

From a unitary on system and environment to Kraus operators

S

Proof (7 steps)

Step 1:
$\mathrm{Tr}_E[M]=\sum_{k=0}^{r-1}\big(I\otimes\langle k\rvert_E\big)\,M\,\big(I\otimes\lvert k\rangle_E\big)$
Definition of the partial trace over $E$ in the orthonormal basis $\{\lvert k\rangle_E\}$ ; here $I\otimes\lvert k\rangle_E$ is the isometry $\lvert v\rangle\mapsto\lvert v\rangle\otimes\lvert k\rangle_E$ from $\mathbb C^d$ to $\mathbb C^d\otimes\mathbb C^r$ .
Step 2:
$\rho\otimes\lvert0\rangle\langle0\rvert_E=\big(I\otimes\lvert0\rangle_E\big)\,\rho\,\big(I\otimes\langle0\rvert_E\big)$
Both sides send $\lvert v\rangle\otimes\lvert e\rangle$ to $\rho\lvert v\rangle\otimes\lvert0\rangle\langle0\vert e\rangle$ , and product vectors span the space.
Step 3:
$\mathcal E(\rho)=\sum_k\big(I\otimes\langle k\rvert\big)U\big(I\otimes\lvert0\rangle\big)\,\rho\,\big(I\otimes\langle0\rvert\big)U^\dagger\big(I\otimes\lvert k\rangle\big)=\sum_kK_k\,\rho\,K_k^\dagger$
Insert the second line into the first; the right factor is the adjoint of $K_k$ because $\big[(I\otimes\langle k\rvert)U(I\otimes\lvert0\rangle)\big]^\dagger=(I\otimes\langle0\rvert)U^\dagger(I\otimes\lvert k\rangle)$ .
Step 4:
$(K_k)_{ij}=\langle i\rvert K_k\lvert j\rangle=\big(\langle i\rvert\otimes\langle k\rvert\big)\,U\,\big(\lvert j\rangle\otimes\lvert0\rangle\big)=U_{(i,k),(j,0)}$
Matrix elements in the product basis $\lvert i,k\rangle=\lvert i\rangle_S\otimes\lvert k\rangle_E$ : the Kraus operator $K_k$ is the $d\times d$ block of $U$ with environment output $k$ and environment input $0$ .
Step 5:
$\Big(\sum_kK_k^\dagger K_k\Big)_{jj'}=\sum_k\sum_{i}\overline{(K_k)_{ij}}\,(K_k)_{ij'}=\sum_{(i,k)}\overline{U_{(i,k),(j,0)}}\,U_{(i,k),(j',0)}$
Matrix product $(K^\dagger K)_{jj'}=\sum_i\overline{K_{ij}}K_{ij'}$ , then the previous line; the double sum over $i$ and $k$ runs over all rows $(i,k)$ of $U$ .
Step 6:
$=\sum_{(i,k)}\big(U^\dagger\big)_{(j,0),(i,k)}\,U_{(i,k),(j',0)}=\big(U^\dagger U\big)_{(j,0),(j',0)}=\delta_{jj'}$
$\overline{U_{ab}}=(U^\dagger)_{ba}$ , the sum is the matrix product of $U^\dagger$ and $U$ , and $U^\dagger U=I$ .
Step 7:
$\sum_kK_k^\dagger K_k=I_d\quad\Longrightarrow\quad\mathrm{Tr}\,\mathcal E(\rho)=\mathrm{Tr}\Big(\sum_kK_k^\dagger K_k\,\rho\Big)=\mathrm{Tr}\,\rho$
The previous line holds for all $j,j'$ ; trace preservation then follows from the cyclicity of the trace.

∎

Theorem

Every Kraus map is completely positive and trace preserving

\sum_kK_k^\dagger K_k=I

Proof (5 steps)

Step 1:
$\mathrm{Tr}\,\mathcal E(A)=\sum_k\mathrm{Tr}\big(K_kAK_k^\dagger\big)=\sum_k\mathrm{Tr}\big(K_k^\dagger K_kA\big)=\mathrm{Tr}\Big(\big(\textstyle\sum_kK_k^\dagger K_k\big)A\Big)=\mathrm{Tr}\,A$
Linearity and cyclicity of the trace, then the completeness relation; this holds for every operator $A$ , not only for density matrices.
Step 2:
$(\mathcal E\otimes\mathrm{id}_m)(A\otimes B)=\mathcal E(A)\otimes B=\sum_k(K_k\otimes I_m)(A\otimes B)(K_k\otimes I_m)^\dagger$
Definition of the extended map on product operators, and $(K\otimes I)(A\otimes B)(K^\dagger\otimes I)=KAK^\dagger\otimes B$ .
Step 3:
$(\mathcal E\otimes\mathrm{id}_m)(\Omega)=\sum_k(K_k\otimes I_m)\,\Omega\,(K_k\otimes I_m)^\dagger\quad\text{for every operator }\Omega$
Product operators span all operators on $\mathbb C^d\otimes\mathbb C^m$ and both sides are linear in $\Omega$ .
Step 4:
$\langle v\rvert(K_k\otimes I)\,\Omega\,(K_k\otimes I)^\dagger\lvert v\rangle=\langle w_k\vert\Omega\vert w_k\rangle\ge0,\qquad\lvert w_k\rangle=(K_k\otimes I)^\dagger\lvert v\rangle$
For $\Omega\succeq0$ each term is positive semidefinite: its expectation value in any vector $\lvert v\rangle$ is an expectation value of $\Omega$ .
Step 5:
$\Omega\succeq0\ \Longrightarrow\ (\mathcal E\otimes\mathrm{id}_m)(\Omega)=\sum_k(\cdots)\succeq0\qquad\text{for all }m$
A sum of positive semidefinite operators is positive semidefinite, and $m$ was arbitrary: $\mathcal E$ is completely positive.

∎

Warning

T(\rho)=\rho^{\mathsf T}

Quick check

Which of these sets of single-qubit operators is a valid Kraus set of a quantum channel?

Amplitude damping and energy relaxation#

\lvert1\rangle

U=\begin{pmatrix}1&0&0&0\\0&\sqrt{1-\gamma}&\sqrt\gamma&0\\0&-\sqrt\gamma&\sqrt{1-\gamma}&0\\0&0&0&1\end{pmatrix}=\begin{pmatrix}1&0&0&0\\0&0.8&0.6&0\\0&-0.6&0.8&0\\0&0&0&1\end{pmatrix}

Excitation-exchanging dilation of amplitude damping, γ = 0.364×4

Unitary verified: $U^\dagger U = I$ (max deviation 0.0e+0)

Derivation

Kraus operators and action of amplitude damping

Step 1:
$K_0\lvert0\rangle=\langle0\rvert_EU\lvert0\rangle\lvert0\rangle_E=\lvert0\rangle,\qquad K_0\lvert1\rangle=\langle0\rvert_EU\lvert1\rangle\lvert0\rangle_E=\sqrt{1-\gamma}\,\lvert1\rangle$
Kraus operator $K_0=\langle0\rvert_EU\lvert0\rangle_E$ of the dilation theorem, evaluated on the two basis states with the defining columns of $U$ .
Step 2:
$K_0=\lvert0\rangle\langle0\rvert+\sqrt{1-\gamma}\,\lvert1\rangle\langle1\rvert=\begin{pmatrix}1&0\\0&\sqrt{1-\gamma}\end{pmatrix}$
A linear map is determined by its action on a basis.
Step 3:
$K_1\lvert0\rangle=\langle1\rvert_EU\lvert0\rangle\lvert0\rangle_E=0,\quad K_1\lvert1\rangle=\langle1\rvert_EU\lvert1\rangle\lvert0\rangle_E=\sqrt\gamma\,\lvert0\rangle\ \Rightarrow\ K_1=\sqrt\gamma\,\lvert0\rangle\langle1\rvert=\begin{pmatrix}0&\sqrt\gamma\\0&0\end{pmatrix}$
Same computation for environment outcome $k=1$ : only the decay branch leaves an excitation in the environment.
Step 4:
$K_0^\dagger K_0+K_1^\dagger K_1=\begin{pmatrix}1&0\\0&1-\gamma\end{pmatrix}+\gamma\,\lvert1\rangle\langle0\vert0\rangle\langle1\rvert=\begin{pmatrix}1&0\\0&1-\gamma\end{pmatrix}+\begin{pmatrix}0&0\\0&\gamma\end{pmatrix}=I$
Completeness, as guaranteed by the dilation theorem: $K_1^\dagger=\sqrt\gamma\lvert1\rangle\langle0\rvert$ .
Step 5:
$K_0\rho K_0^\dagger=\begin{pmatrix}\rho_{00}&\sqrt{1-\gamma}\,\rho_{01}\\\sqrt{1-\gamma}\,\rho_{10}&(1-\gamma)\rho_{11}\end{pmatrix},\qquad K_1\rho K_1^\dagger=\gamma\,\lvert0\rangle\langle1\rvert\rho\lvert1\rangle\langle0\rvert=\begin{pmatrix}\gamma\rho_{11}&0\\0&0\end{pmatrix}$
For a diagonal $K_0=\mathrm{diag}(1,a)$ the entries are $(K_0\rho K_0^\dagger)_{ij}=a_ia_j\rho_{ij}$ with $a_0=1$ , $a_1=\sqrt{1-\gamma}$ ; and $\langle1\rvert\rho\lvert1\rangle=\rho_{11}$ .
Step 6:
$\mathcal A_\gamma(\rho)=\begin{pmatrix}\rho_{00}+\gamma\rho_{11}&\sqrt{1-\gamma}\,\rho_{01}\\\sqrt{1-\gamma}\,\rho_{10}&(1-\gamma)\rho_{11}\end{pmatrix}$
Sum of the two terms. The fraction $\gamma$ of the excited population moves to $\lvert0\rangle$ ; the coherence shrinks by $\sqrt{1-\gamma}$ , not by $1-\gamma$ .

Result

\boxed{\mathcal A_\gamma:\quad\rho_{11}\mapsto(1-\gamma)\rho_{11},\qquad\rho_{01}\mapsto\sqrt{1-\gamma}\,\rho_{01},\qquad\mathcal A_\gamma(\lvert1\rangle\langle1\rvert)=\gamma\lvert0\rangle\langle0\rvert+(1-\gamma)\lvert1\rangle\langle1\rvert}

\gamma

Theorem

Energy relaxation: γ(t) = 1 − e^{−t/T₁} and the semigroup law

\dot\rho_{11}=-\rho_{11}/T_1

Proof (7 steps)

Step 1:
$\frac{\mathrm d}{\mathrm dt}\Big(e^{t/T_1}\rho_{11}(t)\Big)=e^{t/T_1}\Big(\dot\rho_{11}+\frac{\rho_{11}}{T_1}\Big)=0\ \Longrightarrow\ \rho_{11}(t)=e^{-t/T_1}\rho_{11}(0)$
Integrating factor for the linear equation $\dot\rho_{11}=-\rho_{11}/T_1$ .
Step 2:
$\mathcal A_\gamma:\ \rho_{11}\mapsto(1-\gamma)\rho_{11}\quad\Longrightarrow\quad1-\gamma(t)=e^{-t/T_1},\qquad\gamma(t)=1-e^{-t/T_1}$
Match the channel's action on the excited population with the solution of the rate equation for every initial $\rho_{11}(0)$ .
Step 3:
$\rho_{11}\xrightarrow{\mathcal A_{\gamma_1}}(1-\gamma_1)\rho_{11}\xrightarrow{\mathcal A_{\gamma_2}}(1-\gamma_2)(1-\gamma_1)\rho_{11},\qquad\rho_{01}\mapsto\sqrt{1-\gamma_2}\sqrt{1-\gamma_1}\,\rho_{01}$
Apply the action of amplitude damping twice; each step multiplies the excited population by $1-\gamma_i$ and the coherence by $\sqrt{1-\gamma_i}$ .
Step 4:
$\rho_{00}\mapsto\rho_{00}+\gamma_1\rho_{11}+\gamma_2(1-\gamma_1)\rho_{11}=\rho_{00}+\big(1-(1-\gamma_1)(1-\gamma_2)\big)\rho_{11}$
The second channel adds $\gamma_2$ times the excited population left by the first; expand $\gamma_1+\gamma_2-\gamma_1\gamma_2=1-(1-\gamma_1)(1-\gamma_2)$ .
Step 5:
$\mathcal A_{\gamma_2}\circ\mathcal A_{\gamma_1}=\mathcal A_\gamma,\qquad1-\gamma=(1-\gamma_1)(1-\gamma_2)$
The three entries are those of $\mathcal A_\gamma$ , since $\sqrt{1-\gamma_2}\sqrt{1-\gamma_1}=\sqrt{(1-\gamma_1)(1-\gamma_2)}$ for nonnegative factors.
Step 6:
$\big(1-\gamma(t_1)\big)\big(1-\gamma(t_2)\big)=e^{-t_1/T_1}e^{-t_2/T_1}=e^{-(t_1+t_2)/T_1}=1-\gamma(t_1+t_2)$
Insert $\gamma(t)$ : two idle periods compose to one idle period of the total length, the semigroup property of a Markovian process.
Step 7:
$\rho_{01}(t)=\sqrt{1-\gamma(t)}\,\rho_{01}(0)=e^{-t/(2T_1)}\rho_{01}(0)$
Coherence factor of the channel with $\gamma(t)$ inserted: energy relaxation alone destroys coherence at half the population rate.

∎

T_1=80\,\mu\mathrm s

Quick check

T_1=100\,\mu\mathrm s

Pure dephasing and phase flips#

\lvert0\rangle

Derivation

Kraus operators and action of phase damping

Step 1:
$U^\dagger U=\lvert0\rangle\langle0\rvert\otimes I+\lvert1\rangle\langle1\rvert\otimes R_Y(\vartheta)^\dagger R_Y(\vartheta)=I$
$U$ is block diagonal in the system basis with unitary blocks $I$ and $R_Y(\vartheta)$ , and the cross terms vanish because $\langle0\vert1\rangle=0$ .
Step 2:
$K_0=\langle0\rvert_EU\lvert0\rangle_E=\lvert0\rangle\langle0\rvert\,\langle0\vert0\rangle+\lvert1\rangle\langle1\rvert\,\langle0\rvert R_Y(\vartheta)\lvert0\rangle=\begin{pmatrix}1&0\\0&\cos\frac\vartheta2\end{pmatrix}=\begin{pmatrix}1&0\\0&\sqrt{1-\lambda}\end{pmatrix}$
Kraus operator for environment outcome $0$ , with $R_Y(\vartheta)\lvert0\rangle=\cos\frac\vartheta2\lvert0\rangle+\sin\frac\vartheta2\lvert1\rangle$ .
Step 3:
$K_1=\langle1\rvert_EU\lvert0\rangle_E=\lvert1\rangle\langle1\rvert\,\langle1\rvert R_Y(\vartheta)\lvert0\rangle=\begin{pmatrix}0&0\\0&\sin\frac\vartheta2\end{pmatrix}=\begin{pmatrix}0&0\\0&\sqrt\lambda\end{pmatrix}$
Environment outcome $1$ : only the branch in which the system is $\lvert1\rangle$ rotates the environment.
Step 4:
$K_0^\dagger K_0+K_1^\dagger K_1=\mathrm{diag}(1,1-\lambda)+\mathrm{diag}(0,\lambda)=I$
Completeness, as required by the dilation theorem.
Step 5:
$K_0\rho K_0^\dagger=\begin{pmatrix}\rho_{00}&\sqrt{1-\lambda}\,\rho_{01}\\\sqrt{1-\lambda}\,\rho_{10}&(1-\lambda)\rho_{11}\end{pmatrix},\qquad K_1\rho K_1^\dagger=\begin{pmatrix}0&0\\0&\lambda\rho_{11}\end{pmatrix}$
Both Kraus operators are diagonal, so $(K\rho K^\dagger)_{ij}=k_ik_j\rho_{ij}$ with $k=(1,\sqrt{1-\lambda})$ and $k=(0,\sqrt\lambda)$ .
Step 6:
$\mathcal P_\lambda(\rho)=\begin{pmatrix}\rho_{00}&\sqrt{1-\lambda}\,\rho_{01}\\\sqrt{1-\lambda}\,\rho_{10}&\rho_{11}\end{pmatrix}$
Sum of the two terms: $(1-\lambda)\rho_{11}+\lambda\rho_{11}=\rho_{11}$ , so the populations are unchanged and only the coherence is multiplied by $\sqrt{1-\lambda}$ .
Step 7:
$\sqrt{1-\lambda}=\langle e_1\vert e_0\rangle,\qquad\lvert e_0\rangle=\lvert0\rangle_E,\quad\lvert e_1\rangle=R_Y(\vartheta)\lvert0\rangle_E$
Interpretation: $U(\alpha\lvert0\rangle+\beta\lvert1\rangle)\lvert0\rangle_E=\alpha\lvert0\rangle\lvert e_0\rangle+\beta\lvert1\rangle\lvert e_1\rangle$ , whose reduced coherence is $\alpha\bar\beta\,\langle e_1\vert e_0\rangle$ . The more distinguishable the environment states, the smaller the coherence; orthogonal states ( $\lambda=1$ ) erase it.

Result

\boxed{\mathcal P_\lambda:\quad\rho_{00},\rho_{11}\ \text{unchanged},\qquad\rho_{01}\mapsto\sqrt{1-\lambda}\,\rho_{01}}

Circuit · 2 qubits

\lambda=0.64

\lvert\psi\rangle = \tfrac{1}{\sqrt{2}}\lvert 00\rangle + 0.4243\lvert 10\rangle + 0.5657\lvert 11\rangle

Theorem

Phase flip = phase damping (unitary freedom of Kraus operators)

u

Proof (5 steps)

Step 1:
$\sum_jL_j\rho L_j^\dagger=\sum_j\sum_{k,l}u_{jk}\,\overline{u_{jl}}\,K_k\rho K_l^\dagger=\sum_{k,l}\big(u^\dagger u\big)_{lk}\,K_k\rho K_l^\dagger=\sum_kK_k\rho K_k^\dagger$
Expand $L_j$ and $L_j^\dagger=\sum_l\overline{u_{jl}}K_l^\dagger$ ; $\sum_j\overline{u_{jl}}u_{jk}=(u^\dagger u)_{lk}=\delta_{lk}$ because $u$ is unitary. This proves (i).
Step 2:
$Z\rho Z=\begin{pmatrix}\rho_{00}&-\rho_{01}\\-\rho_{10}&\rho_{11}\end{pmatrix}\ \Longrightarrow\ (1-p)\rho+pZ\rho Z=\begin{pmatrix}\rho_{00}&(1-2p)\rho_{01}\\(1-2p)\rho_{10}&\rho_{11}\end{pmatrix}$
Conjugation by $Z=\mathrm{diag}(1,-1)$ flips the sign of the off-diagonal entries; then add the two terms.
Step 3:
$\sqrt{1-\lambda}=1-2p\ \ (0\le p\le\tfrac12)\quad\Longleftrightarrow\quad\lambda=1-(1-2p)^2=4p(1-p)$
Comparison with the action of $\mathcal P_\lambda$ : both channels fix the populations and multiply the coherence by the same factor. The square root is nonnegative, which requires $p\le\tfrac12$ .
Step 4:
$\sqrt{1-p}\,K_0+\sqrt p\,K_1=\mathrm{diag}\Big(\sqrt{1-p},\ \sqrt{1-p}\,(1-2p)+\sqrt p\cdot2\sqrt{p(1-p)}\Big)=\sqrt{1-p}\,I$
First row of $u$ applied to $K_0=\mathrm{diag}(1,1-2p)$ , $K_1=\mathrm{diag}\big(0,2\sqrt{p(1-p)}\big)$ (phase damping with $\sqrt\lambda=2\sqrt{p(1-p)}$ ); the second entry is $\sqrt{1-p}\,(1-2p+2p)$ .
Step 5:
$\sqrt p\,K_0-\sqrt{1-p}\,K_1=\mathrm{diag}\Big(\sqrt p,\ \sqrt p\,(1-2p)-2(1-p)\sqrt p\Big)=\mathrm{diag}\big(\sqrt p,-\sqrt p\big)=\sqrt p\,Z$
Second row of $u$ ; $u$ is real with orthonormal rows, so it is unitary, as (i) requires.

∎

Quick check

Z

Bloch-vector relaxation: T₁, T₂ and T_φ#

\rho=\tfrac12\big(I+r_xX+r_yY+r_zZ\big)

Theorem

Bloch-vector dynamics under amplitude and phase damping

\mathcal A_\gamma

Proof (6 steps)

Step 1:
$\rho=\frac12\begin{pmatrix}1+r_z&r_x-ir_y\\r_x+ir_y&1-r_z\end{pmatrix}$
The Pauli expansion of a qubit density matrix written out entry by entry.
Step 2:
$\frac{1+r_z'}2=\frac{1+r_z}2+\gamma\,\frac{1-r_z}2\quad\Longrightarrow\quad r_z'=r_z+\gamma(1-r_z)=1-(1-\gamma)(1-r_z)$
Amplitude damping maps $\rho_{00}\mapsto\rho_{00}+\gamma\rho_{11}$ ; solve for the new $r_z$ and rearrange.
Step 3:
$\frac{r_x'-ir_y'}2=\sqrt{1-\gamma}\,\frac{r_x-ir_y}2\quad\Longrightarrow\quad r_x'=\sqrt{1-\gamma}\,r_x,\quad r_y'=\sqrt{1-\gamma}\,r_y$
Amplitude damping multiplies $\rho_{01}$ by the real factor $\sqrt{1-\gamma}$ ; compare real and imaginary parts.
Step 4:
$\mathcal P_\lambda:\ \rho_{00},\rho_{11}\ \text{fixed}\Rightarrow r_z'=r_z;\qquad\rho_{01}\mapsto\sqrt{1-\lambda}\,\rho_{01}\Rightarrow r_{x,y}'=\sqrt{1-\lambda}\,r_{x,y}$
Same reading of the phase-damping action.
Step 5:
$r_z(t)=1-e^{-t/T_1}\big(1-r_z(0)\big)$
Insert $1-\gamma(t)=e^{-t/T_1}$ ; phase damping does not change $r_z$ . The fixed point is $r_z=1$ , the ground state $\lvert0\rangle$ .
Step 6:
$r_{x,y}(t)=e^{-t/(2T_1)}\,e^{-t/T_\varphi}\,r_{x,y}(0)=e^{-t\left(\frac1{2T_1}+\frac1{T_\varphi}\right)}r_{x,y}(0)=e^{-t/T_2}\,r_{x,y}(0)$
The transverse factors of the two channels multiply; $\sqrt{1-\gamma(t)}=e^{-t/(2T_1)}$ by the semigroup theorem; the exponents add, which defines $T_2$ .

∎

Theorem

Amplitude and phase damping commute, and T₂ ≤ 2T₁

\mathcal P_\lambda\circ\mathcal A_\gamma=\mathcal A_\gamma\circ\mathcal P_\lambda

Proof (5 steps)

Step 1:
$\mathcal P_\lambda\big(\mathcal A_\gamma(\rho)\big)=\begin{pmatrix}\rho_{00}+\gamma\rho_{11}&\sqrt{1-\lambda}\sqrt{1-\gamma}\,\rho_{01}\\\sqrt{1-\lambda}\sqrt{1-\gamma}\,\rho_{10}&(1-\gamma)\rho_{11}\end{pmatrix}$
Apply amplitude damping, then phase damping, which keeps the populations and rescales the coherence.
Step 2:
$\mathcal A_\gamma\big(\mathcal P_\lambda(\rho)\big)=\begin{pmatrix}\rho_{00}+\gamma\rho_{11}&\sqrt{1-\gamma}\sqrt{1-\lambda}\,\rho_{01}\\\sqrt{1-\gamma}\sqrt{1-\lambda}\,\rho_{10}&(1-\gamma)\rho_{11}\end{pmatrix}$
Phase damping first leaves the populations unchanged, and amplitude damping acts on the populations through $\rho_{11}$ alone; the coherence factors are the same numbers in the opposite order.
Step 3:
$\frac1{T_2}=\frac1{2T_1}+\frac1{T_\varphi},\qquad\frac1{T_\varphi}\ge0$
Definition of $T_2$ from the Bloch theorem; a dephasing rate cannot be negative because $\sqrt{1-\lambda}\le1$ .
Step 4:
$\frac1{T_2}\ge\frac1{2T_1}\quad\Longleftrightarrow\quad T_2\le2T_1,\qquad T_2=2T_1\iff\frac1{T_\varphi}=0$
Take reciprocals of positive numbers, which reverses the inequality.
Step 5:
$\frac1{T_\varphi}=\frac1{T_2}-\frac1{2T_1}=\frac{2T_1-T_2}{2T_1T_2}\ \ge0\quad\text{iff}\quad T_2\le2T_1$
Solve the defining relation for the pure-dephasing rate; a measurement with $T_2>2T_1$ would contradict the model and indicates a non-Markovian or miscalibrated experiment.

∎

Derivation

The thermal-relaxation channel of the simulator

Step 1:
$\gamma=1-e^{-t/T_1},\qquad\lambda=1-e^{-2t/T_\varphi}=1-e^{-2t\left(\frac1{T_2}-\frac1{2T_1}\right)}\quad\Longrightarrow\quad\sqrt{1-\lambda}=e^{-t/T_\varphi}$
The engine's parameters for thermalRelaxation $(T_1,T_2,t)$ with $0<T_2\le2T_1$ : amplitude damping for the population decay, phase damping for the remaining coherence decay.
Step 2:
$P_0A_0=\begin{pmatrix}1&0\\0&\sqrt{1-\lambda}\sqrt{1-\gamma}\end{pmatrix}=\begin{pmatrix}1&0\\0&e^{-t/T_2}\end{pmatrix},\qquad P_0A_1=\begin{pmatrix}1&0\\0&\sqrt{1-\lambda}\end{pmatrix}\begin{pmatrix}0&\sqrt\gamma\\0&0\end{pmatrix}=\begin{pmatrix}0&\sqrt\gamma\\0&0\end{pmatrix}$
Products of the phase-damping Kraus operators $P_i$ and amplitude-damping Kraus operators $A_j$ (the Kraus set of a composition); $e^{-t/T_\varphi}e^{-t/(2T_1)}=e^{-t/T_2}$ .
Step 3:
$P_1A_0=\begin{pmatrix}0&0\\0&\sqrt{\lambda(1-\gamma)}\end{pmatrix},\qquad P_1A_1=\begin{pmatrix}0&0\\0&\sqrt\lambda\end{pmatrix}\begin{pmatrix}0&\sqrt\gamma\\0&0\end{pmatrix}=0$
The last product vanishes and is dropped, so the channel has three Kraus operators when $\lambda>0$ and two when $T_2=2T_1$ (then $\lambda=0$ and $P_1A_0=0$ too).
Step 4:
$\sum K^\dagger K=\mathrm{diag}\big(1,(1-\lambda)(1-\gamma)\big)+\mathrm{diag}(0,\gamma)+\mathrm{diag}\big(0,\lambda(1-\gamma)\big)=I$
Completeness: $(1-\gamma)(1-\lambda+\lambda)+\gamma=1$ .
Step 5:
$\rho\mapsto\begin{pmatrix}\rho_{00}+\gamma\rho_{11}&e^{-t/T_2}\rho_{01}\\e^{-t/T_2}\rho_{10}&e^{-t/T_1}\rho_{11}\end{pmatrix}$
Sum of the three conjugations: the excited population collects $(1-\lambda)(1-\gamma)\rho_{11}+\lambda(1-\gamma)\rho_{11}=(1-\gamma)\rho_{11}$ , the ground population gains $\gamma\rho_{11}$ , and the coherence decays with $T_2$ , exactly the Bloch laws.

Result

\boxed{\text{thermalRelaxation}(T_1,T_2,t)=\mathcal P_{\lambda(t)}\circ\mathcal A_{\gamma(t)}:\qquad\rho_{11}(t)=e^{-t/T_1}\rho_{11}(0),\quad\rho_{01}(t)=e^{-t/T_2}\rho_{01}(0)}

Decoherence on the Bloch sphereDecoherenceBloch

Bloch-vector trajectories under T₁ relaxation and T₂ dephasing with sliders, together with the density-matrix entries over time.

r_z(t) = 1 - \bigl(1 - r_z(0)\bigr)e^{-t/T_1},\qquad r_\perp(t) = r_\perp(0)\,e^{-t/T_2}

\varphi(t) = \varphi_0 + \Delta t,\qquad \frac{1}{T_\varphi} = \frac{1}{T_2} - \frac{1}{2T_1} \ge 0

x–z plane (populations)

x–y plane (coherence)

○ initial state · ● state at t · line: trajectory on [0, 150] μs. Every trajectory ends at |0⟩ (r = (0, 0, 1)): T₁ pulls r_z up, T₂ shrinks the transverse part.

Components against time

$r_x$
$r_y$
$r_z$
$|\vec r|$
$e^{-t/T_1}$
$e^{-t/T_2}$

t: 0.000 μs
r(t) = (x, y, z): (1.000, 0.000, 0.000)
|r(t)|: 1.0000
P(1) = (1 − r_z)/2: 0.5000
γ(t) = 1 − e^(−t/T₁): 0.0000
T_φ = (1/T₂ − 1/(2T₁))⁻¹: 42.857 μs
|ρ₀₁| = r_⊥/2: 0.5000
purity Tr ρ²: 1.0000
Kraus channel vs closed form: …

\rho(t) = \begin{pmatrix} 0.5000 & 0.5000 \\ 0.5000 & 0.5000 \end{pmatrix}

Populations ρ₀₀ = 0.5000, ρ₁₁ = 0.5000; coherence ρ₀₁ = (r_x − i r_y)/2 decays as e^(−t/T₂).

Relaxation times and time

T₁ (energy relaxation)50.00 μs

T₂ (coherence)30.00 μs

T₂ ≤ 2T₁: relaxation alone already destroys coherence at rate 1/(2T₁), so 1/T₂ = 1/(2T₁) + 1/T_φ cannot be smaller.

Time t0.00 μs

Detuning Δ0.0000 rad/μs

Initial state

Polar angle θ1.571

Azimuth φ₀0.000

r(0) = (1.000, 0.000, 0.000)

T_1=50\,\mu\mathrm s

Decoherence on the Bloch sphereDecoherenceBloch

Bloch-vector trajectories under T₁ relaxation and T₂ dephasing with sliders, together with the density-matrix entries over time.

r_z(t) = 1 - \bigl(1 - r_z(0)\bigr)e^{-t/T_1},\qquad r_\perp(t) = r_\perp(0)\,e^{-t/T_2}

\varphi(t) = \varphi_0 + \Delta t,\qquad \frac{1}{T_\varphi} = \frac{1}{T_2} - \frac{1}{2T_1} \ge 0

x–z plane (populations)

x–y plane (coherence)

○ initial state · ● state at t · line: trajectory on [0, 150] μs. Every trajectory ends at |0⟩ (r = (0, 0, 1)): T₁ pulls r_z up, T₂ shrinks the transverse part.

Components against time

$r_x$
$r_y$
$r_z$
$|\vec r|$
$e^{-t/T_1}$
$e^{-t/T_2}$

t: 0.000 μs
r(t) = (x, y, z): (1.000, 0.000, 0.000)
|r(t)|: 1.0000
P(1) = (1 − r_z)/2: 0.5000
γ(t) = 1 − e^(−t/T₁): 0.0000
T_φ = (1/T₂ − 1/(2T₁))⁻¹: ∞ (T₂ = 2T₁)
|ρ₀₁| = r_⊥/2: 0.5000
purity Tr ρ²: 1.0000
Kraus channel vs closed form: …

\rho(t) = \begin{pmatrix} 0.5000 & 0.5000 \\ 0.5000 & 0.5000 \end{pmatrix}

Populations ρ₀₀ = 0.5000, ρ₁₁ = 0.5000; coherence ρ₀₁ = (r_x − i r_y)/2 decays as e^(−t/T₂).

Relaxation times and time

T₁ (energy relaxation)50.00 μs

T₂ (coherence)100.00 μs

T₂ ≤ 2T₁: relaxation alone already destroys coherence at rate 1/(2T₁), so 1/T₂ = 1/(2T₁) + 1/T_φ cannot be smaller.

Time t0.00 μs

Detuning Δ0.0000 rad/μs

Initial state

Polar angle θ1.571

Azimuth φ₀0.000

r(0) = (1.000, 0.000, 0.000)

T_2=20\,\mu\mathrm s

Decoherence on the Bloch sphereDecoherenceBloch

Bloch-vector trajectories under T₁ relaxation and T₂ dephasing with sliders, together with the density-matrix entries over time.

r_z(t) = 1 - \bigl(1 - r_z(0)\bigr)e^{-t/T_1},\qquad r_\perp(t) = r_\perp(0)\,e^{-t/T_2}

\varphi(t) = \varphi_0 + \Delta t,\qquad \frac{1}{T_\varphi} = \frac{1}{T_2} - \frac{1}{2T_1} \ge 0

x–z plane (populations)

x–y plane (coherence)

○ initial state · ● state at t · line: trajectory on [0, 150] μs. Every trajectory ends at |0⟩ (r = (0, 0, 1)): T₁ pulls r_z up, T₂ shrinks the transverse part.

Components against time

$r_x$
$r_y$
$r_z$
$|\vec r|$
$e^{-t/T_1}$
$e^{-t/T_2}$

t: 0.000 μs
r(t) = (x, y, z): (1.000, 0.000, 0.000)
|r(t)|: 1.0000
P(1) = (1 − r_z)/2: 0.5000
γ(t) = 1 − e^(−t/T₁): 0.0000
T_φ = (1/T₂ − 1/(2T₁))⁻¹: 25.000 μs
|ρ₀₁| = r_⊥/2: 0.5000
purity Tr ρ²: 1.0000

\rho(t) = \begin{pmatrix} 0.5000 & 0.5000 \\ 0.5000 & 0.5000 \end{pmatrix}

Populations ρ₀₀ = 0.5000, ρ₁₁ = 0.5000; coherence ρ₀₁ = (r_x − i r_y)/2 decays as e^(−t/T₂).

Relaxation times and time

T₁ (energy relaxation)50.00 μs

T₂ (coherence)20.00 μs

T₂ ≤ 2T₁: relaxation alone already destroys coherence at rate 1/(2T₁), so 1/T₂ = 1/(2T₁) + 1/T_φ cannot be smaller.

Time t0.00 μs

Detuning Δ0.2000 rad/μs

Initial state

Polar angle θ1.571

Azimuth φ₀0.000

r(0) = (1.000, 0.000, 0.000)

On hardware

T_1

Quick check

T_1=100\,\mu\mathrm s

The depolarizing channel#

X

Lemma

Pauli twirl: Σ_P PAP = 2 Tr(A) I

2\times2

Proof (5 steps)

Step 1:
$A=\begin{pmatrix}a&b\\c&d\end{pmatrix},\qquad\mathrm{Tr}A=a+d$
An arbitrary complex $2\times2$ matrix.
Step 2:
$XAX=\begin{pmatrix}0&1\\1&0\end{pmatrix}\begin{pmatrix}a&b\\c&d\end{pmatrix}\begin{pmatrix}0&1\\1&0\end{pmatrix}=\begin{pmatrix}d&c\\b&a\end{pmatrix}$
Left multiplication by $X$ swaps the rows, right multiplication swaps the columns.
Step 3:
$YA=\begin{pmatrix}-ic&-id\\ia&ib\end{pmatrix},\qquad YAY=\begin{pmatrix}-ic&-id\\ia&ib\end{pmatrix}\begin{pmatrix}0&-i\\i&0\end{pmatrix}=\begin{pmatrix}d&-c\\-b&a\end{pmatrix}$
Direct multiplication with $Y=\begin{pmatrix}0&-i\\i&0\end{pmatrix}$ , using $(-i)(i)=1$ and $i\cdot i=-1$ .
Step 4:
$ZAZ=\begin{pmatrix}1&0\\0&-1\end{pmatrix}\begin{pmatrix}a&b\\c&d\end{pmatrix}\begin{pmatrix}1&0\\0&-1\end{pmatrix}=\begin{pmatrix}a&-b\\-c&d\end{pmatrix}$
Conjugation by $Z$ flips the sign of the off-diagonal entries.
Step 5:
$A+XAX+YAY+ZAZ=\begin{pmatrix}a+d+d+a&b+c-c-b\\c+b-b-c&d+a+a+d\end{pmatrix}=2(a+d)\,I=2\,\mathrm{Tr}(A)\,I$
Add the four matrices entry by entry: the off-diagonal entries cancel in pairs.

∎

Derivation

The depolarizing channel shrinks the Bloch vector by 1 − 4p/3

Step 1:
$X\rho X+Y\rho Y+Z\rho Z=2\,\mathrm{Tr}(\rho)\,I-\rho=2I-\rho$
The Pauli-twirl lemma with $A=\rho$ , minus the term $I\rho I$ ; $\mathrm{Tr}\rho=1$ .
Step 2:
$\mathcal D_p(\rho)=(1-p)\rho+\frac p3\,(2I-\rho)=\Big(1-\frac{4p}3\Big)\rho+\frac{2p}3\,I=\Big(1-\frac{4p}3\Big)\rho+\frac{4p}3\,\frac I2$
Substitute into the Kraus form and collect the terms in $\rho$ : $1-p-\tfrac p3=1-\tfrac{4p}3$ .
Step 3:
$\mathcal D_p\Big(\tfrac12(I+\vec r\cdot\vec\sigma)\Big)=\tfrac12\Big(I+\big(1-\tfrac{4p}3\big)\,\vec r\cdot\vec\sigma\Big)\quad\Longrightarrow\quad\vec r\mapsto f\,\vec r,\quad f=1-\frac{4p}3$
Insert the Bloch form: the identity part is preserved and every Bloch component is multiplied by the same factor, so the ball shrinks uniformly towards its centre.
Step 4:
$\text{Qiskit: }\ \mathcal E_\lambda(\rho)=(1-\lambda)\rho+\lambda\,\frac I2\quad\Longrightarrow\quad\lambda=\frac{4p}3$
Qiskit's depolarizing_error with parameter $\lambda$ mixes with the maximally mixed state with weight $\lambda$ ; equate the coefficients of $\rho$ . The same physical error has $\lambda=\tfrac43p$ , e.g. $p=0.3$ is $\lambda=0.4$ .
Step 5:
$p=\tfrac34\ (\lambda=1):\ \mathcal D_{3/4}(\rho)=\frac I2\ \text{for every }\rho;\qquad p=1:\ \vec r\mapsto-\tfrac13\vec r$
At $p=\tfrac34$ the factor vanishes and every input is mapped to the maximally mixed state (complete depolarisation). For $p>\tfrac34$ the factor is negative: the channel is still completely positive (it has a Kraus form), but it inverts the Bloch vector while shrinking it.

Result

\boxed{\mathcal D_p(A)=f\,A+(1-f)\,\mathrm{Tr}(A)\,\frac I2\ \ \text{for every operator }A,\qquad f=1-\frac{4p}3}

\lvert+\rangle

Quick check

\lambda=0.64

Noisy entanglement: the Bell-state fidelity#

\lvert\Phi^+\rangle=\tfrac1{\sqrt2}(\lvert00\rangle+\lvert11\rangle)

Lemma

Pauli expansion of the Bell projector

\lvert\Phi^+\rangle\langle\Phi^+\rvert=\tfrac14\big(I\otimes I+X\otimes X-Y\otimes Y+Z\otimes Z\big)

Proof (5 steps)

Step 1:
$\lvert\Phi^+\rangle\langle\Phi^+\rvert=\tfrac12\big(\lvert00\rangle\langle00\rvert+\lvert00\rangle\langle11\rvert+\lvert11\rangle\langle00\rvert+\lvert11\rangle\langle11\rvert\big)=\frac12\begin{pmatrix}1&0&0&1\\0&0&0&0\\0&0&0&0\\1&0&0&1\end{pmatrix}$
Expand the outer product in the basis $\lvert00\rangle,\lvert01\rangle,\lvert10\rangle,\lvert11\rangle$ .
Step 2:
$Z\otimes Z=\mathrm{diag}(1,-1,-1,1)\ \Longrightarrow\ I\otimes I+Z\otimes Z=\mathrm{diag}(2,0,0,2)$
$Z\otimes Z\lvert ab\rangle=(-1)^{a+b}\lvert ab\rangle$ .
Step 3:
$X\otimes X\lvert ab\rangle=\lvert\bar a\bar b\rangle,\qquad Y\otimes Y\lvert00\rangle=-\lvert11\rangle,\ \ Y\otimes Y\lvert11\rangle=-\lvert00\rangle,\ \ Y\otimes Y\lvert01\rangle=\lvert10\rangle,\ \ Y\otimes Y\lvert10\rangle=\lvert01\rangle$
$Y\lvert0\rangle=i\lvert1\rangle$ and $Y\lvert1\rangle=-i\lvert0\rangle$ , so the phases are $i\cdot i=-1$ , $(-i)(-i)=-1$ and $i(-i)=1$ .
Step 4:
$X\otimes X-Y\otimes Y=2\big(\lvert00\rangle\langle11\rvert+\lvert11\rangle\langle00\rvert\big)$
On the pair $\lvert00\rangle,\lvert11\rangle$ the matrix elements are $1-(-1)=2$ ; on the pair $\lvert01\rangle,\lvert10\rangle$ they are $1-1=0$ .
Step 5:
$\tfrac14\big(II+ZZ+XX-YY\big)=\tfrac14\Big(2\lvert00\rangle\langle00\rvert+2\lvert11\rangle\langle11\rvert+2\lvert00\rangle\langle11\rvert+2\lvert11\rangle\langle00\rvert\Big)=\lvert\Phi^+\rangle\langle\Phi^+\rvert$
Add the two previous results and compare with the first line.

∎

Theorem

Bell pair under local depolarizing noise

f=1-\tfrac{4p}3

Proof (8 steps)

Step 1:
$\mathcal D_p(I)=f\,I+(1-f)\cdot2\cdot\frac I2=I,\qquad\mathcal D_p(\sigma_j)=f\,\sigma_j+(1-f)\,\mathrm{Tr}(\sigma_j)\,\frac I2=f\,\sigma_j$
The operator form $\mathcal D_p(A)=fA+(1-f)\mathrm{Tr}(A)\tfrac I2$ of the depolarizing derivation, with $\mathrm{Tr}\,I=2$ and $\mathrm{Tr}\,\sigma_j=0$ .
Step 2:
$(\mathcal D_p\otimes\mathcal D_p)(A\otimes B)=\sum_{a,b}(K_a\otimes K_b)(A\otimes B)(K_a\otimes K_b)^\dagger=\mathcal D_p(A)\otimes\mathcal D_p(B)$
The Kraus operators of the product channel are $K_a\otimes K_b$ ; the double sum factorises.
Step 3:
$II\mapsto II,\qquad XX\mapsto f^2XX,\qquad YY\mapsto f^2YY,\qquad ZZ\mapsto f^2ZZ$
Each non-identity factor contributes one factor $f$ .
Step 4:
$\rho=\tfrac14\big(II+f^2(XX-YY+ZZ)\big)=f^2\cdot\tfrac14(II+XX-YY+ZZ)+(1-f^2)\cdot\tfrac14II$
Linearity applied to the Pauli expansion of the Bell projector; regroup the identity term.
Step 5:
$\rho=f^2\,\lvert\Phi^+\rangle\langle\Phi^+\rvert+(1-f^2)\,\frac I4$
The lemma again: the noisy pair is a mixture of the ideal pair and white noise, a Werner state with weight $w=f^2$ .
Step 6:
$\rho_{00,00}=\rho_{11,11}=\frac{f^2}2+\frac{1-f^2}4=\frac{1+f^2}4,\qquad\rho_{01,01}=\rho_{10,10}=\frac{1-f^2}4,\qquad\rho_{00,11}=\rho_{11,00}=\frac{f^2}2$
Read off the entries from $f^2$ times the Bell matrix plus $\tfrac{1-f^2}4$ on the diagonal.
Step 7:
$F=\mathrm{Tr}\big(\rho\,\lvert\Phi^+\rangle\langle\Phi^+\rvert\big)=\frac1{16}\,\mathrm{Tr}\Big[\big(II+f^2(XX-YY+ZZ)\big)\big(II+XX-YY+ZZ\big)\Big]=\frac{4+12f^2}{16}=\frac{1+3f^2}4$
$\mathrm{Tr}\big[(\sigma_a\otimes\sigma_b)(\sigma_c\otimes\sigma_d)\big]=\mathrm{Tr}(\sigma_a\sigma_c)\,\mathrm{Tr}(\sigma_b\sigma_d)=4\,\delta_{ac}\delta_{bd}$ , so only the four diagonal pairs survive: $4+f^2(4+4+4)$ ; note $(-1)(-1)=1$ for $YY$ .
Step 8:
$p=0.1:\ f=1-\tfrac{0.4}3=\tfrac{13}{15},\quad f^2=\tfrac{169}{225},\quad F=\tfrac14\Big(1+\tfrac{507}{225}\Big)=\tfrac{732}{900}=\tfrac{61}{75}\approx0.813333$
Exact arithmetic. Each qubit alone keeps the fidelity $\tfrac12(1+f)\approx0.9333$ with its own pure input, but the pair loses more, because the correlations $XX,YY,ZZ$ decay with $f^2$ .

∎

\rho=\begin{pmatrix}\frac{1+f^2}4&0&0&\frac{f^2}2\\0&\frac{1-f^2}4&0&0\\0&0&\frac{1-f^2}4&0\\\frac{f^2}2&0&0&\frac{1+f^2}4\end{pmatrix}\approx\begin{pmatrix}0.437778&0&0&0.375556\\0&0.062222&0&0\\0&0&0.062222&0\\0.375556&0&0&0.437778\end{pmatrix}

Bell pair after local depolarizing noise, p = 0.14×4

Density-matrix viewDensityMatrixView

Heatmap of the entries of ρ (populations on the diagonal, coherences off it) together with purity Tr ρ².

Tr ρ = 1.0000Purity Tr ρ² = 0.6731 (mixed)Colour: |ρ_ij| (hue = phase), full opacity at 0.4378

\rho=\begin{pmatrix} 0.4378 & 0.0000 & 0.0000 & 0.3756 \\ 0.0000 & 0.0622 & 0.0000 & 0.0000 \\ 0.0000 & 0.0000 & 0.0622 & 0.0000 \\ 0.3756 & 0.0000 & 0.0000 & 0.4378 \end{pmatrix}

\operatorname{Tr}\rho=1.0000,\qquad \gamma=\operatorname{Tr}\rho^2=\sum_{r,c}\lvert\rho_{rc}\rvert^2=0.6731

\operatorname{spec}(\rho)=\{0.0622,\ 0.0622,\ 0.0622,\ 0.8133\}

Hermitian ρ = ρ†: yes (max |ρ − ρ†| = 0.0e+0)
Positive semidefinite: yes (λ_min = 0.0622)
Unit trace: yes (Tr ρ = 1.0000)

ρ is a valid density operator (mixed state, γ < 1).

Corollary

When does the noisy Bell pair stay entangled?

\rho_w=w\,\lvert\Phi^+\rangle\langle\Phi^+\rvert+(1-w)\frac I4

Proof (7 steps)

Step 1:
$\big(\lvert\Phi^+\rangle\langle\Phi^+\rvert\big)^{T_B}=\tfrac12\sum_{i,j}\lvert i\rangle\langle j\rvert\otimes\big(\lvert i\rangle\langle j\rvert\big)^{\mathsf T}=\tfrac12\sum_{i,j}\lvert i\rangle\langle j\rvert\otimes\lvert j\rangle\langle i\rvert=\tfrac12\,\mathrm{SWAP}$
Partial transpose on the second qubit; $\sum_{ij}\lvert i\rangle\langle j\rvert\otimes\lvert j\rangle\langle i\rvert$ maps $\lvert ab\rangle$ to $\lvert ba\rangle$ .
Step 2:
$\rho_w^{T_B}=\frac w2\,\mathrm{SWAP}+\frac{1-w}4\,I,\qquad\mathrm{spec}=\Big\{\tfrac{1+w}4\ (\times3),\ \tfrac{1-3w}4\Big\}$
$\mathrm{SWAP}$ has eigenvalue $+1$ on the symmetric states $\lvert00\rangle,\lvert11\rangle,\lvert\Psi^+\rangle$ and $-1$ on the singlet $\lvert\Psi^-\rangle$ ; $I^{T_B}=I$ .
Step 3:
$\rho=\sum_kq_k\,\alpha_k\otimes\beta_k\ \Longrightarrow\ \rho^{T_B}=\sum_kq_k\,\alpha_k\otimes\beta_k^{\mathsf T}\succeq0$
Peres criterion: a separable state has a positive semidefinite partial transpose, since $\beta_k^{\mathsf T}$ is again a density matrix. Hence $w>\tfrac13$ , where $\tfrac{1-3w}4<0$ , implies entanglement.
Step 4:
$\rho_{1/3}=\tfrac14\Big(II+\tfrac13(XX-YY+ZZ)\Big)=\tfrac13\Big[\tfrac14(II+XX)+\tfrac14(II-YY)+\tfrac14(II+ZZ)\Big]$
Converse at $w=\tfrac13$ : regroup the Pauli expansion into three terms.
Step 5:
$\tfrac14(II+XX)=\tfrac12\big(\lvert{+}{+}\rangle\langle{+}{+}\rvert+\lvert{-}{-}\rangle\langle{-}{-}\rvert\big),\ \ \tfrac14(II-YY)=\tfrac12\big(\lvert{+i},{-i}\rangle\langle\cdot\rvert+\lvert{-i},{+i}\rangle\langle\cdot\rvert\big),\ \ \tfrac14(II+ZZ)=\tfrac12\big(\lvert00\rangle\langle00\rvert+\lvert11\rangle\langle11\rvert\big)$
Each term is a mixture of product states: e.g. $\lvert{\pm}{\pm}\rangle\langle{\pm}{\pm}\rvert=\tfrac14(I\pm X)\otimes(I\pm X)$ , and the cross terms $\pm(XI+IX)$ cancel in the average; likewise with the $Y$ eigenstates $\lvert{\pm i}\rangle$ , for which $\tfrac14(I+Y)\otimes(I-Y)$ and $\tfrac14(I-Y)\otimes(I+Y)$ average to $\tfrac14(II-YY)$ .
Step 6:
$w\le\tfrac13:\quad\rho_w=3w\,\rho_{1/3}+(1-3w)\,\frac I4$
A convex combination of the separable $\rho_{1/3}$ and the product state $\tfrac I2\otimes\tfrac I2$ is separable: $3w\big(\tfrac13\Phi^++\tfrac23\tfrac I4\big)+(1-3w)\tfrac I4=w\Phi^++(1-w)\tfrac I4$ .
Step 7:
$F=\frac{1+3w}4>\frac12\iff w>\frac13\iff f>\frac1{\sqrt3}\iff p<\frac34\Big(1-\frac1{\sqrt3}\Big)\approx0.316987$
Insert $w=f^2$ , $f=1-\tfrac{4p}3$ , and solve; for $p=0.1$ , $w=\tfrac{169}{225}>\tfrac13$ , so the pair is still entangled (the partial transpose has the eigenvalue $\tfrac{1-3w}4\approx-0.3133$ ).

∎

p

Noisy circuits on the platform#

\{\text{kind: channel},\ \text{channel},\ \text{targets},\ \text{params}\}

TypeScriptrunnable

import { DensityMatrix, amplitudeDamping, blochVector, runCircuit } from '@/lib/quantum-simulator';

// Circuit-IR helpers: a gate op and a noise-channel op (channels run on the density-matrix backend).
const gate = (name: string, targets: number[], params?: number[]) => ({ kind: 'gate', gate: name, targets, ...(params ? { params } : {}) });
const channel = (name: string, target: number, params: number[]) => ({ kind: 'channel', channel: name, targets: [target], params });
const run = (numQubits: number, ops: object[]) =>
  new DensityMatrix(numQubits, runCircuit({ version: 1, numQubits, numClbits: 0, ops }, { backend: 'density', seed: 0 }).densityMatrix!);
const f6 = (v: number) => (Math.abs(v) < 5e-13 ? 0 : v).toFixed(6);
const show = (rho: DensityMatrix) =>
  '[' + Array.from({ length: rho.dim }, (_, r) => '[' + Array.from({ length: rho.dim }, (_, c) => f6(rho.get(r, c).re)).join(', ') + ']').join(', ') + ']';

// (1) Same output on |+>, different channels: phase damping 0.64 and depolarizing 0.3. On |0> they differ.
console.log('H; phaseDamping(0.64)  -> ' + show(run(1, [gate('H', [0]), channel('phaseDamping', 0, [0.64])])));
console.log('H; depolarizing(0.3)   -> ' + show(run(1, [gate('H', [0]), channel('depolarizing', 0, [0.3])])));
console.log('|0>: phaseDamping(0.64) -> ' + show(run(1, [channel('phaseDamping', 0, [0.64])])) + ', depolarizing(0.3) -> ' + show(run(1, [channel('depolarizing', 0, [0.3])])));

// (2) Amplitude damping is a semigroup: idling 10 us then 20 us equals idling 30 us (T1 = 50 us).
const T1 = 50, T2 = 30;
const gammaOf = (t: number) => 1 - Math.exp(-t / T1);
const excited = () => run(1, [gate('X', [0])]);
const twoSteps = excited().applyKraus(amplitudeDamping(gammaOf(10)), [0]).applyKraus(amplitudeDamping(gammaOf(20)), [0]);
const oneStep = excited().applyKraus(amplitudeDamping(gammaOf(30)), [0]);
console.log('P(1) after 10 + 20 us = ' + f6(twoSteps.get(1, 1).re) + ', after 30 us = ' + f6(oneStep.get(1, 1).re) + ', e^{-30/50} = ' + f6(Math.exp(-30 / T1)));

// (3) thermalRelaxation(T1, T2, t) on RY(2 pi/3)|0>: the Bloch vector follows the T1 / T2 formulas exactly.
const r0 = blochVector(run(1, [gate('RY', [0], [(2 * Math.PI) / 3])]));
for (const t of [0, 10, 30, 60, 150]) {
  const r = blochVector(run(1, [gate('RY', [0], [(2 * Math.PI) / 3]), channel('thermalRelaxation', 0, [T1, T2, t])]));
  const rx = r0.x * Math.exp(-t / T2);
  const rz = 1 - (1 - r0.z) * Math.exp(-t / T1);
  console.log('t = ' + String(t).padStart(3) + ': r_x = ' + f6(r.x) + ' (' + f6(rx) + '), r_z = ' + f6(r.z) + ' (' + f6(rz) + ')');
}

// (4) Bell state with local depolarizing noise p = 0.1 on both qubits: F = (1 + 3 f^2)/4 with f = 1 - 4p/3.
const p = 0.1, f = 1 - (4 * p) / 3;
const noisyBell = run(2, [gate('H', [0]), gate('CNOT', [0, 1]), channel('depolarizing', 0, [p]), channel('depolarizing', 1, [p])]);
const F = 0.5 * (noisyBell.get(0, 0).re + noisyBell.get(0, 3).re + noisyBell.get(3, 0).re + noisyBell.get(3, 3).re);
console.log('rho = ' + show(noisyBell));
console.log('F = <Phi+|rho|Phi+> = ' + F.toFixed(12) + ', (1 + 3 f^2)/4 = ' + ((1 + 3 * f * f) / 4).toFixed(12) + ', purity = ' + f6(noisyBell.purity()));

Expected output

H; phaseDamping(0.64)  -> [[0.500000, 0.300000], [0.300000, 0.500000]]
H; depolarizing(0.3)   -> [[0.500000, 0.300000], [0.300000, 0.500000]]
|0>: phaseDamping(0.64) -> [[1.000000, 0.000000], [0.000000, 0.000000]], depolarizing(0.3) -> [[0.800000, 0.000000], [0.000000, 0.200000]]
P(1) after 10 + 20 us = 0.548812, after 30 us = 0.548812, e^{-30/50} = 0.548812
t =   0: r_x = 0.866025 (0.866025), r_z = -0.500000 (-0.500000)
t =  10: r_x = 0.620534 (0.620534), r_z = -0.228096 (-0.228096)
t =  30: r_x = 0.318593 (0.318593), r_z = 0.176783 (0.176783)
t =  60: r_x = 0.117204 (0.117204), r_z = 0.548209 (0.548209)
t = 150: r_x = 0.005835 (0.005835), r_z = 0.925319 (0.925319)
rho = [[0.437778, 0.000000, 0.000000, 0.375556], [0.000000, 0.062222, 0.000000, 0.000000], [0.000000, 0.000000, 0.062222, 0.000000], [0.375556, 0.000000, 0.000000, 0.437778]]
F = <Phi+|rho|Phi+> = 0.813333333333, (1 + 3 f^2)/4 = 0.813333333333, purity = 0.673126

0.0058

Pythonrunnable

import numpy as np

I2 = np.eye(2, dtype=complex)
X = np.array([[0, 1], [1, 0]], dtype=complex)
Y = np.array([[0, -1j], [1j, 0]])
Z = np.diag([1.0 + 0j, -1.0])

def apply(kraus, rho):
    return sum(K @ rho @ K.conj().T for K in kraus)

def complete(kraus):
    return np.allclose(sum(K.conj().T @ K for K in kraus), np.eye(kraus[0].shape[0]), atol=1e-12)

def amp_damp(g):
    return [np.diag([1, np.sqrt(1 - g)]).astype(complex), np.array([[0, np.sqrt(g)], [0, 0]], dtype=complex)]

def phase_damp(lam):
    return [np.diag([1, np.sqrt(1 - lam)]).astype(complex), np.diag([0, np.sqrt(lam)]).astype(complex)]

def depolarizing(p):
    return [np.sqrt(1 - p) * I2] + [np.sqrt(p / 3) * P for P in (X, Y, Z)]

def show(m):
    m = np.where(np.abs(m) < 5e-13, 0, m)
    cell = lambda z: f"{z.real:.6f}" if abs(z.imag) < 5e-13 else f"{z.real:.6f}{z.imag:+.6f}j"
    return "[[" + "], [".join(", ".join(cell(z) for z in row) for row in m) + "]]"

# (1) Stinespring -> Kraus: K_k = <k|_E U |0>_E for the amplitude-damping unitary on (system, environment).
g = 0.36
c, s = np.sqrt(1 - g), np.sqrt(g)
U = np.array([[1, 0, 0, 0], [0, c, s, 0], [0, -s, c, 0], [0, 0, 0, 1]], dtype=complex)  # basis |system env>
print("U unitary:", np.allclose(U.conj().T @ U, np.eye(4)))
K = [U.reshape(2, 2, 2, 2)[:, k, :, 0] for k in range(2)]  # K_k[i, j] = <i k| U |j 0>
print("K_0 =", show(K[0]), " K_1 =", show(K[1]), " complete:", complete(K))
rho = np.array([[0.3, 0.2 - 0.1j], [0.2 + 0.1j, 0.7]])
dilated = U @ np.kron(rho, np.diag([1, 0])) @ U.conj().T
reduced = np.einsum("iaja->ij", dilated.reshape(2, 2, 2, 2))  # partial trace over the environment
print("Tr_E[U(rho x |0><0|)U+] =", show(reduced), " = sum_k K rho K+:", np.allclose(reduced, apply(K, rho)))

# (2)-(4) T1 / T2 relaxation of the Bloch vector of RY(2 pi / 3)|0>, r(0) = (sin 2pi/3, 0, cos 2pi/3).
T1, T2 = 50.0, 30.0
Tphi = 1 / (1 / T2 - 1 / (2 * T1))
print(f"T1 = {T1}, T2 = {T2}: T_phi = {Tphi:.6f}, 2*T1 = {2 * T1}")
th = 2 * np.pi / 3
psi = np.array([np.cos(th / 2), np.sin(th / 2)], dtype=complex)
rho0 = np.outer(psi, psi.conj())
bloch = lambda r: np.real([np.trace(r @ P) for P in (X, Y, Z)])
rx0, rz0 = bloch(rho0)[0], bloch(rho0)[2]
print("   t    r_x       r_z     | r_x(0)e^{-t/T2}  1-(1-r_z(0))e^{-t/T1}")
for t in (0, 10, 30, 60, 150):
    gamma, lam = 1 - np.exp(-t / T1), 1 - np.exp(-2 * t / Tphi)
    r = bloch(apply(phase_damp(lam), apply(amp_damp(gamma), rho0)))
    print(f"{t:4d}  {r[0]:+.6f} {r[2]:+.6f} |  {rx0 * np.exp(-t / T2):+.6f}         {1 - (1 - rz0) * np.exp(-t / T1):+.6f}")

# Semigroup: idle periods of 10 us and 20 us compose to one of 30 us, since 1 - gamma = (1 - g1)(1 - g2).
g1, g2, g3 = (1 - np.exp(-t / T1) for t in (10, 20, 30))
r1 = apply(amp_damp(g2), apply(amp_damp(g1), rho0))
print(f"AD(g2) AD(g1) = AD(g): {np.allclose(r1, apply(amp_damp(g3), rho0))}, (1-g1)(1-g2) = {(1 - g1) * (1 - g2):.6f} = 1-g = {1 - g3:.6f}")

# (5) Depolarizing: sum_P P A P = 2 Tr(A) I for any 2x2 A, and the Bloch vector shrinks by 1 - 4p/3.
A = np.array([[1 + 2j, -0.5], [3j, 4 - 1j]])
print("sum_P P A P = 2 Tr(A) I:", np.allclose(sum(P @ A @ P for P in (I2, X, Y, Z)), 2 * np.trace(A) * I2))
for p in (0.3, 0.75, 1.0):
    out = apply(depolarizing(p), rho)
    print(f"p = {p:.2f}: r_out = {1 - 4 * p / 3:+.6f} * r_in:", np.allclose(bloch(out), (1 - 4 * p / 3) * bloch(rho)))

# Same output on |+>, different channels: phase damping 0.64 and depolarizing 0.3; they differ on |0>.
plus = np.full((2, 2), 0.5, dtype=complex)
zero = np.diag([1, 0]).astype(complex)
print("|+>: PD(0.64) ->", show(apply(phase_damp(0.64), plus)), " DEP(0.3) ->", show(apply(depolarizing(0.3), plus)))
print("|0>: PD(0.64) ->", show(apply(phase_damp(0.64), zero)), " DEP(0.3) ->", show(apply(depolarizing(0.3), zero)))

# (6) Bell state |Phi+> with local depolarizing p = 0.1 on both qubits.
p = 0.1
f = 1 - 4 * p / 3
phi = np.array([1, 0, 0, 1], dtype=complex) / np.sqrt(2)
bell = np.outer(phi, phi.conj())
print("|Phi+><Phi+| = (II + XX - YY + ZZ)/4:", np.allclose(bell, (np.eye(4) + np.kron(X, X) - np.kron(Y, Y) + np.kron(Z, Z)) / 4))
noisy = apply([np.kron(Ka, Kb) for Ka in depolarizing(p) for Kb in depolarizing(p)], bell)
F = np.real(phi.conj() @ noisy @ phi)
print(f"f = {f:.10f}, F = {F:.10f}, (1 + 3 f^2)/4 = {(1 + 3 * f ** 2) / 4:.10f}")
print("rho =", show(noisy))
print(f"(1 + f^2)/4 = {(1 + f ** 2) / 4:.6f}, (1 - f^2)/4 = {(1 - f ** 2) / 4:.6f}, f^2/2 = {f ** 2 / 2:.6f}")
pt = noisy.reshape(2, 2, 2, 2).transpose(0, 3, 2, 1).reshape(4, 4)  # partial transpose on qubit 1
print(f"min eigenvalue of the partial transpose = {np.linalg.eigvalsh(pt).min():.10f} = (1 - 3 f^2)/4 = {(1 - 3 * f ** 2) / 4:.10f}")

Expected output

U unitary: True
K_0 = [[1.000000, 0.000000], [0.000000, 0.800000]]  K_1 = [[0.000000, 0.600000], [0.000000, 0.000000]]  complete: True
Tr_E[U(rho x |0><0|)U+] = [[0.552000, 0.160000-0.080000j], [0.160000+0.080000j, 0.448000]]  = sum_k K rho K+: True
T1 = 50.0, T2 = 30.0: T_phi = 42.857143, 2*T1 = 100.0
   t    r_x       r_z     | r_x(0)e^{-t/T2}  1-(1-r_z(0))e^{-t/T1}
   0  +0.866025 -0.500000 |  +0.866025         -0.500000
  10  +0.620534 -0.228096 |  +0.620534         -0.228096
  30  +0.318593 +0.176783 |  +0.318593         +0.176783
  60  +0.117204 +0.548209 |  +0.117204         +0.548209
 150  +0.005835 +0.925319 |  +0.005835         +0.925319
AD(g2) AD(g1) = AD(g): True, (1-g1)(1-g2) = 0.548812 = 1-g = 0.548812
sum_P P A P = 2 Tr(A) I: True
p = 0.30: r_out = +0.600000 * r_in: True
p = 0.75: r_out = +0.000000 * r_in: True
p = 1.00: r_out = -0.333333 * r_in: True
|+>: PD(0.64) -> [[0.500000, 0.300000], [0.300000, 0.500000]]  DEP(0.3) -> [[0.500000, 0.300000], [0.300000, 0.500000]]
|0>: PD(0.64) -> [[1.000000, 0.000000], [0.000000, 0.000000]]  DEP(0.3) -> [[0.800000, 0.000000], [0.000000, 0.200000]]
|Phi+><Phi+| = (II + XX - YY + ZZ)/4: True
f = 0.8666666667, F = 0.8133333333, (1 + 3 f^2)/4 = 0.8133333333
rho = [[0.437778, 0.000000, 0.000000, 0.375556], [0.000000, 0.062222, 0.000000, 0.000000], [0.000000, 0.000000, 0.062222, 0.000000], [0.375556, 0.000000, 0.000000, 0.437778]]
(1 + f^2)/4 = 0.437778, (1 - f^2)/4 = 0.062222, f^2/2 = 0.375556
min eigenvalue of the partial transpose = -0.3133333333 = (1 - 3 f^2)/4 = -0.3133333333

10\,\mu\mathrm s

References

M. A. Nielsen and I. L. Chuang, Quantum Computation and Quantum Information, 10th anniversary edition, Cambridge University Press (2010), chapter 8: Quantum noise and quantum operations.doi:10.1017/CBO9780511976667
P. Krantz, M. Kjaergaard, F. Yan, T. P. Orlando, S. Gustavsson and W. D. Oliver, A quantum engineer's guide to superconducting qubits, Applied Physics Reviews 6, 021318 (2019).doi:10.1063/1.5089550
W. F. Stinespring, Positive functions on C*-algebras, Proceedings of the American Mathematical Society 6, 211–216 (1955).doi:10.1090/S0002-9939-1955-0069403-4
M.-D. Choi, Completely positive linear maps on complex matrices, Linear Algebra and its Applications 10, 285–290 (1975).doi:10.1016/0024-3795(75)90075-0
A. Peres, Separability criterion for density matrices, Physical Review Letters 77, 1413–1415 (1996).doi:10.1103/PhysRevLett.77.1413
R. F. Werner, Quantum states with Einstein-Podolsky-Rosen correlations admitting a hidden-variable model, Physical Review A 40, 4277–4281 (1989).doi:10.1103/PhysRevA.40.4277
PennyLane documentation: qml.DepolarizingChannel (Kraus operators √(1−p) I, √(p/3) X, √(p/3) Y, √(p/3) Z).docs.pennylane.ai
PennyLane documentation: qml.AmplitudeDamping, qml.PhaseDamping and qml.ThermalRelaxationError.docs.pennylane.ai
Qiskit Aer documentation: noise.depolarizing_error (E(ρ) = (1 − λ)ρ + λ I/2ⁿ) and noise.thermal_relaxation_error.qiskit.github.io

Noise, Decoherence, $T_1$ and $T_2$

Open systems and quantum channels#

Quantum channel (CPTP map) and Kraus representation

From a unitary on system and environment to Kraus operators

Every Kraus map is completely positive and trace preserving

Amplitude damping and energy relaxation#

Kraus operators and action of amplitude damping

Energy relaxation: γ(t) = 1 − e^{−t/T₁} and the semigroup law

Pure dephasing and phase flips#

Kraus operators and action of phase damping

Phase flip = phase damping (unitary freedom of Kraus operators)

Bloch-vector relaxation: T₁, T₂ and T_φ#

Bloch-vector dynamics under amplitude and phase damping

Amplitude and phase damping commute, and T₂ ≤ 2T₁

The thermal-relaxation channel of the simulator

The depolarizing channel#

Pauli twirl: Σ_P PAP = 2 Tr(A) I

The depolarizing channel shrinks the Bloch vector by 1 − 4p/3

Noisy entanglement: the Bell-state fidelity#

Pauli expansion of the Bell projector

Bell pair under local depolarizing noise

When does the noisy Bell pair stay entangled?

Noisy circuits on the platform#

References

Hands-on labs

Amplitude damping, T1 decay and engineered dephasing