From 3a2d10900aba939ce51bccd5b1c94327a36e6884 Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Thu, 1 Aug 2024 11:58:49 -0700 Subject: [PATCH 1/8] Add Vlasov equation --- .../bloqs/arithmetic/multiplication.ipynb | 4 +- qualtran/bloqs/arithmetic/multiplication.py | 39 +- ...nian_simulation_by_matrix_arithmetic.ipynb | 522 ++++++++++++++++++ .../block_encoding/vlasov_equation.ipynb | 232 ++++++++ .../bloqs/block_encoding/vlasov_equation.py | 276 +++++++++ .../block_encoding/vlasov_equation_test.py | 19 + 6 files changed, 1074 insertions(+), 18 deletions(-) create mode 100644 qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb create mode 100644 qualtran/bloqs/block_encoding/vlasov_equation.ipynb create mode 100644 qualtran/bloqs/block_encoding/vlasov_equation.py create mode 100644 qualtran/bloqs/block_encoding/vlasov_equation_test.py diff --git a/qualtran/bloqs/arithmetic/multiplication.ipynb b/qualtran/bloqs/arithmetic/multiplication.ipynb index 0909cd39a..beca475aa 100644 --- a/qualtran/bloqs/arithmetic/multiplication.ipynb +++ b/qualtran/bloqs/arithmetic/multiplication.ipynb @@ -858,7 +858,9 @@ }, "outputs": [], "source": [ - "invert_real_number = InvertRealNumber(bitsize=10, num_frac=7)" + "from qualtran import QFxp\n", + "\n", + "invert_real_number = InvertRealNumber(QFxp(bitsize=6, num_frac=3), QFxp(bitsize=10, num_frac=7))" ] }, { diff --git a/qualtran/bloqs/arithmetic/multiplication.py b/qualtran/bloqs/arithmetic/multiplication.py index a8f558f14..37baa6dd8 100644 --- a/qualtran/bloqs/arithmetic/multiplication.py +++ b/qualtran/bloqs/arithmetic/multiplication.py @@ -33,7 +33,7 @@ from qualtran.bloqs.arithmetic.subtraction import Subtract from qualtran.bloqs.basic_gates import CNOT, TGate, Toffoli, XGate from qualtran.bloqs.mcmt import MultiControlPauli -from qualtran.symbolics import HasLength, is_symbolic, smax, SymbolicInt +from qualtran.symbolics import ceil, HasLength, is_symbolic, log2, smax, SymbolicInt if TYPE_CHECKING: import quimb.tensor as qtn @@ -539,31 +539,34 @@ class InvertRealNumber(Bloq): [Quantum Algorithms and Circuits for Scientific Computing](https://arxiv.org/pdf/1511.08253). Section 2.1. """ - bitsize: int - num_frac: int + a_dtype: QFxp + b_dtype: QFxp def __attrs_post_init__(self): - if self.num_frac == self.bitsize: + if ( + is_symbolic(self.a_dtype.num_frac) + or is_symbolic(self.a_dtype.bitsize) + or is_symbolic(self.b_dtype.bitsize) + ): + return + if self.a_dtype.num_frac == self.a_dtype.bitsize: raise ValueError("num_frac must be < bitsize since a >= 1.") + if self.b_dtype.bitsize < self.a_dtype.bitsize: + raise ValueError("b_dtype too small") @property def signature(self): - return Signature( - [ - Register("a", QFxp(self.bitsize, self.num_frac)), - Register("result", QFxp(self.bitsize, self.num_frac)), - ] - ) + return Signature([Register("a", self.a_dtype), Register("result", self.b_dtype)]) def pretty_name(self) -> str: return "1/a" def build_call_graph(self, ssa: 'SympySymbolAllocator') -> Set['BloqCountT']: # initial approximation: Figure 4 - num_int = self.bitsize - self.num_frac + num_int = self.b_dtype.bitsize - self.b_dtype.num_frac # Newton-Raphson: Eq. (1) # x' = -a * x^2 + 2 * x - num_iters = int(np.ceil(np.log2(self.bitsize))) + num_iters = ceil(log2(self.b_dtype.bitsize)) # TODO: When decomposing we will potentially need to use larger registers. # Related issue: https://github.com/quantumlib/Qualtran/issues/655 return { @@ -571,16 +574,18 @@ def build_call_graph(self, ssa: 'SympySymbolAllocator') -> Set['BloqCountT']: (CNOT(), 2 + num_int - 1), (MultiControlPauli(cvs=HasLength(num_int), target_gate=cirq.X), 1), (XGate(), 1), - (SquareRealNumber(self.bitsize), num_iters), # x^2 - (MultiplyTwoReals(self.bitsize), num_iters), # a * x^2 - (ScaleIntByReal(self.bitsize, 2), num_iters), # 2 * x - (Subtract(QUInt(self.bitsize)), num_iters), # 2 * x - a * x^2 + (SquareRealNumber(self.b_dtype.bitsize), num_iters), # x^2 + (MultiplyTwoReals(self.b_dtype.bitsize), num_iters), # a * x^2 + (ScaleIntByReal(self.b_dtype.bitsize, 2), num_iters), # 2 * x + (Subtract(QUInt(self.b_dtype.bitsize)), num_iters), # 2 * x - a * x^2 } @bloq_example def _invert_real_number() -> InvertRealNumber: - invert_real_number = InvertRealNumber(bitsize=10, num_frac=7) + from qualtran import QFxp + + invert_real_number = InvertRealNumber(QFxp(bitsize=6, num_frac=3), QFxp(bitsize=10, num_frac=7)) return invert_real_number diff --git a/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb b/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb new file mode 100644 index 000000000..918eb1364 --- /dev/null +++ b/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb @@ -0,0 +1,522 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Hamiltonian Simulation by Matrix Arithmetic" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In this notebook, we will demonstrate how we can use Qualtran's block encoding library to encode and transform data as part of a complex quantum algorithm. \n", + "\n", + "Specifically, we will study [A quantum algorithm for the linear Vlasov equation with collions](https://arxiv.org/abs/2303.03450) by Ameri et al. (2023), which simulates the dynamics of plasmas." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Note on Speedup\n", + "\n", + "As stated by the authors, this \"quantum algorithm ... can yield a quadratic speedup compared to classical algorithms that solve the same system of equations. An exponential speedup, however, does not seem possible with currently known methods due to the large norm of the matrices involved.\"\n", + "\n", + "As such, one goal of this notebook will be to perform a concrete resource estimate of the quantum algorithm to help us ascertain its speedup in practice." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import functools\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "from scipy.special import jv\n", + "\n", + "from qualtran.bloqs.basic_gates import Identity\n", + "from qualtran.bloqs.block_encoding import (\n", + " BlockEncoding,\n", + " ScaledChebyshevPolynomial,\n", + " LinearCombination,\n", + " Phase,\n", + " Product,\n", + " TensorProduct,\n", + " Unitary,\n", + ")\n", + "from qualtran.bloqs.block_encoding.sparse_matrix import (\n", + " SymmetricBandedRowColumnOracle,\n", + " ExplicitEntryOracle,\n", + " SparseMatrix,\n", + ")\n", + "from qualtran.bloqs.block_encoding.rewrites import collect_like_terms\n", + "from qualtran.bloqs.block_encoding.vlasov_equation import VlasovEntryOracle\n", + "from qualtran.resource_counting.generalizers import ignore_split_join\n", + "from qualtran.resource_counting.t_counts_from_sigma import t_counts_from_sigma" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Mathematical Formulation\n", + "\n", + "The algorithm solves the collisional Vlasov-Poisson equation (Eq. 1 and 2 in the reference). We focus on the core of the algorithm, which solves the following ODE:\n", + "\n", + "$$ \\frac{d\\mathbf{g}_k}{dt} = A\\mathbf{g}_k $$\n", + "\n", + "with solution\n", + "\n", + "$$ \\mathbf{g}_k(t) = e^{At} \\mathbf{g}_k(t = 0) $$\n", + "\n", + "where $\\mathbf{g}_k(0)$ is the initial condition and $A$ is a 3-sparse matrix:\n", + "\n", + "$$ A = -ikH + \\nu\\ \\mathrm{diag}([0, 0, 2, \\ldots, M - 1, M]) $$\n", + "\n", + "with\n", + "\n", + "$$ H = \\begin{bmatrix}\n", + "0 & \\sqrt{\\frac{1 + \\alpha}{2}} \\\\\n", + "\\sqrt{\\frac{1 + \\alpha}{2}} & 0 & 1 \\\\\n", + " & 1 & 0 & \\sqrt{\\frac{3}{2}} \\\\\n", + " & & \\ddots & \\ddots & \\ddots \\\\\n", + "& & & \\sqrt{\\frac{M - 1}{2}} & 0 & \\sqrt{\\frac{M}{2}} \\\\\n", + "& & & & \\sqrt{\\frac{M}{2}} & 0\n", + "\\end{bmatrix} $$\n", + "\n", + "where $k, \\nu, \\alpha$ are parameters of the physical system and $n$ is the number of qubits in the system, with the system size $M = 2^n - 1$.\n", + "\n", + "For simplicity, we focus on the case where $\\nu = 0, k = 1$, in which case the solution is:\n", + "\n", + "$$ \\mathbf{g}_k(t) = e^{-iHt} \\mathbf{g}_k(t = 0) $$\n", + "\n", + "which takes the familiar form of Hamiltonian simulation." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Implementation\n", + "\n", + "In principle, there are numerous methods to realize the unitary operator $e^{-iHt}$ in this problem. Our approach will be to decompose $e^{-iHt} = \\cos(Ht) + i\\sin(Ht)$, and then use the Jacobi-Anger decomposition:\n", + "\n", + "$$\\cos(Xt) = J_0(t) + 2 \\sum_{k=1}^\\infty (-1)^k J_{2k}(t)T_{2k}(X)$$\n", + "\n", + "$$\\sin(Xt) = 2 \\sum_{k=0}^\\infty (-1)^k J_{2k+1}(t)T_{2k+1}(X)$$\n", + "\n", + "where $J_n$ denotes the $n$-th Bessel function of the first kind and $T_n$ denotes the $n$-th Chebyshev polynomial of the first kind.\n", + "\n", + "Using Qualtran's block encoding library, we will produce a block encoding of $H$ as given by the matrix above, and complete each of the above decomposition steps via matrix arithmetic on block encodings." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Cost Model" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We will use $T$-complexity as the cost metric. However, following [Sunderhauf et al.](https://arxiv.org/abs/2302.10949), we will multiply the $T$-complexity of each block encoding $\\mathcal{B}[A/\\alpha]$ by its normalization factor $\\alpha$. This because $\\alpha$ is proportional to the number of repetitions we expect to need to successfully post-select $A$ when given the block encoding $\\mathcal{B}[A/\\alpha]$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def t_count(bloq: BlockEncoding):\n", + " _, sigma = bloq.call_graph(generalizer=ignore_split_join)\n", + " return int(t_counts_from_sigma(sigma) * bloq.alpha)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Matrix Powers\n", + "\n", + "Given a block encoding $\\mathcal{B}[X]$ of a matrix $X$, Qualtran makes it simple to compute a block encoding of its powers $\\mathcal{B}[X^n]$:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def identity(x: BlockEncoding) -> BlockEncoding:\n", + " '''Given B[A], return B[A^0]'''\n", + " return TensorProduct(tuple(Unitary(Identity()) for _ in range(x.system_bitsize)))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "@functools.cache\n", + "def power(n: int, x: BlockEncoding) -> BlockEncoding:\n", + " '''Given B[A], return B[A^n].'''\n", + " if n == 0:\n", + " return identity(x)\n", + " return Product.of(x, power(n - 1, x))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Chebyshev Polynomials\n", + "\n", + "We can next define a block encoding of the Chebyshev polynomial $\\mathcal{B}[T_n(X)]$ of a matrix. The definition of $T_n(x)$ is:\n", + "\n", + "$$ T_0(x) = 1 $$\n", + "$$ T_1(x) = x $$\n", + "$$ T_{n+1}(x) = 2x T_n(x) - T_{n-1}(x) $$\n", + "\n", + "We can translate this textbook recursive definition into code:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "@functools.cache\n", + "def chebyshev_exp(n: int, x: BlockEncoding) -> BlockEncoding:\n", + " '''Implementation of B[T_n(x)] that scales as O(2^n).'''\n", + " if n == 0:\n", + " return identity(x)\n", + " if n == 1:\n", + " return x\n", + " return LinearCombination.of_terms(\n", + " (2.0, Product.of(x, chebyshev_exp(n - 1, x))), (-1.0, chebyshev_exp(n - 2, x))\n", + " )" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Though simple and correct, the above implementation has a gate complexity that scales as $O(2^n)$ due to the two recursive calls. We can optimize it to $O(n^2)$ by precomputing the Chebyshev coefficients and multiplying each by the appropriate power of the matrix:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def chebyshev_quad(n: int, x: BlockEncoding) -> BlockEncoding:\n", + " '''Implementation of B[T_n(x)] that scales as O(n^2).'''\n", + " if n == 0:\n", + " return identity(x)\n", + " if n == 1:\n", + " return x\n", + " coeffs = np.polynomial.chebyshev.cheb2poly([0] * n + [1])\n", + " terms = ((coeff, power(i, x)) for i, coeff in enumerate(coeffs) if coeff != 0)\n", + " return LinearCombination.of_terms(*terms)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In principle, we could further improve the complexity to $O(n)$ with the qubitization technique of [Low and Chuang (2018)](https://arxiv.org/abs/1610.06546). This technique is implemented by Qualtran's `ChebyshevPolynomial` bloq.\n", + "\n", + "There is a technical wrinkle, however: the qubitization technique only works correctly when given a block encoding with $\\alpha = 1$, which will not be the case here. \n", + "\n", + "Our solution is to use Qualtran's `ScaledChebyshevPolynomial` bloq, which assembles a linear combination of Chebyshev polynomials that works correctly with $\\alpha \\neq 1$. By itself, this step increases the complexity back to $O(n^2)$. However, we will use Qualtran's `collect_like_terms` optimizer to fuse this linear combination with other terms in the decomposition of sine and cosine, essentially amortizing out this cost." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def chebyshev_fast(n: int, x: BlockEncoding) -> BlockEncoding:\n", + " '''Implementation of B[T_n(x)] that scales as O(n).'''\n", + "\n", + " return ScaledChebyshevPolynomial(x, n).linear_combination" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "chebyshev = chebyshev_fast" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Jacobi-Anger Approximation" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We can then construct the Jacobi-Anger approximation of sine and cosine:\n", + "\n", + "$$\\cos(Xt) = J_0(t) + 2 \\sum_{k=1}^\\infty (-1)^k J_{2k}(t)T_{2k}(X)$$\n", + "\n", + "$$\\sin(Xt) = 2 \\sum_{k=0}^\\infty (-1)^k J_{2k+1}(t)T_{2k+1}(X)$$\n", + "\n", + "where we truncate $\\infty$ to some $K$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def jacobi_anger_cosine(t: float, h: BlockEncoding, K: int) -> BlockEncoding:\n", + " return LinearCombination(\n", + " (identity(h),) + tuple(chebyshev(2 * k, h) for k in range(1, K)),\n", + " (jv(0, t),) + tuple(2 * (-1) ** k * jv(2 * k, t) for k in range(1, K)),\n", + " lambd_bits=9,\n", + " )" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def jacobi_anger_sine(t: float, h: BlockEncoding, K: int) -> BlockEncoding:\n", + " return LinearCombination(\n", + " tuple(chebyshev(2 * k + 1, h) for k in range(K)),\n", + " tuple(2 * (-1) ** k * jv(2 * k + 1, t) for k in range(K)),\n", + " lambd_bits=9,\n", + " )" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### System Block Encoding" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let us define the physical parameters of the system:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "wavenumber = 2\n", + "alpha = 2 / wavenumber**2\n", + "t = 0.5" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now we are ready to construct a block-encoding of $H$. \n", + "\n", + "We use an oracle to perform arithmetic (subtraction, division, square root) to compute each element $H_{ij}$ of the matrix given its index $(i, j)$. This oracle is implemented by the `VlasovEntryOracle` bloq in Qualtran.\n", + "\n", + "We also use an oracle to specify the positions of the 3 nonzero elements in each row/column $i$, which are $i-1$, $i$, and $i+1$. This oracle is implemented by `SymmetricBandedRowColumnOracle` in Qualtran. Putting it all together:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def vlasov_bloq(n, K, entry_oracle):\n", + " row_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + " col_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + " h = SparseMatrix(row_oracle, col_oracle, entry_oracle, eps=0)\n", + " cos = jacobi_anger_cosine(t * h.alpha, h, K)\n", + " isin = Phase(jacobi_anger_sine(t * h.alpha, h, K), phi=np.pi, eps=1e-11)\n", + " return collect_like_terms(LinearCombination.of_terms((1.0, cos), (1.0, isin)))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def vlasov_sparse(n, K):\n", + " return vlasov_bloq(n, K, VlasovEntryOracle(system_bitsize=n, entry_bitsize=n + 3, alpha=alpha))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "For comparison, we can also implement an alternative block encoding of $H$ from a classically pre-computed matrix, which we load via QROM:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def gen_vlasov_hamiltonian(n: int, alpha: float):\n", + " data = np.zeros((2**n, 2**n))\n", + " data[0][1] = data[1][0] = np.sqrt((1 + alpha) / 2)\n", + " for i in range(2, 2**n):\n", + " data[i - 1][i] = data[i][i - 1] = np.sqrt(i / 2)\n", + " data /= np.max(data)\n", + " return data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def vlasov_explicit(n, K):\n", + " data = gen_vlasov_hamiltonian(n, alpha)\n", + " return vlasov_bloq(n, K, ExplicitEntryOracle(system_bitsize=n, data=data, entry_bitsize=n + 3))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Cost Comparison\n", + "\n", + "We can compare the cost of using QROM to load the explicit matrix for $H$ versus computing the block encoding using quantum arithmetic, in terms of the system size $n$ and the degree $K$ to which the Jacobi-Anger decomposition is truncated." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "chebyshev = chebyshev_fast\n", + "explicit_n = list(t_count(vlasov_explicit(n, 5)) for n in range(2, 7))\n", + "sparse_n = list(t_count(vlasov_sparse(n, 5)) for n in range(2, 14))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.ylabel('Number of T gates')\n", + "plt.ylim(0, 5e6)\n", + "plt.xlabel('dimension n of Vlasov equation Hamiltonian')\n", + "plt.plot(range(2, 7), explicit_n, marker='o', label='explicit matrix from QROM')\n", + "plt.plot(range(2, 14), sparse_n, marker='o', label='sparse matrix from oracle')\n", + "plt.legend(loc=\"upper left\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The plot shows that the complexity of the block encoding of $e^{-iHt}$ scales with $O(2^n)$ when using QROM and $O(n)$ with quantum arithmetic." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We can also compare the costs of the different ways of computing the Chebyshev polynomial." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "chebyshev = chebyshev_fast\n", + "sparse_k = list(t_count(vlasov_sparse(5, K)) for K in range(2, 14))\n", + "\n", + "chebyshev = chebyshev_exp\n", + "sparse0_k = list(t_count(vlasov_sparse(5, K)) for K in range(2, 7))\n", + "\n", + "chebyshev = chebyshev_quad\n", + "sparse1_k = list(t_count(vlasov_sparse(5, K)) for K in range(2, 10))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.ylabel('Number of T gates')\n", + "plt.ylim(0, 3e7)\n", + "plt.xlabel('degree K of polynomial in Jacobi-Anger approximation')\n", + "plt.plot(range(2, 7), sparse0_k, marker='o', label='chebyshev_exp')\n", + "plt.plot(range(2, 10), sparse1_k, marker='o', label='chebyshev_quad')\n", + "plt.plot(range(2, 14), sparse_k, marker='o', label='chebyshev_fast')\n", + "plt.legend(loc=\"upper left\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Empirically, when using the efficient method to compute the Chebyshev polynomial, the overall complexity of the block encoding of $e^{-iHt}$ scales as $O(nK^2)$. The quadratic scaling in $K$ is due to the computation of the linear combination of $K$ Chebyshev polynomials, each of which has $O(K)$ amortized complexity, in the Jacobi-Anger decomposition." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "qualtran", + "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.11.2" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/qualtran/bloqs/block_encoding/vlasov_equation.ipynb b/qualtran/bloqs/block_encoding/vlasov_equation.ipynb new file mode 100644 index 000000000..f3f93dd95 --- /dev/null +++ b/qualtran/bloqs/block_encoding/vlasov_equation.ipynb @@ -0,0 +1,232 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import sympy\n", + "\n", + "from qualtran.bloqs.block_encoding.sparse_matrix import (\n", + " ExplicitEntryOracle,\n", + " SymmetricBandedRowColumnOracle,\n", + " SparseMatrix,\n", + ")\n", + "from qualtran.bloqs.block_encoding.vlasov_equation import VlasovEntryOracle\n", + "from qualtran.drawing import show_counts_sigma\n", + "from qualtran.resource_counting.generalizers import ignore_split_join\n", + "from qualtran.resource_counting.t_counts_from_sigma import t_counts_from_sigma" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "n = sympy.Symbol('n')\n", + "row_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + "col_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + "entry_oracle = VlasovEntryOracle(n, 7, alpha=0.2)\n", + "bloq = SparseMatrix(row_oracle, col_oracle, entry_oracle, eps=0)" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "#### Counts totals:\n", + " - `Allocate`: 1\n", + " - `Allocate`: 1\n", + " - `ArbitraryClifford(n=2)`: $\\displaystyle 221 n + 557$\n", + " - `CNOT`: $\\displaystyle 98 n + 223$\n", + " - `C[2][TwoBitSwap]`: 14\n", + " - `C[2][XGate]`: 7\n", + " - `C[Ry(0.0024867959858108648π)]`: 7\n", + " - `C[Ry(0.435905783151025π)]`: 1\n", + " - `C[TwoBitSwap]`: $\\displaystyle n$\n", + " - `Free`: 1\n", + " - `Free`: 1\n", + " - `H`: 12\n", + " - `Rz(-0.3918265520306073π)`: 2\n", + " - `Rz(0.3918265520306073π)`: 2\n", + " - `T`: $\\displaystyle 8716 n - 32 \\left(n - 2\\right) \\left(30 n - 30\\right) - 72 \\left(n - 1\\right)^{2} - 32 \\left(n - 1\\right) \\left(6 n - 12\\right) + 8 \\left\\lceil{480 n - 60 \\left(n - 1\\right)^{2} + 1486.0}\\right\\rceil + 34084$\n", + " - `TwoBitSwap`: $\\displaystyle n$\n", + " - `XGate`: $\\displaystyle 10 n + 165$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "g, sigma = bloq.call_graph(generalizer=ignore_split_join)\n", + "show_counts_sigma(sigma)" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "8716*n - 32*(n - 2)*(30*n - 30) - 72*(n - 1)**2 - 32*(n - 1)*(6*n - 12) + 8*ceiling(480*n - 60*(n - 1)**2 + 1486.0) + 34292\n" + ] + } + ], + "source": [ + "ts = t_counts_from_sigma(sigma)\n", + "print(ts)" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "text/latex": [ + "$\\displaystyle 25416$" + ], + "text/plain": [ + "25416" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "ts.subs(n, 11)" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "def gen_vlasov_hamiltonian(n, alpha, m):\n", + " data = np.zeros((2**n, 2**n))\n", + " data[0][1] = data[1][0] = np.sqrt((1 + alpha) / 2)\n", + " for i in range(2, m + 1):\n", + " data[i - 1][i] = data[i][i - 1] = np.sqrt(i / 2)\n", + " data /= np.max(data)\n", + " return data" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "n = 11\n", + "k = 2\n", + "alpha = 2 / k**2\n", + "data = gen_vlasov_hamiltonian(n, alpha, m=(2**n - 1))\n", + "row_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + "col_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + "entry_oracle = ExplicitEntryOracle(system_bitsize=n, data=data, entry_bitsize=7)\n", + "bloq = SparseMatrix(row_oracle, col_oracle, entry_oracle, eps=0)" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "#### Counts totals:\n", + " - `Allocate`: 1\n", + " - `Allocate`: 1\n", + " - `Allocate`: 1\n", + " - `ArbitraryClifford(n=2)`: 860046\n", + " - `CNOT`: 131200\n", + " - `C[Ry(0.015625π)]`: 1\n", + " - `C[Ry(0.03125π)]`: 1\n", + " - `C[Ry(0.0625π)]`: 1\n", + " - `C[Ry(0.125π)]`: 1\n", + " - `C[Ry(0.25π)]`: 1\n", + " - `C[Ry(0.5π)]`: 1\n", + " - `C[Ry(π)]`: 1\n", + " - `Free`: 1\n", + " - `Free`: 1\n", + " - `Free`: 1\n", + " - `H`: 12\n", + " - `Rz(-0.3918265520306073π)`: 2\n", + " - `Rz(0.3918265520306073π)`: 2\n", + " - `T`: 229368\n", + " - `TwoBitSwap`: 11\n", + " - `XGate`: 50" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "g, sigma = bloq.call_graph(generalizer=ignore_split_join)\n", + "show_counts_sigma(sigma)" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "229576\n" + ] + } + ], + "source": [ + "ts = t_counts_from_sigma(sigma)\n", + "print(ts)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "qualtran", + "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.11.2" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/qualtran/bloqs/block_encoding/vlasov_equation.py b/qualtran/bloqs/block_encoding/vlasov_equation.py new file mode 100644 index 000000000..02c2e1fba --- /dev/null +++ b/qualtran/bloqs/block_encoding/vlasov_equation.py @@ -0,0 +1,276 @@ +# Copyright 2024 Google LLC +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# https://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +from typing import cast, Dict, Set + +import cirq +import numpy as np +from attrs import frozen +from numpy.typing import NDArray + +from qualtran import ( + Bloq, + bloq_example, + BloqBuilder, + DecomposeTypeError, + QFxp, + Register, + Signature, + Soquet, + SoquetT, +) +from qualtran._infra.controlled import CtrlSpec +from qualtran.bloqs.arithmetic.comparison import GreaterThan +from qualtran.bloqs.arithmetic.multiplication import InvertRealNumber +from qualtran.bloqs.arithmetic.trigonometric import ArcSin +from qualtran.bloqs.basic_gates import Ry, Swap, Toffoli +from qualtran.bloqs.basic_gates.x_basis import XGate +from qualtran.bloqs.block_encoding.sparse_matrix import EntryOracle, SparseMatrix +from qualtran.bloqs.mcmt.and_bloq import And +from qualtran.bloqs.mcmt.multi_control_multi_target_pauli import MultiControlPauli +from qualtran.bloqs.rotations.phase_gradient import _fxp +from qualtran.resource_counting import BloqCountT, SympySymbolAllocator +from qualtran.simulation.classical_sim import ClassicalValT +from qualtran.symbolics import ceil, HasLength, is_symbolic, SymbolicInt + + +@frozen +class InverseSquareRoot(Bloq): + r"""Compute the inverse square root of a fixed-point number. + + Args: + bitsize: Number of bits used to represent the number. + num_frac: Number of fraction bits in the number. + num_iters: Number of Newton-Raphson iterations. + Registers: + x: `bitsize`-sized input register. + result: `bitsize`-sized output register. + References: + [Optimizing Quantum Circuits for Arithmetic](https://arxiv.org/abs/1805.12445). Appendix C. + """ + + bitsize: SymbolicInt + num_frac: SymbolicInt + num_iters: SymbolicInt = 4 # reference studies 3, 4, or 5 iterations + + def __attrs_post_init__(self): + if ( + not is_symbolic(self.num_frac) + and not is_symbolic(self.bitsize) + and self.num_frac > self.bitsize + ): + raise ValueError("num_frac must be < bitsize.") + + @property + def signature(self): + return Signature( + [ + Register("x", QFxp(self.bitsize, self.num_frac)), + Register("result", QFxp(self.bitsize, self.num_frac)), + ] + ) + + def pretty_name(self) -> str: + return "1/sqrt(x)" + + def on_classical_vals( + self, x: ClassicalValT, result: ClassicalValT + ) -> Dict[str, ClassicalValT]: + if is_symbolic(self.bitsize): + raise ValueError(f"Symbolic bitsize {self.bitsize} not supported") + x_fxp: float = _fxp(x / 2**self.bitsize, self.bitsize).astype(float) + result ^= int(1 / np.sqrt(x_fxp) * 2**self.bitsize) + return {'x': x, 'result': result} + + def build_call_graph(self, ssa: 'SympySymbolAllocator') -> Set['BloqCountT']: + n = self.bitsize + p = self.bitsize - self.num_frac + m = self.num_iters + # directly copied from T_invsqrt on page 10 of reference + ts = ceil( + n**2 * (15 * m / 2 + 3) + + 15 * n * p * m + + n * (23 * m / 2 + 5) + - 15 * p**2 * m + + 15 * p * m + - 2 * m + ) + return {(Toffoli(), ts)} + + +@frozen +class VlasovEntryOracle(EntryOracle): + r"""Oracle specifying the entries of the Hamiltonian for the Vlasov equation. + + Args: + system_bitsize: The number of bits used to represent an index. The value of `M` in the + referenced equation must be equal to `2**system_bitsize - 1`. + entry_bitsize: The number of bits of precision to represent the arcsin of each entry. + alpha: The physical parameter $\alpha$ in the referenced equation. + + Registers: + q: The flag qubit that is rotated proportionally to the value of the entry. + i: The row index. + j: The column index. + + References: + [A quantum algorithm for the linear Vlasov equation with collisions](https://arxiv.org/pdf/2303.03450) (2022). Eq. (27). + """ + + system_bitsize: SymbolicInt + entry_bitsize: SymbolicInt + alpha: float + + def build_call_graph(self, ssa: 'SympySymbolAllocator') -> Set['BloqCountT']: + mcx = MultiControlPauli(cvs=HasLength(self.system_bitsize), target_gate=cirq.X) + gt = GreaterThan(self.system_bitsize, self.system_bitsize) + num_frac = self.entry_bitsize - (self.system_bitsize - 1) + inv = InvertRealNumber(QFxp(self.system_bitsize, 1), QFxp(self.entry_bitsize, num_frac)) + invsqrt = InverseSquareRoot(self.entry_bitsize, num_frac) + arcsin = ArcSin(self.entry_bitsize, num_frac) + swap = Swap(1).controlled(CtrlSpec(cvs=(1, 0))) + x = XGate().controlled(CtrlSpec(cvs=(1, 0))) + ry = Ry(2**-self.entry_bitsize).controlled() + return { + (mcx, 1), + (And(), 1), + (Ry(2 * np.arccos(np.sqrt((1 + self.alpha) / 2))).controlled(), 1), + (gt, 2), + (Swap(self.system_bitsize).controlled(), 1), + (inv, 1), + (invsqrt, 1), + (arcsin, 1), + (And(cv1=0, cv2=0), 1), + (XGate(), 1), + (swap, 2 * self.entry_bitsize), + (x, self.entry_bitsize), + (ry, self.entry_bitsize), + (arcsin.adjoint(), 1), + (invsqrt.adjoint(), 1), + (inv.adjoint(), 1), + (Swap(self.system_bitsize).controlled(), 1), + (XGate(), 1), + (And(cv1=0, cv2=0, uncompute=True), 1), + (gt.adjoint(), 2), + (And().adjoint(), 1), + (mcx.adjoint(), 2), + } + + def build_composite_bloq( + self, bb: BloqBuilder, q: SoquetT, i: SoquetT, j: SoquetT + ) -> Dict[str, SoquetT]: + if is_symbolic(self.system_bitsize) or is_symbolic(self.entry_bitsize): + raise DecomposeTypeError(f"Cannot decompose symbolic {self=}") + + mcx = MultiControlPauli(cvs=(0,) * self.system_bitsize, target_gate=cirq.X) + i_bits, i_zero = bb.add_t(mcx, controls=bb.split(cast(Soquet, i)), target=bb.allocate(1)) + j_bits, j_zero = bb.add_t(mcx, controls=bb.split(cast(Soquet, j)), target=bb.allocate(1)) + i = bb.join(cast(NDArray, i_bits)) + j = bb.join(cast(NDArray, j_bits)) + i_zero_j_zero, i_or_j_zero = bb.add_t(And(), ctrl=np.array([i_zero, j_zero])) + + # case 1: i = 0 or j = 0, entry is sqrt((1 + alpha) / 2) + i_or_j_zero, q = bb.add_t( + Ry(2 * np.arccos(np.sqrt((1 + self.alpha) / 2))).controlled(), ctrl=i_or_j_zero, q=q + ) + + gt = GreaterThan(self.system_bitsize, self.system_bitsize) + i_gt_j = bb.allocate(1) + j_gt_i = bb.allocate(1) + i, j, i_gt_j = bb.add_t(gt, a=i, b=j, target=i_gt_j) + i, j, j_gt_i = bb.add_t(gt, a=j, b=i, target=j_gt_i) + + # case 2: i > j, entry is sqrt(i / 2) + # swap i and j such that we fall into case 3 + i_gt_j, i, j = bb.add_t(Swap(self.system_bitsize).controlled(), ctrl=i_gt_j, x=i, y=j) + + # case 3: i < j, entry is sqrt(j / 2) + num_frac = self.entry_bitsize - (self.system_bitsize - 1) + entry0 = bb.allocate(dtype=QFxp(self.entry_bitsize, num_frac)) + # compute entry = 1 / (j / 2) + inv = InvertRealNumber(QFxp(self.system_bitsize, 1), QFxp(self.entry_bitsize, num_frac)) + j, entry0 = bb.add_t(inv, a=j, result=entry0) + # compute entry' = 1 / sqrt(1 / (j / 2)) = sqrt(j / 2) + invsqrt = InverseSquareRoot(self.entry_bitsize, num_frac) + entry1 = bb.allocate(dtype=QFxp(self.entry_bitsize, num_frac)) + entry0, entry1 = bb.add_t(invsqrt, x=entry0, result=entry1) + # compute entry' = arcsin(entry) + arcsin = ArcSin(self.entry_bitsize, num_frac) + entry = bb.allocate(dtype=QFxp(self.entry_bitsize, num_frac)) + entry1, entry = bb.add_t(arcsin, x=entry1, result=entry) + + i_gt_j_j_gt_i, i_neq_j = bb.add_t(And(cv1=0, cv2=0), ctrl=np.array([i_gt_j, j_gt_i])) + i_gt_j, j_gt_i = cast(NDArray, i_gt_j_j_gt_i) + i_neq_j = cast(Soquet, bb.add(XGate(), q=i_neq_j)) + + # if i_neq_j and not i_or_j_zero, rotate by entry + entry_bits = bb.split(cast(Soquet, entry)) + zero = bb.allocate(1) + ctrls = np.array([i_neq_j, i_or_j_zero]) + swap = Swap(1).controlled(CtrlSpec(cvs=(1, 0))) + for k in range(len(entry_bits)): + # selectively swap entry bit with 0 so we only rotate if i_neq_j and not i_or_j_zero + ctrls, entry_bits[k], zero = bb.add_t(swap, ctrl=ctrls, x=entry_bits[k], y=zero) + # flip q because entry is arcsin, not arccos + ctrls, q = bb.add_t(XGate().controlled(CtrlSpec(cvs=(1, 0))), ctrl=ctrls, q=q) + entry_bits[k], q = bb.add_t(Ry(2**-k).controlled(), ctrl=entry_bits[k], q=q) + ctrls, entry_bits[k], zero = bb.add_t(swap, ctrl=ctrls, x=entry_bits[k], y=zero) + i_neq_j, i_or_j_zero = cast(NDArray, ctrls) + bb.free(cast(Soquet, zero)) + entry = bb.join(entry_bits) + + entry1, entry = bb.add_t(arcsin.adjoint(), x=entry1, result=entry) + bb.free(entry) + entry0, entry1 = bb.add_t(invsqrt.adjoint(), x=entry0, result=entry1) + bb.free(entry1) + j, entry0 = bb.add_t(inv.adjoint(), a=j, result=entry0) + bb.free(cast(Soquet, entry0)) + + i_gt_j, i, j = bb.add_t(Swap(self.system_bitsize).controlled(), ctrl=i_gt_j, x=i, y=j) + + i_neq_j = bb.add(XGate(), q=i_neq_j) + i_gt_j, j_gt_i = cast( + NDArray, + bb.add( + And(cv1=0, cv2=0, uncompute=True), ctrl=np.array([i_gt_j, j_gt_i]), target=i_neq_j + ), + ) + + i, j, j_gt_i = bb.add_t(gt.adjoint(), a=j, b=i, target=j_gt_i) + i, j, i_gt_j = bb.add_t(gt.adjoint(), a=i, b=j, target=i_gt_j) + bb.free(cast(Soquet, i_gt_j)) + bb.free(cast(Soquet, j_gt_i)) + + i_zero, j_zero = cast( + NDArray, bb.add(And().adjoint(), ctrl=i_zero_j_zero, target=i_or_j_zero) + ) + j_bits, j_zero = bb.add_t(mcx.adjoint(), controls=bb.split(cast(Soquet, j)), target=j_zero) + i_bits, i_zero = bb.add_t(mcx.adjoint(), controls=bb.split(cast(Soquet, i)), target=i_zero) + bb.free(cast(Soquet, j_zero)) + bb.free(cast(Soquet, i_zero)) + + return {"q": q, "i": bb.join(cast(NDArray, i_bits)), "j": bb.join(cast(NDArray, j_bits))} + + +@bloq_example +def _vlasov_block_encoding() -> SparseMatrix: + from qualtran.bloqs.block_encoding.sparse_matrix import SymmetricBandedRowColumnOracle + + row_oracle = SymmetricBandedRowColumnOracle(3, bandsize=1) + col_oracle = SymmetricBandedRowColumnOracle(3, bandsize=1) + entry_oracle = VlasovEntryOracle(3, 7, alpha=0.2) + symmetric_banded_matrix_block_encoding = SparseMatrix( + row_oracle, col_oracle, entry_oracle, eps=0 + ) + return symmetric_banded_matrix_block_encoding diff --git a/qualtran/bloqs/block_encoding/vlasov_equation_test.py b/qualtran/bloqs/block_encoding/vlasov_equation_test.py new file mode 100644 index 000000000..b5b31d453 --- /dev/null +++ b/qualtran/bloqs/block_encoding/vlasov_equation_test.py @@ -0,0 +1,19 @@ +# Copyright 2024 Google LLC +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# https://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. + +from qualtran.bloqs.block_encoding.vlasov_equation import _vlasov_block_encoding + + +def test_vlasov_block_encoding(bloq_autotester): + bloq_autotester(_vlasov_block_encoding) From 281624136be41e190438246c7ddb8ea03d6dc2be Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Thu, 1 Aug 2024 18:36:34 -0700 Subject: [PATCH 2/8] Update notebook --- ...nian_simulation_by_matrix_arithmetic.ipynb | 48 ++++++++++++++----- 1 file changed, 37 insertions(+), 11 deletions(-) diff --git a/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb b/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb index 918eb1364..5aa349e39 100644 --- a/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb +++ b/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb @@ -353,8 +353,7 @@ "outputs": [], "source": [ "def vlasov_bloq(n, K, entry_oracle):\n", - " row_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", - " col_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", + " row_oracle = col_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", " h = SparseMatrix(row_oracle, col_oracle, entry_oracle, eps=0)\n", " cos = jacobi_anger_cosine(t * h.alpha, h, K)\n", " isin = Phase(jacobi_anger_sine(t * h.alpha, h, K), phi=np.pi, eps=1e-11)\n", @@ -430,12 +429,25 @@ "metadata": {}, "outputs": [], "source": [ - "plt.ylabel('Number of T gates')\n", + "#import seaborn as sns\n", + "#sns.set_theme()\n", + "#plt.rcParams[\"font.family\"] = \"Google Sans\"" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.ylabel('Number of $T$ gates')\n", "plt.ylim(0, 5e6)\n", - "plt.xlabel('dimension n of Vlasov equation Hamiltonian')\n", + "plt.xlabel('dimension $n$ of Vlasov equation Hamiltonian')\n", "plt.plot(range(2, 7), explicit_n, marker='o', label='explicit matrix from QROM')\n", "plt.plot(range(2, 14), sparse_n, marker='o', label='sparse matrix from oracle')\n", - "plt.legend(loc=\"upper left\")" + "plt.legend(loc=\"lower right\")\n", + "plt.title('QROM vs Quantum Arithmetic for Vlasov Matrix')\n", + "plt.savefig('scaling-n.png', dpi=300)" ] }, { @@ -474,13 +486,15 @@ "metadata": {}, "outputs": [], "source": [ - "plt.ylabel('Number of T gates')\n", + "plt.ylabel('Number of $T$ gates')\n", "plt.ylim(0, 3e7)\n", - "plt.xlabel('degree K of polynomial in Jacobi-Anger approximation')\n", - "plt.plot(range(2, 7), sparse0_k, marker='o', label='chebyshev_exp')\n", - "plt.plot(range(2, 10), sparse1_k, marker='o', label='chebyshev_quad')\n", - "plt.plot(range(2, 14), sparse_k, marker='o', label='chebyshev_fast')\n", - "plt.legend(loc=\"upper left\")" + "plt.xlabel('degree $K$ of polynomial in Jacobi-Anger approximation')\n", + "plt.plot(range(2, 7), sparse0_k, marker='o', label='naive recursive')\n", + "plt.plot(range(2, 10), sparse1_k, marker='o', label='iterative')\n", + "plt.plot(range(2, 14), sparse_k, marker='o', label='qubitization')\n", + "plt.legend(loc=\"upper right\")\n", + "plt.title('Comparison of Chebyshev Implementations')\n", + "plt.savefig('chebyshev.png', dpi=300)" ] }, { @@ -490,6 +504,18 @@ "Empirically, when using the efficient method to compute the Chebyshev polynomial, the overall complexity of the block encoding of $e^{-iHt}$ scales as $O(nK^2)$. The quadratic scaling in $K$ is due to the computation of the linear combination of $K$ Chebyshev polynomials, each of which has $O(K)$ amortized complexity, in the Jacobi-Anger decomposition." ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# from qualtran.drawing import show_bloq_via_qpic\n", + "# from qualtran.bloqs.basic_gates import Hadamard\n", + "# bloq = chebyshev(5, Unitary(Hadamard())).decompose_bloq().flatten()\n", + "# show_bloq_via_qpic(bloq)" + ] + }, { "cell_type": "code", "execution_count": null, From 4e33084bc370416221b5b38479dcba1bd3418740 Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Fri, 2 Aug 2024 11:53:19 -0700 Subject: [PATCH 3/8] Fix sign --- .../hamiltonian_simulation_by_matrix_arithmetic.ipynb | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb b/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb index 5aa349e39..307015ae3 100644 --- a/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb +++ b/qualtran/bloqs/block_encoding/hamiltonian_simulation_by_matrix_arithmetic.ipynb @@ -103,7 +103,7 @@ "source": [ "## Implementation\n", "\n", - "In principle, there are numerous methods to realize the unitary operator $e^{-iHt}$ in this problem. Our approach will be to decompose $e^{-iHt} = \\cos(Ht) + i\\sin(Ht)$, and then use the Jacobi-Anger decomposition:\n", + "In principle, there are numerous methods to realize the unitary operator $e^{-iHt}$ in this problem. Our approach will be to decompose $e^{-iHt} = \\cos(Ht) - i\\sin(Ht)$, and then use the Jacobi-Anger decomposition:\n", "\n", "$$\\cos(Xt) = J_0(t) + 2 \\sum_{k=1}^\\infty (-1)^k J_{2k}(t)T_{2k}(X)$$\n", "\n", @@ -356,8 +356,8 @@ " row_oracle = col_oracle = SymmetricBandedRowColumnOracle(n, bandsize=1)\n", " h = SparseMatrix(row_oracle, col_oracle, entry_oracle, eps=0)\n", " cos = jacobi_anger_cosine(t * h.alpha, h, K)\n", - " isin = Phase(jacobi_anger_sine(t * h.alpha, h, K), phi=np.pi, eps=1e-11)\n", - " return collect_like_terms(LinearCombination.of_terms((1.0, cos), (1.0, isin)))" + " isin = Phase(jacobi_anger_sine(t * h.alpha, h, K), phi=1 / 2, eps=1e-11)\n", + " return collect_like_terms(LinearCombination.of_terms((1.0, cos), (-1.0, isin)))" ] }, { From 36baa593c0fa567002d4df2095d6d8a32f0624a9 Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Wed, 7 Aug 2024 11:13:13 -0700 Subject: [PATCH 4/8] Fix mypy --- .../bloqs/block_encoding/vlasov_equation.py | 46 +++++++++---------- 1 file changed, 23 insertions(+), 23 deletions(-) diff --git a/qualtran/bloqs/block_encoding/vlasov_equation.py b/qualtran/bloqs/block_encoding/vlasov_equation.py index 02c2e1fba..e7a09a332 100644 --- a/qualtran/bloqs/block_encoding/vlasov_equation.py +++ b/qualtran/bloqs/block_encoding/vlasov_equation.py @@ -174,43 +174,43 @@ def build_composite_bloq( raise DecomposeTypeError(f"Cannot decompose symbolic {self=}") mcx = MultiControlPauli(cvs=(0,) * self.system_bitsize, target_gate=cirq.X) - i_bits, i_zero = bb.add_t(mcx, controls=bb.split(cast(Soquet, i)), target=bb.allocate(1)) - j_bits, j_zero = bb.add_t(mcx, controls=bb.split(cast(Soquet, j)), target=bb.allocate(1)) + i_bits, i_zero = bb.add(mcx, controls=bb.split(cast(Soquet, i)), target=bb.allocate(1)) + j_bits, j_zero = bb.add(mcx, controls=bb.split(cast(Soquet, j)), target=bb.allocate(1)) i = bb.join(cast(NDArray, i_bits)) j = bb.join(cast(NDArray, j_bits)) - i_zero_j_zero, i_or_j_zero = bb.add_t(And(), ctrl=np.array([i_zero, j_zero])) + i_zero_j_zero, i_or_j_zero = bb.add(And(), ctrl=np.array([i_zero, j_zero])) # case 1: i = 0 or j = 0, entry is sqrt((1 + alpha) / 2) - i_or_j_zero, q = bb.add_t( + i_or_j_zero, q = bb.add( Ry(2 * np.arccos(np.sqrt((1 + self.alpha) / 2))).controlled(), ctrl=i_or_j_zero, q=q ) gt = GreaterThan(self.system_bitsize, self.system_bitsize) i_gt_j = bb.allocate(1) j_gt_i = bb.allocate(1) - i, j, i_gt_j = bb.add_t(gt, a=i, b=j, target=i_gt_j) - i, j, j_gt_i = bb.add_t(gt, a=j, b=i, target=j_gt_i) + i, j, i_gt_j = bb.add(gt, a=i, b=j, target=i_gt_j) + i, j, j_gt_i = bb.add(gt, a=j, b=i, target=j_gt_i) # case 2: i > j, entry is sqrt(i / 2) # swap i and j such that we fall into case 3 - i_gt_j, i, j = bb.add_t(Swap(self.system_bitsize).controlled(), ctrl=i_gt_j, x=i, y=j) + i_gt_j, i, j = bb.add(Swap(self.system_bitsize).controlled(), ctrl=i_gt_j, x=i, y=j) # case 3: i < j, entry is sqrt(j / 2) num_frac = self.entry_bitsize - (self.system_bitsize - 1) entry0 = bb.allocate(dtype=QFxp(self.entry_bitsize, num_frac)) # compute entry = 1 / (j / 2) inv = InvertRealNumber(QFxp(self.system_bitsize, 1), QFxp(self.entry_bitsize, num_frac)) - j, entry0 = bb.add_t(inv, a=j, result=entry0) + j, entry0 = bb.add(inv, a=j, result=entry0) # compute entry' = 1 / sqrt(1 / (j / 2)) = sqrt(j / 2) invsqrt = InverseSquareRoot(self.entry_bitsize, num_frac) entry1 = bb.allocate(dtype=QFxp(self.entry_bitsize, num_frac)) - entry0, entry1 = bb.add_t(invsqrt, x=entry0, result=entry1) + entry0, entry1 = bb.add(invsqrt, x=entry0, result=entry1) # compute entry' = arcsin(entry) arcsin = ArcSin(self.entry_bitsize, num_frac) entry = bb.allocate(dtype=QFxp(self.entry_bitsize, num_frac)) - entry1, entry = bb.add_t(arcsin, x=entry1, result=entry) + entry1, entry = bb.add(arcsin, x=entry1, result=entry) - i_gt_j_j_gt_i, i_neq_j = bb.add_t(And(cv1=0, cv2=0), ctrl=np.array([i_gt_j, j_gt_i])) + i_gt_j_j_gt_i, i_neq_j = bb.add(And(cv1=0, cv2=0), ctrl=np.array([i_gt_j, j_gt_i])) i_gt_j, j_gt_i = cast(NDArray, i_gt_j_j_gt_i) i_neq_j = cast(Soquet, bb.add(XGate(), q=i_neq_j)) @@ -221,23 +221,23 @@ def build_composite_bloq( swap = Swap(1).controlled(CtrlSpec(cvs=(1, 0))) for k in range(len(entry_bits)): # selectively swap entry bit with 0 so we only rotate if i_neq_j and not i_or_j_zero - ctrls, entry_bits[k], zero = bb.add_t(swap, ctrl=ctrls, x=entry_bits[k], y=zero) + ctrls, entry_bits[k], zero = bb.add(swap, ctrl=ctrls, x=entry_bits[k], y=zero) # flip q because entry is arcsin, not arccos - ctrls, q = bb.add_t(XGate().controlled(CtrlSpec(cvs=(1, 0))), ctrl=ctrls, q=q) - entry_bits[k], q = bb.add_t(Ry(2**-k).controlled(), ctrl=entry_bits[k], q=q) - ctrls, entry_bits[k], zero = bb.add_t(swap, ctrl=ctrls, x=entry_bits[k], y=zero) + ctrls, q = bb.add(XGate().controlled(CtrlSpec(cvs=(1, 0))), ctrl=ctrls, q=q) + entry_bits[k], q = bb.add(Ry(2**-k).controlled(), ctrl=entry_bits[k], q=q) + ctrls, entry_bits[k], zero = bb.add(swap, ctrl=ctrls, x=entry_bits[k], y=zero) i_neq_j, i_or_j_zero = cast(NDArray, ctrls) bb.free(cast(Soquet, zero)) entry = bb.join(entry_bits) - entry1, entry = bb.add_t(arcsin.adjoint(), x=entry1, result=entry) + entry1, entry = bb.add(arcsin.adjoint(), x=entry1, result=entry) bb.free(entry) - entry0, entry1 = bb.add_t(invsqrt.adjoint(), x=entry0, result=entry1) + entry0, entry1 = bb.add(invsqrt.adjoint(), x=entry0, result=entry1) bb.free(entry1) - j, entry0 = bb.add_t(inv.adjoint(), a=j, result=entry0) + j, entry0 = bb.add(inv.adjoint(), a=j, result=entry0) bb.free(cast(Soquet, entry0)) - i_gt_j, i, j = bb.add_t(Swap(self.system_bitsize).controlled(), ctrl=i_gt_j, x=i, y=j) + i_gt_j, i, j = bb.add(Swap(self.system_bitsize).controlled(), ctrl=i_gt_j, x=i, y=j) i_neq_j = bb.add(XGate(), q=i_neq_j) i_gt_j, j_gt_i = cast( @@ -247,16 +247,16 @@ def build_composite_bloq( ), ) - i, j, j_gt_i = bb.add_t(gt.adjoint(), a=j, b=i, target=j_gt_i) - i, j, i_gt_j = bb.add_t(gt.adjoint(), a=i, b=j, target=i_gt_j) + i, j, j_gt_i = bb.add(gt.adjoint(), a=j, b=i, target=j_gt_i) + i, j, i_gt_j = bb.add(gt.adjoint(), a=i, b=j, target=i_gt_j) bb.free(cast(Soquet, i_gt_j)) bb.free(cast(Soquet, j_gt_i)) i_zero, j_zero = cast( NDArray, bb.add(And().adjoint(), ctrl=i_zero_j_zero, target=i_or_j_zero) ) - j_bits, j_zero = bb.add_t(mcx.adjoint(), controls=bb.split(cast(Soquet, j)), target=j_zero) - i_bits, i_zero = bb.add_t(mcx.adjoint(), controls=bb.split(cast(Soquet, i)), target=i_zero) + j_bits, j_zero = bb.add(mcx.adjoint(), controls=bb.split(cast(Soquet, j)), target=j_zero) + i_bits, i_zero = bb.add(mcx.adjoint(), controls=bb.split(cast(Soquet, i)), target=i_zero) bb.free(cast(Soquet, j_zero)) bb.free(cast(Soquet, i_zero)) From 1bef37b312aea0e0eb43dc6caf7587962f72f73d Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Wed, 7 Aug 2024 11:20:09 -0700 Subject: [PATCH 5/8] Update conftest --- qualtran/conftest.py | 1 + 1 file changed, 1 insertion(+) diff --git a/qualtran/conftest.py b/qualtran/conftest.py index 6545294bf..f66c61bff 100644 --- a/qualtran/conftest.py +++ b/qualtran/conftest.py @@ -115,6 +115,7 @@ def assert_bloq_example_serializes_for_pytest(bloq_ex: BloqExample): 'scaled_chebyshev_poly_odd', 'black_box_select', # cannot serialize AutoPartition 'black_box_prepare', # cannot serialize AutoPartition + 'vlasov_block_encoding' ]: pytest.xfail("Skipping serialization test for bloq examples that cannot yet be serialized.") From ed58e3244fc0770219a9123215cbed7a0d7ce0d1 Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Wed, 7 Aug 2024 12:31:27 -0700 Subject: [PATCH 6/8] Fix format --- qualtran/conftest.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/qualtran/conftest.py b/qualtran/conftest.py index f66c61bff..40d91cc1a 100644 --- a/qualtran/conftest.py +++ b/qualtran/conftest.py @@ -115,7 +115,7 @@ def assert_bloq_example_serializes_for_pytest(bloq_ex: BloqExample): 'scaled_chebyshev_poly_odd', 'black_box_select', # cannot serialize AutoPartition 'black_box_prepare', # cannot serialize AutoPartition - 'vlasov_block_encoding' + 'vlasov_block_encoding', ]: pytest.xfail("Skipping serialization test for bloq examples that cannot yet be serialized.") From 29b8b67d2c5dc4111c0608a601da3b7334118bae Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Mon, 26 Aug 2024 13:22:28 -0700 Subject: [PATCH 7/8] Fix for merge --- .../block_encoding/vlasov_equation.ipynb | 122 ++---------------- .../bloqs/block_encoding/vlasov_equation.py | 2 +- .../resource_counting/t_counts_from_sigma.py | 16 ++- 3 files changed, 25 insertions(+), 115 deletions(-) diff --git a/qualtran/bloqs/block_encoding/vlasov_equation.ipynb b/qualtran/bloqs/block_encoding/vlasov_equation.ipynb index f3f93dd95..bf4838fab 100644 --- a/qualtran/bloqs/block_encoding/vlasov_equation.ipynb +++ b/qualtran/bloqs/block_encoding/vlasov_equation.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -22,7 +22,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -35,39 +35,9 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/markdown": [ - "#### Counts totals:\n", - " - `Allocate`: 1\n", - " - `Allocate`: 1\n", - " - `ArbitraryClifford(n=2)`: $\\displaystyle 221 n + 557$\n", - " - `CNOT`: $\\displaystyle 98 n + 223$\n", - " - `C[2][TwoBitSwap]`: 14\n", - " - `C[2][XGate]`: 7\n", - " - `C[Ry(0.0024867959858108648π)]`: 7\n", - " - `C[Ry(0.435905783151025π)]`: 1\n", - " - `C[TwoBitSwap]`: $\\displaystyle n$\n", - " - `Free`: 1\n", - " - `Free`: 1\n", - " - `H`: 12\n", - " - `Rz(-0.3918265520306073π)`: 2\n", - " - `Rz(0.3918265520306073π)`: 2\n", - " - `T`: $\\displaystyle 8716 n - 32 \\left(n - 2\\right) \\left(30 n - 30\\right) - 72 \\left(n - 1\\right)^{2} - 32 \\left(n - 1\\right) \\left(6 n - 12\\right) + 8 \\left\\lceil{480 n - 60 \\left(n - 1\\right)^{2} + 1486.0}\\right\\rceil + 34084$\n", - " - `TwoBitSwap`: $\\displaystyle n$\n", - " - `XGate`: $\\displaystyle 10 n + 165$" - ], - "text/plain": [ - "" - ] - }, - "metadata": {}, - "output_type": "display_data" - } - ], + "outputs": [], "source": [ "g, sigma = bloq.call_graph(generalizer=ignore_split_join)\n", "show_counts_sigma(sigma)" @@ -75,17 +45,9 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "8716*n - 32*(n - 2)*(30*n - 30) - 72*(n - 1)**2 - 32*(n - 1)*(6*n - 12) + 8*ceiling(480*n - 60*(n - 1)**2 + 1486.0) + 34292\n" - ] - } - ], + "outputs": [], "source": [ "ts = t_counts_from_sigma(sigma)\n", "print(ts)" @@ -93,30 +55,16 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/latex": [ - "$\\displaystyle 25416$" - ], - "text/plain": [ - "25416" - ] - }, - "execution_count": 5, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ "ts.subs(n, 11)" ] }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -131,7 +79,7 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -147,43 +95,9 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/markdown": [ - "#### Counts totals:\n", - " - `Allocate`: 1\n", - " - `Allocate`: 1\n", - " - `Allocate`: 1\n", - " - `ArbitraryClifford(n=2)`: 860046\n", - " - `CNOT`: 131200\n", - " - `C[Ry(0.015625π)]`: 1\n", - " - `C[Ry(0.03125π)]`: 1\n", - " - `C[Ry(0.0625π)]`: 1\n", - " - `C[Ry(0.125π)]`: 1\n", - " - `C[Ry(0.25π)]`: 1\n", - " - `C[Ry(0.5π)]`: 1\n", - " - `C[Ry(π)]`: 1\n", - " - `Free`: 1\n", - " - `Free`: 1\n", - " - `Free`: 1\n", - " - `H`: 12\n", - " - `Rz(-0.3918265520306073π)`: 2\n", - " - `Rz(0.3918265520306073π)`: 2\n", - " - `T`: 229368\n", - " - `TwoBitSwap`: 11\n", - " - `XGate`: 50" - ], - "text/plain": [ - "" - ] - }, - "metadata": {}, - "output_type": "display_data" - } - ], + "outputs": [], "source": [ "g, sigma = bloq.call_graph(generalizer=ignore_split_join)\n", "show_counts_sigma(sigma)" @@ -191,17 +105,9 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "229576\n" - ] - } - ], + "outputs": [], "source": [ "ts = t_counts_from_sigma(sigma)\n", "print(ts)" diff --git a/qualtran/bloqs/block_encoding/vlasov_equation.py b/qualtran/bloqs/block_encoding/vlasov_equation.py index e7a09a332..2bdb33722 100644 --- a/qualtran/bloqs/block_encoding/vlasov_equation.py +++ b/qualtran/bloqs/block_encoding/vlasov_equation.py @@ -38,7 +38,7 @@ from qualtran.bloqs.basic_gates.x_basis import XGate from qualtran.bloqs.block_encoding.sparse_matrix import EntryOracle, SparseMatrix from qualtran.bloqs.mcmt.and_bloq import And -from qualtran.bloqs.mcmt.multi_control_multi_target_pauli import MultiControlPauli +from qualtran.bloqs.mcmt import MultiControlPauli from qualtran.bloqs.rotations.phase_gradient import _fxp from qualtran.resource_counting import BloqCountT, SympySymbolAllocator from qualtran.simulation.classical_sim import ClassicalValT diff --git a/qualtran/resource_counting/t_counts_from_sigma.py b/qualtran/resource_counting/t_counts_from_sigma.py index 07ffc8f8e..0ffd7f8cb 100644 --- a/qualtran/resource_counting/t_counts_from_sigma.py +++ b/qualtran/resource_counting/t_counts_from_sigma.py @@ -21,17 +21,21 @@ def t_counts_from_sigma(sigma: Mapping['Bloq', SymbolicInt]) -> SymbolicInt: """Aggregates T-counts from a sigma dictionary by summing T-costs for all rotation bloqs.""" - from qualtran.bloqs.basic_gates import TGate + from qualtran.bloqs.basic_gates import SGate, TGate from qualtran.cirq_interop.t_complexity_protocol import TComplexity from qualtran.resource_counting.classify_bloqs import bloq_is_rotation ret = sigma.get(TGate(), 0) + sigma.get(TGate().adjoint(), 0) for bloq, counts in sigma.items(): - if bloq_is_rotation(bloq) and not cirq.has_stabilizer_effect(bloq): - if isinstance(bloq, Controlled): - # TODO native controlled rotation bloqs missing (CRz, CRy etc.) - # https://github.com/quantumlib/Qualtran/issues/878 - bloq = bloq.subbloq + while isinstance(bloq, Controlled): + # TODO native controlled rotation bloqs missing (CRz, CRy etc.) + # https://github.com/quantumlib/Qualtran/issues/878 + bloq = bloq.subbloq + if ( + not isinstance(bloq, (TGate, SGate)) + and bloq_is_rotation(bloq) + and not cirq.has_stabilizer_effect(bloq) + ): assert hasattr(bloq, 'eps') ret += ceil(TComplexity.rotation_cost(bloq.eps)) * counts return ret From f2f518ff3f5a64fbe0764fc75a9168aae844c4bd Mon Sep 17 00:00:00 2001 From: Charles Yuan Date: Mon, 26 Aug 2024 15:23:58 -0700 Subject: [PATCH 8/8] Fix format --- qualtran/bloqs/block_encoding/vlasov_equation.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/qualtran/bloqs/block_encoding/vlasov_equation.py b/qualtran/bloqs/block_encoding/vlasov_equation.py index 2bdb33722..0d4f64e1a 100644 --- a/qualtran/bloqs/block_encoding/vlasov_equation.py +++ b/qualtran/bloqs/block_encoding/vlasov_equation.py @@ -37,8 +37,8 @@ from qualtran.bloqs.basic_gates import Ry, Swap, Toffoli from qualtran.bloqs.basic_gates.x_basis import XGate from qualtran.bloqs.block_encoding.sparse_matrix import EntryOracle, SparseMatrix -from qualtran.bloqs.mcmt.and_bloq import And from qualtran.bloqs.mcmt import MultiControlPauli +from qualtran.bloqs.mcmt.and_bloq import And from qualtran.bloqs.rotations.phase_gradient import _fxp from qualtran.resource_counting import BloqCountT, SympySymbolAllocator from qualtran.simulation.classical_sim import ClassicalValT