Theorem 9 assumes the multistep method is of order $p$ and $A(\alpha)$-stable, that $A$ is diagonalisable as $A = VDV^{-1}$, and that its eigenvalues all satisfy $|\arg(-\lambda_i)| \leq \alpha$ — so the eigenvalues must sit in a wedge in the left half-plane, and where $A$ has an eigenvalue with positive real part the paper's remedy is to "simply subtract a multiple of the identity, and rescale the solution". Two further hypotheses: $\max_{t \in [t_0, t_0+\Delta t]} ||x(t)|| = O(||x(t_0+\Delta t)||)$, which rules out a solution that decays away from its own final value, and $\varepsilon = o(||x_{in}||)$. $A(\alpha)$-stability rather than $A$-stability is forced rather than chosen: "the second Dahlquist barrier is that an $A$-stable multistep method must be of order $p \leq 2$", and the point of the paper is $p > 2$. **The rigorous result is for constant coefficients**: the Conclusions state that the approach "can also be used to solve linear differential equations with time-dependent coefficients, though the error analysis will be more difficult", so unlike the Dyson and LCHS routes this one does not meet the slot's $A(t)$ contract with a proof behind it.

A block-encoding of $A(t)$ with a normalization $\alpha_A \ge \max_t \lVert A(t)\rVert$, preparation unitaries for $u_0$ and $b$, the evolution time $T$, and an error tolerance $\varepsilon$.

1. Linear ODE system → Linear system Ax = bChoose a time discretization or propagator approximation

   What crosses here is a system of rows, and the register that indexes them. The step is Eq. (9), $\sum_{\ell=0}^{k}\alpha_{\ell}x_{j+\ell} = h\sum_{\ell=0}^{k}\beta_{\ell}[A(t_{j+\ell})x_{j+\ell} + b(t_{j+\ell})]$, and the trick that makes one linear system out of a whole evolution is a "Feynman's clock": an added register encoding the time, so that $|\psi\rangle := \sum_{j=0}^{N_t}|t_j\rangle|x_j\rangle$ holds the solution at every step rather than at the last one. approximation: The method has order $p$ exactly when $\rho(e^h) - h\sigma(e^h) = O(h^{p+1})$ for the generating polynomials $\rho(\zeta) = \sum_j\alpha_j\zeta^j$, $\sigma(\zeta) = \sum_j\beta_j\zeta^j$; going high-order is the whole point, because Euler at order one forces $N_t = O(\Delta t^2/\varepsilon)$ and "a scaling of the complexity of at least $O(\Delta t^4)$". assumption: The method is $A(\alpha)$-stable — $S \supset S_{\alpha} = \{\mu; |\arg(-\mu)| < \alpha, \mu \neq 0\}$, Definition 2 — rather than $A$-stable, because the second Dahlquist barrier caps an $A$-stable method at order two; and every eigenvalue of $A$ lies in that wedge, with a shift by a multiple of the identity where one does not. **The horizon is doubled on purpose.** Reading the answer off the clock succeeds with probability only $1/(N_t+1)$, so $x$ is held constant over $[t_0+\Delta t, t_0+2\Delta t]$ with $N_t = 2\Delta t/h$ — "then any measurement result for the time in this interval will give the state corresponding to the solution" — which costs nothing in the scaling of $N_t$.

   approximationassumption

   Open this step

2. Linear system Ax = b → Solution as a stateQuantum linear solve

   The assembled system $\mathcal{A}\vec{x} = \vec{b}$ of Eq. (7) goes to a quantum linear system algorithm — Harrow, Hassidim and Lloyd's, whose runtime Berry quotes as $\tilde{O}(\log(N)s^4\kappa^2/\varepsilon_L)$. assumption: The condition number of the assembled system is not the condition number of $A$: Theorem 7 gives $\kappa = O(N_t\kappa_V)$, so the clock register's length enters the solve directly, and $\kappa_V$ — the condition number of the eigenvector matrix $V$ in $A = VDV^{-1}$ — enters in place of any conditioning of $A$ itself. approximation: The solver's error is a trace distance on the state it returns, and the two error budgets are tied by $\varepsilon_L = \Theta(\varepsilon)$; the argument for that is Eqs. (60)-(70), which bound the post-measurement state against the exact one by $O(\varepsilon_L p_{ex}^{-1})$ and then show $p_{ex} = \Omega(1)$. What makes $p_{ex}$ a constant rather than $1/(N_t+1)$ is the doubled horizon of the previous hop, so the two hops are not independent: the discretization pays for the readout. Preparing the right-hand side costs $O(\sqrt{s})$ by Lemma 8, and the superposition over times $O(\log N_t)$ "simply by choosing $N_t$ to be a power of two".

   approximationassumption

   Open this step

