docs/source/whitepaper/5-density-matrix_cn.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "f08fed11",
   "metadata": {},
   "source": [
    "# 密度矩阵和混态演化"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8ca168b6",
   "metadata": {},
   "source": [
    "## 概述\n",
    "TensorCircuit 提供了两种含噪声、混态量子演化的方法。\n",
    "$n$ 量子比特的全密度矩阵模拟是通过使用 ``tc.DMCircuit(n)`` 提供的，然后将量子操作——包括幺正门以及由 Kraus 算子指定的一般量子操作——添加到电路中。\n",
    "相对于通过 ``tc.Circuit`` 对 $n$ 个量子比特进行纯态模拟，全密度矩阵模拟会占据两倍内存，因此可模拟的最大系统大小将是纯态情况下可以模拟的一半。\n",
    "内存需求较小的选项是使用标准的 ``tc.Circuit(n)`` 对象并通过蒙特卡罗轨迹方法随机模拟开放系统演化。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "36cfaa32",
   "metadata": {},
   "source": [
    "## 设置"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "a7c8f51c",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import tensorcircuit as tc\n",
    "\n",
    "K = tc.set_backend(\"tensorflow\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0afa9c53",
   "metadata": {},
   "source": [
    "## 使用 ``tc.DMCircuit`` 进行密度矩阵模拟"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "68872b6a",
   "metadata": {},
   "source": [
    "我们在下面通过考虑单个量子比特上的简单电路来说明这种方法，该电路将对应于 $\\vert{0}\\rangle$ 状态和最大混合状态的概率混合的混合状态作为输入\n",
    "$\\rho(\\alpha) = \\alpha\\vert 0\\rangle \\langle 0\\vert + (1-\\alpha)I/2。$\n",
    "\n",
    "然后这个状态通过一个应用 $X$ 门的电路，然后是对应于带有参数 $\\gamma$ 的振幅阻尼通道 $\\mathcal{E}_\\gamma$ 的量子操作。\n",
    "这有 Kraus 运算符$K_0 = \\begin{pmatrix}\n",
    "1 & 0 \\\\ 0 & \\sqrt{1-\\gamma}\n",
    "\\end{pmatrix}, \\quad K_1 = \\begin{pmatrix}\n",
    "0 & \\sqrt{\\gamma} \\\\ 0 & 0\n",
    "\\end{pmatrix}$\n",
    "因此，该电路导致演化\n",
    "$\\rho(\\alpha) \\xrightarrow[]{X} X\\rho(\\alpha)X\\xrightarrow[]{\\mathcal{E}_\\gamma}\\sum_{i=0}^1 K_i X\\rho(\\alpha)X K_i^\\dagger$\n",
    "\n",
    "\n",
    "为了在 TensorCircuit 中模拟这一点，我们首先创建一个 ``tc.DMCircuit``（密度矩阵电路）对象并使用 ``dminputs`` 可选参数设置输入状态\n",
    "（请注意，如果将纯状态输入提供给 `` tc.DMCircuit ，这应该通过 ``inputs`` 可选参数来完成)。\n",
    "\n",
    "$\\rho(\\alpha)$ 有矩阵形式\n",
    "$\\rho(\\alpha) = \\begin{pmatrix}\n",
    "\\frac{1+\\alpha}{2} & \\\\ & \\frac{1-\\alpha}{2}\n",
    "\\end{pmatrix},$\n",
    "因此（取 $\\alpha=0.6$）我们如下初始化密度矩阵电路。\n",
    "\n",
    "为了实现诸如振幅阻尼通道之类的通用量子操作，我们使用了 ``general_kraus``，并提供了相应的 Kraus 运算符列表。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "332fd15c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def rho(alpha):\n",
    "    return np.array([[(1 + alpha) / 2, 0], [0, (1 - alpha) / 2]])\n",
    "\n",
    "\n",
    "input_state = rho(0.6)\n",
    "dmc = tc.DMCircuit(1, dminputs=input_state)\n",
    "\n",
    "dmc.x(0)\n",
    "\n",
    "\n",
    "def amp_damp_kraus(gamma):\n",
    "    K0 = np.array([[1, 0], [0, np.sqrt(1 - gamma)]])\n",
    "    K1 = np.array([[0, np.sqrt(gamma)], [0, 0]])\n",
    "    return K0, K1\n",
    "\n",
    "\n",
    "K0, K1 = amp_damp_kraus(0.3)\n",
    "dmc.general_kraus([K0, K1], 0)  # 将具有 Kraus 算子 [K0,K1] 的通道应用于 qubit 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "17aa5310",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor: shape=(2, 2), dtype=complex64, numpy=\n",
       "array([[0.44+0.j, 0.  +0.j],\n",
       "       [0.  +0.j, 0.56+0.j]], dtype=complex64)>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 得到输出密度矩阵\n",
    "dmc.state()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "e37cf64d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tf.Tensor((-0.11999999+0j), shape=(), dtype=complex64) (<tf.Tensor: shape=(1,), dtype=float32, numpy=array([1.], dtype=float32)>, -1.0)\n"
     ]
    }
   ],
   "source": [
    "# 将期望作为电路对象评估\n",
    "print(dmc.expectation_ps(z=[0]), dmc.measure(0))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6fa9e7c6",
   "metadata": {},
   "source": [
    "在上面的例子中，我们手动输入振幅阻尼通道的 Kraus 算子，以说明实现常见量子通道的一般方法。\n",
    "事实上，TensorCircuit 包含用于返回许多公共通道的 Kraus 算子的内置方法，包括振幅阻尼、去极化、相位阻尼和复位通道。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "9205dac6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "amplitudedampingchannel\n",
      "depolarizingchannel\n",
      "phasedampingchannel\n",
      "resetchannel\n"
     ]
    }
   ],
   "source": [
    "# 一组内置量子通道\n",
    "\n",
    "for k in dir(tc.channels):\n",
    "    if k.endswith(\"channel\"):\n",
    "        print(k)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "58c04997",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor: shape=(4, 4), dtype=complex64, numpy=\n",
       "array([[0.49999997+0.j, 0.        +0.j, 0.4472136 +0.j, 0.        +0.j],\n",
       "       [0.        +0.j, 0.        +0.j, 0.        +0.j, 0.        +0.j],\n",
       "       [0.4472136 +0.j, 0.        +0.j, 0.49999994+0.j, 0.        +0.j],\n",
       "       [0.        +0.j, 0.        +0.j, 0.        +0.j, 0.        +0.j]],\n",
       "      dtype=complex64)>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dmc = tc.DMCircuit(2)\n",
    "dmc.h(0)\n",
    "gamma = 0.2\n",
    "K0, K1 = tc.channels.phasedampingchannel(gamma)\n",
    "dmc.general_kraus([K0, K1], 0)\n",
    "dmc.state()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "251224a1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor: shape=(4, 4), dtype=complex64, numpy=\n",
       "array([[0.49999997+0.j, 0.        +0.j, 0.4472136 +0.j, 0.        +0.j],\n",
       "       [0.        +0.j, 0.        +0.j, 0.        +0.j, 0.        +0.j],\n",
       "       [0.4472136 +0.j, 0.        +0.j, 0.49999994+0.j, 0.        +0.j],\n",
       "       [0.        +0.j, 0.        +0.j, 0.        +0.j, 0.        +0.j]],\n",
       "      dtype=complex64)>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 或者我们可以直接使用下面的API进行速记\n",
    "\n",
    "dmc = tc.DMCircuit(2)\n",
    "dmc.h(0)\n",
    "gamma = 0.2\n",
    "dmc.phasedamping(0, gamma=0.2)\n",
    "dmc.state()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f6cf564a",
   "metadata": {},
   "source": [
    "### 自动微分、即时编译兼容性\n",
    "\n",
    "``tc.DMCircuit`` 和 ``tc.Circuit`` 一样也兼容 ML 范式，例如 自动微分、即时编译和 vmap。请参见下面的示例。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "bdd81210",
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 3\n",
    "nbatch = 2\n",
    "\n",
    "\n",
    "def loss(params, noisep):\n",
    "    c = tc.DMCircuit(n)\n",
    "    for i in range(n):\n",
    "        c.rx(i, theta=params[i])\n",
    "    for i in range(n):\n",
    "        c.depolarizing(i, px=noisep, py=noisep, pz=noisep)\n",
    "    return K.real(K.sum([c.expectation_ps(z=[i]) for i in range(n)]))\n",
    "\n",
    "\n",
    "loss_vvg = K.jit(\n",
    "    K.vectorized_value_and_grad(loss, argnums=(0, 1), vectorized_argnums=(0))\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "db09f435",
   "metadata": {},
   "outputs": [],
   "source": [
    "vs, (gparams, gnoisep) = loss_vvg(0.1 * K.ones([nbatch, n]), 0.1 * K.ones([]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "c0514b0d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(TensorShape([2]), TensorShape([2, 3]), TensorShape([]))"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "vs.shape, gparams.shape, gnoisep.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "88a993f3",
   "metadata": {},
   "source": [
    "注意噪声参数也可以被微分和即时编译！"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d37b6c4f",
   "metadata": {},
   "source": [
    "## 使用 ``tc.Circuit`` 进行蒙特卡罗噪声模拟"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fb05308d",
   "metadata": {},
   "source": [
    "对于纯态输入，蒙特卡洛方法可用于使用 ``tc.Circuit`` 而不是 ``tc.DMCircuit`` 对嘈杂的量子演化进行采样，其中混合状态是用纯状态的集合有效模拟的。\n",
    "\n",
    "至于密度矩阵模拟，可以通过提供相关的 Kraus 算子 $\\{K_i\\}$ 的列表将量子通道 $\\mathcal{E}$ 添加到电路对象中。API 与全密度矩阵模拟相同。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "772e634c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor: shape=(2,), dtype=complex64, numpy=array([0.+0.j, 1.+0.j], dtype=complex64)>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "input_state = np.array([1, 1] / np.sqrt(2))\n",
    "c = tc.Circuit(1, inputs=input_state)\n",
    "c.general_kraus(tc.channels.phasedampingchannel(0.5), 0)\n",
    "c.state()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9de97bcc",
   "metadata": {},
   "source": [
    "不过，在这个框架中，作用于 $\\vert{\\psi}\\rangle$ 的通道的输出，即\n",
    "$\\mathcal{E} ( \\vert{\\psi}\\rangle\\langle{\\psi}\\vert) = \\sum_i K_i \\vert{\\psi}\\rangle\\langle{\\psi}\\vert K_i^ \\dagger$\n",
    "被视为状态的集合 $\\frac{K_i\\vert{\\psi}\\rangle}{\\sqrt{\\langle{\\psi}\\vert K_i^\\dagger K_i \\vert{\\psi}\\rangle}}$\n",
    "每个发生的概率为 $p_i = \\langle{\\psi}\\vert K_i^\\dagger K_i \\vert{\\psi}\\rangle$.\n",
    "因此，上面的代码随机产生在状态 $\\vert{\\psi}\\rangle=\\frac{\\vert{0}\\rangle+\\vert{1}\\rangle}{\\sqrt{2}  }$\n",
    "中初始化的单个量子比特的输出通过参数 $\\gamma=0.5$ 的相位阻尼通道。\n",
    "通过使用 ``unitary_kraus`` 而不是 ``general_kraus``，可以更有效地处理 Kraus 算子都是幺正矩阵(可相差一个常数因子)的通道的蒙特卡罗模拟。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "f414e0b9",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor: shape=(), dtype=int32, numpy=3>"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "px, py, pz = 0.1, 0.2, 0.3\n",
    "c.unitary_kraus(tc.channels.depolarizingchannel(px, py, pz), 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "84a47c16",
   "metadata": {},
   "source": [
    "请注意，上面返回的 int 张量指示在此轨迹中，在电路上应用了哪个算子。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "330b5d24",
   "metadata": {},
   "source": [
    "### 外化随机性\n",
    "\n",
    "上面的 ``general_kraus`` 和 ``unitary\\_kraus`` 示例都从各自的方法内部处理随机性生成。\n",
    "也就是说，当将 Kraus 运算符的列表 $[K_0, K_1, \\ldots, K_{m-1}]$ 提供给 ``general_kraus``或`` unitary_kraus``时，\n",
    "该方法将区间 $[0, 1]$ 到 $m$ 连续区间 $[0,1] = I_0 \\cup I_1 \\cup \\ldots I_{m-1}$ 其中$I_i$的长度等于获得结果$i$的相对概率。\n",
    "然后从方法内部生成$[0,1]$中的均匀随机变量$x$，并根据$x$所在的区间选择结果$i$。\n",
    "\n",
    "在 TensorCircuit 中，我们拥有用于随机数生成和管理的完整的后端不可知基础设施。\n",
    "但是，如果我们依赖这些方法中的随机数生成，jit、随机数和后端切换之间的相互作用通常是微妙的。\n",
    "有关详细信息，请参阅 [advance.html\\#randoms-jit-backend-agnostic-and-their-interplay](../advance.html#randoms-jit-backend-agnostic-and-their-interplay)。\n",
    "\n",
    "在某些情况下，最好先从方法外部生成随机变量，然后将生成的值传递给 ``general_kraus`` 或 ``unitary_kraus`` 这可以通过可选的 ``status`` 参数来完成："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "4ce0ccf9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tf.Tensor(1, shape=(), dtype=int32)\n",
      "tf.Tensor(3, shape=(), dtype=int32)\n"
     ]
    }
   ],
   "source": [
    "px, py, pz = 0.1, 0.2, 0.3\n",
    "x = 0.5\n",
    "print(c.unitary_kraus(tc.channels.depolarizingchannel(px, py, pz), 0, status=x))\n",
    "x = 0.8\n",
    "print(c.unitary_kraus(tc.channels.depolarizingchannel(px, py, pz), 0, status=x))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6bdf395b",
   "metadata": {},
   "source": [
    "这很有用，例如，当希望使用 ``vmap`` 批量计算蒙特卡罗模拟的多次运行时。这在下面的示例中进行了说明，其中 ``vmap`` 用于并行计算 10 次模拟运行。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "c9ccd608",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<tf.Tensor: shape=(10,), dtype=complex64, numpy=\n",
       "array([ 0.99999994+0.j,  0.99999994+0.j,  0.99999994+0.j, -0.99999994+0.j,\n",
       "        0.99999994+0.j,  0.99999994+0.j,  0.99999994+0.j,  0.99999994+0.j,\n",
       "       -0.99999994+0.j,  0.99999994+0.j], dtype=complex64)>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def f(x):\n",
    "    c = tc.Circuit(1)\n",
    "    c.h(0)\n",
    "    c.unitary_kraus(tc.channels.depolarizingchannel(0.1, 0.2, 0.3), 0, status=x)\n",
    "    return c.expectation_ps(x=[0])\n",
    "\n",
    "\n",
    "f_vmap = K.vmap(f, vectorized_argnums=0)\n",
    "X = K.implicit_randn(10)\n",
    "f_vmap(X)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}