A state proportional to $u(T)$, or a history state, together with separately stated matrix-query and initial-state-query complexity.

Read in full: 1010.2745 carries no numerics at all. It is a construction and a complexity analysis end to end — Definitions 1 and 2, the stability polynomials, the block system of Eq. (8), and Theorem 9 — with no instance, no simulation and no figure. The register row records that as a full-text read rather than an abstract one.

given  A (sparse, constant coefficients), b, x_in, horizon dt, error budget e

choose a linear k-step method of order p that is A(alpha)-stable
    # A-stability would cap p at 2 (the second Dahlquist barrier), and the
    # whole point of going high-order is p > 2
    # alpha must be wide enough that every eigenvalue of A has |arg(-lambda_i)| <= alpha

set N_t = 2*dt/h  time steps
    # twice the horizon on purpose: x is held constant on [dt, 2dt], so that
    # measuring the clock register lands on the answer with probability Omega(1)
    # instead of 1/(N_t + 1)

assemble the sparse linear system  script_A x_vec = b_vec  whose rows are:
    row 0:            x_0 = x_in
    stepping rows:    sum_l alpha_l x_{j+l}  =  h * sum_l beta_l ( A x_{j+l} + b )
    holding rows:     x_{j+1} - x_j = 0        # the constant tail

hand the system to the layer below
    # the solution comes back as one state  sum_j |t_j> |x_j>  over a clock
    # register: the whole history at once, not the final time alone

measure the clock register; keep the outcomes in [dt, 2dt]
    # what is left is a state encoding x(dt), to within trace distance e

Theorem 9, the main result: $\tilde{O}( \log(N_x) s^{9/2} (||A||\Delta t)^{2+2/p} \kappa_V^5 (||x_{in}|| + ||b||/||A||)/\varepsilon^2 )$ calls to the oracles for $A$, $b$ and $x_{in}$, where $s$ is the sparsity, $N_x$ the dimension and $\kappa_V$ the condition number of the eigenvector matrix. Dropping the term for the error in *starting* the multistep method improves the other quantities without changing the $||A||\Delta t$ scaling: $\tilde{O}( \log(N_x) s^{9/2} (||A||\Delta t)^{2+2/p} \kappa_V^{2+4/p} (||x_{in}|| + ||b||/||A||)^{1/p}/\varepsilon^{1+2/p} )$. The Conclusions quote the headline in $||A||$ and $\Delta t$ alone: $\tilde{O}( (||A||\Delta t)^2 )$. The paper says outright that this is probably not the last word — "it is likely that this is suboptimal, because the lower bound is linear scaling", the no-fast-forwarding bound $O(||A||\Delta t)$ — and every later route on this slot is an attempt on that gap.

• The clock-register linear system $\mathcal{A}\vec{x} = \vec{b}$ (sections III and IV)

  • High-order quantum algorithm for solving linear differential equations

    Dominic W. Berry · 2010

  About

  Why a linear system rather than a simulation, in the paper's own order. Section II sets up the alternative first and rejects it: for a general $A$, split into Hermitian and antiHermitian parts and combine the two evolutions by a Lie-Trotter formula, and "it appears to give a complexity that increases exponentially with the time interval $\Delta t$" — because the amplitude "must therefore decay exponentially, because we must allow for the maximum eigenvalue of $A_H$". That approach also "does not enable simulation of inhomogeneous differential equations". Section III's answer is the trick the rest of the construction rests on: a "Feynman's clock", an additional register encoding the time, so that one state $|\psi\rangle := \sum_{j=0}^{N_t}|t_j\rangle|x_j\rangle$ holds the solution at all times rather than the final time alone.

  Methods

  The system is $\mathcal{A}\vec{x} = \vec{b}$ of Eq. (7), where "each entry of $\mathcal{A}$ is a block of the dimension of $A$, and each entry of $\vec{x}$ and $\vec{b}$ is a block of the dimension of $x$". Eq. (8) writes it out for the Euler case with $A$ and $b$ time-independent: the first block row is $\mathbb{1}$ against $x_{in}$ and "sets the initial value, $x_0 = x_{in}$"; the stepping rows are $-(\mathbb{1}+Ah)$ on the subdiagonal against $\mathbb{1}$ on the diagonal, so that "$x_{j+1} - (x_j + Ax_jh) = bh$"; and "the final rows indicate equations where $x_{j+1} - x_j = 0$", the tail that holds $x$ constant. Section IV replaces the stepping rows with a general linear multistep method, Eq. (9): $\sum_{\ell=0}^{k}\alpha_{\ell}x_{j+\ell} = h\sum_{\ell=0}^{k}\beta_{\ell}[A(t_{j+\ell})x_{j+\ell} + b(t_{j+\ell})]$, whose stability is read off the generating polynomials $\rho(\zeta) = \sum_j \alpha_j\zeta^j$ and $\sigma(\zeta) = \sum_j \beta_j\zeta^j$ through the roots of $\rho(\zeta) - \mu\sigma(\zeta) = 0$: the stability domain is $S := \{\mu \in \mathbb{C};\ \text{all roots } \zeta_j(\mu) \text{ satisfy } |\zeta_j(\mu)| \leq 1,\ \text{multiple roots satisfy } |\zeta_j(\mu)| < 1\}$, and the method has order $p$ if and only if $\rho(e^h) - h\sigma(e^h) = O(h^{p+1})$. Definition 2 is the property the algorithm needs: $A(\alpha)$-stable for $0 < \alpha < \pi/2$ means $S \supset S_{\alpha} = \{\mu;\ |\arg(-\mu)| < \alpha,\ \mu \neq 0\}$. **The doubling is load-bearing and easy to miss.** Reading the answer off the clock register succeeds with probability only $1/(N_t+1)$, so the paper solves to $t_0 + 2\Delta t$ with $x$ held constant on the second half — "then any measurement result for the time in this interval will give the state corresponding to the solution" — which is what makes $p_{ex} = \Omega(1)$ in Eq. (69) and costs nothing in the scaling of $N_t$.

  Data

  No dataset and no instantiated matrix. $A$, $b$ and $x_{in}$ stay symbolic throughout, reached only through oracles, and the sole concrete display is the five-block-row Euler example of Eq. (8) — drawn to make the block structure legible, with $A$ and $b$ still symbolic. The one quantitative choice the construction fixes is $N_t = 2\Delta t/h$.

  Code

  None found yet.

  Results

  None found yet.

• Quantum algorithms for linear differential equations · Qiskit

  From the repository — run, not written up from a paper · unsupported

  About

  Given a first-order linear differential equation d/dt x = A(t)x + b(t) with N-dimensional vectors x and b and an N×N matrix A, and given an initial condition x(0), produce the solution x(t) at a later time t to precision ε, in the sense that the normalized vector x(t)/‖x(t)‖ returned is at distance at most ε from the exact solution.

  Methods

  None found yet.

  Data

  None found yet.

  Code

  Qiskit

  Results

  Literature record · problem statement and speedup class checked against the Quantum Algorithm Zoo entry · primary reference checked against its arXiv abs page

• Quantum algorithms for linear differential equations

  Given a first-order linear differential equation d/dt x = A(t)x + b(t) with N-dimensional vectors x and b and an N×N matrix A, and given an initial condition x(0), produce the solution x(t) at a later time t to precision ε, in the sense that the normalized vector x(t)/‖x(t)‖ returned is at distance at most ε from the exact solution.