Stated for systems of possibly inhomogeneous linear ODEs with constant coefficients. Berry, Childs, Ostrander and Wang describe the encoding as a sparse, well-conditioned linear system, and state that unlike with finite difference methods their approach does not require additional hypotheses to ensure numerical stability. It reduces to a quantum linear solve; it does not remove that layer.

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

   Over one step of length $h$ the exact solution $\vec{x}(h) = e^{Ah}\vec{x}(0) + (e^{Ah} - I)A^{-1}\vec{b}$ is replaced by $T_k(Ah)\vec{x}(0) + S_k(Ah)h\vec{b}$, where $T_k(z) = \sum_{j=0}^{k} z^j/j!$ and $S_k(z) = \sum_{j=1}^{k} z^{j-1}/j!$. approximation: $T_k(Ah)$ stands in for $e^{Ah}$ and $S_k(Ah)h$ for $(e^{Ah}-I)A^{-1}$. Lemmas 10 and 12 bound the two replacements in their unscaled form, $\lvert T_k(z) - e^{z}\rvert \le 1/(k+1)!$ and $\lvert S_k(z) - (e^{z}-1)z^{-1}\rvert \le 1/(k+1)!$, so the $\vec{b}$ term picks up a factor $h$ when the bound is carried to the step above. Each Taylor term becomes one row of $C_{m,k,p}(Ah)$ built from the row before it, $x_{i,j} = (Ah/j)\,x_{i,j-1}$ for $2 \le j \le k$ with $x_{i,1} = Ah\,x_{i,0} + h\,b$, and the row that closes a step sums them, $x_{i,0} = \sum_{j=0}^{k} x_{i-1,j}$. assumption: The step is short — $h \le 1/\lVert A\rVert$, so $\lVert Ah\rVert \le 1$ — and the lemmas hold only for $\lvert z\rvert \le 1$ with $\operatorname{Re}(z) \le 0$, which is why Theorem 6 is stated for a diagonalizable $A = VDV^{-1}$ whose every eigenvalue satisfies $\operatorname{Re}(\lambda_i) \le 0$, with $k \ge 5$ and $(k+1)! \ge 2m$; its accumulated bound carries $\kappa_V = \lVert V\rVert\,\lVert V^{-1}\rVert$.

   approximationassumption

   Open this step

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

   The system handed down is $C_{m,k,p}(Ah)\,x = \lvert 0\rangle\lvert x_{\mathrm{in}}\rangle + h\sum_{i=0}^{m-1} \lvert i(k+1)+1\rangle\lvert b\rangle$, a $(d+1)N \times (d+1)N$ matrix with $d = m(k+1)+p$ and, for $A$ of sparsity $s$, $O(ks)$ nonzero entries in any row or column; it is nonsingular because it is lower triangular with nonzero diagonal entries. Theorem 5 bounds its condition number by $\kappa_C \le 6\kappa_V k(m+p)$, which at the algorithm's parameter choice $m = p = \lceil T\lVert A\rVert\rceil$ is $O(\kappa_V k m)$. assumption: $A = VDV^{-1}$ is diagonalizable with $\lVert A\rVert \le 1$ and $\operatorname{Re}(\lambda_i) \le 0$ for every eigenvalue, and $k \ge 5$ with $(k+1)! \ge 2m$; the bound is in terms of $\kappa_V = \lVert V\rVert\,\lVert V^{-1}\rVert$, the condition number of the diagonalizing $V$. approximation: The solver returns a state only $\delta$-close to the normalised solution, with $\delta = \varepsilon/(25\sqrt{m}\,g)$ chosen so that the extracted state is $\varepsilon$-close.

   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.

Berry, Childs, Ostrander and Wang report no run of any kind. A full-text read (2026-08-12, which moved the register row for arxiv:1701.03684 from an abstract basis to a full-text one) found no figure, no table, no worked instance with numbers, no dataset, no simulator and no hardware anywhere in the paper — the only occurrences of "simulation" are in the sense of Hamiltonian simulation. There is therefore nothing to write up, and there is no neighbouring numerical work to describe either.

given  A (constant coefficients, diagonalizable A = V D V^-1, Re(lambda_i) <= 0),
       b, x_in with known norms ||x_in||, ||b||, horizon T, error budget e <= 1/2,
       oracles O_A for the entries of A, and controlled O_x, O_b

set    h     = T / ceil(T ||A||)
       m = p = T/h = ceil(T ||A||)
       g     = max_{t in [0,T]} ||x(t)|| / ||x(T)||
       delta = e / (25 sqrt(m) g)                   # the solver's error budget
       Omega = 70 g kappa_V m^{3/2} (||x_in|| + T ||b||) / (e ||x(T)||)
       k     = ceil( 2 log(Omega) / log log(Omega) )   # ensures (k+1)! >= Omega

# --- time-discretization, via truncated-taylor-propagator -------------------
assemble  C_{m,k,p}(Ah) x = |0>|x_in> + h sum_{i=0..m-1} |i(k+1)+1>|b>
    Taylor-term rows :  x_{i,j} = (Ah/j) x_{i,j-1},   2 <= j <= k
                        x_{i,1} = Ah x_{i,0} + h b
    step-closing row :  x_{i,0} = sum_{j=0..k} x_{i-1,j}
    padding rows     :  x_{m,j} = x_{m,j-1},          1 <= j <= p
    # sparse and well-conditioned: O(ks) nonzeros per row or column, and
    # condition number O(kappa_V k m) by Theorem 5.  Unlike with finite
    # difference methods, no additional hypothesis is needed to ensure
    # numerical stability.

repeat O(g) times, coherently -- the amplitude amplification rounds:

    # --- own work: prepare the right-hand side ---
    prepare |0>|x_in> + h sum_{i=0..m-1} |i(k+1)+1>|b>
        # a constant number of calls to O_x and O_b, poly(log(mk)) gates

    # --- quantum-linear-solve ---
    solve C_{m,k,p}(Ah) x = (that state) with a QLSA, to within delta of the
        normalised solution
        # this route reduces to a quantum linear solve; it does not remove
        # that layer

    # --- own work: extract the answer ---
    measure the first register in the standard basis
    accept iff the outcome lies in S = {m(k+1), ..., m(k+1)+p}
        # the last p+1 blocks all hold x_{m,0}, the state at time T
        # one round accepts with probability at least 1/78 g^2  (p = m)

return the second register -- a state e-close to x(T)/||x(T)||, with a flag
       indicating success

The complexity is polynomial in the logarithm of the inverse error, an exponential improvement over previous quantum algorithms for this problem. No matrix-query or state-preparation query count appears in the abstract; a $\kappa_V$-dependent expression often attached to this method in secondary summaries belongs to the spectral-method row of a later comparison table, so it is not reproduced here.

• Krovi's reanalysis of the all-at-once encoding

  Reanalyses the all-at-once propagator encoding and shows that the norm of the matrix exponential, rather than the eigenvector condition number, characterizes the run time. It still forms a global linear system and still calls a quantum linear solver.

  What would earn it a path of its own: The paper's contribution is analytic — the norm of the matrix exponential, not the eigenvector condition number, characterizes the run time — and it chooses no step of its own. Nothing in the walk differs, so a drawn path would need the map to represent analysis-level distinctions (which bound governs a lane), which it does not draw today.

• The sparse-oracle construction of Theorem 9 (sections 2, 6 and 7)

  • Quantum algorithm for linear differential equations with exponentially improved dependence on precision

    Dominic W. Berry, Andrew M. Childs, Aaron Ostrander, Guoming Wang · 2017

  About

  The algorithm the paper states and proves as its main result. Section 1 gives the reason it was built: the earlier approach handed a linear multistep method to a quantum linear systems algorithm, and its complexity stayed $\mathrm{poly}(1/\varepsilon)$ even when the improved solver was substituted, because the multistep method is itself a significant source of error. Direct use of a linear combination of unitaries was considered and rejected — for general non-anti-Hermitian $A$ the best implementation the authors are aware of has an exponentially small success probability — so a truncation of the Taylor series of the propagator $\exp(At)$ is encoded into a linear system instead, which the paper describes as conceptually similar to a linear combination of unitaries but with significantly better performance because the combination is built stepwise through a system of linear equations. The problem statement is Theorem 9's: $A = VDV^{-1}$ is an $N \times N$ diagonalizable matrix whose eigenvalues all satisfy $\operatorname{Re}(\lambda_j) \le 0$, with at most $s$ nonzero entries in any row and column and an oracle $O_A$ computing them; $\vec{x}_{\mathrm{in}}$ and $\vec{b}$ have known norms and controlled oracles $O_x$ and $O_b$ preparing states proportional to them; the output is a state $\varepsilon$-close in $\ell_2$ norm to $\vec{x}(T)/\lVert\vec{x}(T)\rVert$ with a flag indicating success, obtained with probability $\Omega(1)$.

  Methods

  The parameters are all fixed in one place, the Statement of the Algorithm in section 7: $h = T/\lceil T\lVert A\rVert\rceil$, $m = p = T/h = \lceil T\lVert A\rVert\rceil$, $\delta = \varepsilon/(25\sqrt{m}\,g)$ with $\varepsilon \le 1/2$, and $k = \lceil 2\log\Omega/\log\log\Omega\rceil$ where $\Omega = 2.8\kappa_V m(\lVert x_{\mathrm{in}}\rVert + mh\lVert b\rVert)/(\delta\lVert x(T)\rVert) = 70g\kappa_V m^{3/2}(\lVert x_{\mathrm{in}}\rVert + T\lVert b\rVert)/(\varepsilon\lVert x(T)\rVert)$. That choice of $k$ ensures $(k+1)! \ge \Omega$, and because $\Omega \ge 70$ it also ensures $k \ge 5$, which is the hypothesis the conditioning and success-probability theorems need. The system built is $C_{m,k,p}(Ah)\lvert x\rangle = \lvert 0\rangle\lvert x_{\mathrm{in}}\rangle + h\sum_{i=0}^{m-1}\lvert i(k+1)+1\rangle\lvert b\rangle$. Its right-hand side is prepared by Lemma 8 of section 6, in three moves: a rotation on the block register that weights $\lVert\vec{x}_{\mathrm{in}}\rVert$ against $\sqrt{m}h\lVert\vec{b}\rVert$, then $O_x$ and $O_b$ applied in arbitrary order, then a unitary that fixes $\lvert 0\rangle$ and sends $\lvert 1\rangle$ to $\frac{1}{\sqrt{m}}\sum_{j=0}^{m-1}\lvert j(k+1)+1\rangle$, spreading the inhomogeneous part over the time-step blocks — a constant number of calls to $O_x$ and $O_b$ plus $\mathrm{poly}(\log(mk))$ elementary gates. The system is solved with the linear systems algorithm of Childs, Kothari and Somma, which section 1 singles out because it avoids phase estimation and inverts the system directly by a linear combination of unitaries, giving a state within $\delta$ of the normalised solution. The first register is then measured in the standard basis and the outcome accepted when it lies in $S = \{m(k+1), m(k+1)+1, \ldots, m(k+1)+p\}$, the block labels that all carry $\lvert x_{m,0}\rangle$. Theorem 7 with $p = m$ puts one round's success probability at $1/(78g^2)$ for the exact solution vector, and Eq. (127) at $1/(121g^2)$ once the solver's own error $\delta$ is carried through. $O(g)$ rounds of amplitude amplification raise it to $\Omega(1)$; section 2 records classical repetition as the alternative to amplification.

  Data

  No dataset. $A$ enters only through the oracle $O_A$ that computes its at most $s$ nonzero entries per row and column, and the two vectors only through controlled oracles $O_x$ and $O_b$ preparing $\vec{x}_{\mathrm{in}}/\lVert\vec{x}_{\mathrm{in}}\rVert$ and $\vec{b}/\lVert\vec{b}\rVert$, with both norms assumed known. Section 6 marks this input model as a deliberate change from the multistep algorithm it is compared against, which assumes $\vec{x}_{\mathrm{in}}$ and $\vec{b}$ are sparse vectors whose entries are given by oracles: instead of assuming sparsity, this paper simply assumes controlled oracles that produce states proportional to the two vectors.

  Code

  None found yet.

  Results

  Counted, never run. The solve alone costs $O(\kappa_V k^2 m s \cdot \mathrm{poly}(\log(\kappa_V kms/\delta)))$ queries, using Theorem 5's $O(\kappa_V km)$ bound on the condition number of $C_{m,k,p}(Ah)$ at these parameters and the fact that the matrix is $(d+1)N \times (d+1)N$ with $O(ks)$ nonzero entries in any row or column. With the $O(g)$ amplification rounds, Theorem 9 gives the whole algorithm $O(\kappa_V sgT\lVert A\rVert \cdot \mathrm{poly}(\log(\kappa_V sg\beta T\lVert A\rVert/\varepsilon)))$ queries to $O_A$, $O_x$ and $O_b$, where $\beta = (\lVert x_{\mathrm{in}}\rVert + T\lVert b\rVert)/\lVert x(T)\rVert$, and a gate complexity larger by a factor $\mathrm{poly}(\log(\kappa_V sg\beta T\lVert A\rVert N/\varepsilon))$. Section 8 argues that neither of the two parameters this is charged in can be improved much: evolution in which the state decays can be used to implement postselection, and the ability to postselect on an exponentially small amplitude would imply BQP = PP, so the dependence on $g$ is close to necessary, while Hamiltonian simulation is the special case $g = 1$ and the no-fast-forwarding theorem then makes the dependence on $T$ optimal up to logarithmic factors. Nothing here was measured. The paper reports no numerical experiment, contains no figure and no table, and names no hardware and no simulator, so this implementation exists as a specified and analysed algorithm rather than an executed one.

• Non-diagonalizable $A$ through a nearby diagonalizable matrix (section 8)

  • Quantum algorithm for linear differential equations with exponentially improved dependence on precision

    Dominic W. Berry, Andrew M. Childs, Aaron Ostrander, Guoming Wang · 2017

  About

  Theorem 9 assumes $A$ is diagonalizable, and section 8 states that the algorithm can nonetheless produce approximate solutions when it is not. This is the paper's own variant for that case, offered as a remark rather than as a second theorem.

  Methods

  Diagonalizable matrices are dense within the set of all complex matrices: for any non-diagonalizable $A$ and any $\delta > 0$ there is a diagonalizable $B$ with $\lVert A - B\rVert < \delta$. The construction is then run unchanged with $B$ in place of $A$, so that Eq. (1) is simulated approximately.

  Data

  None found yet.

  Code

  None found yet.

  Results

  The paper states the price in the same sentence it offers the variant: this approach can yield a matrix $B$ whose diagonalizing transformation has a condition number polynomial in $1/\delta$, so the complexity would no longer be $\mathrm{poly}(\log(1/\varepsilon))$. The variant keeps the algorithm applicable and gives up the property the paper was written for. No bound is proved for it and no numerics are reported.

• Hamiltonian simulation as the $g = 1$ specialisation (sections 1 and 8)

  • Quantum algorithm for linear differential equations with exponentially improved dependence on precision

    Dominic W. Berry, Andrew M. Childs, Aaron Ostrander, Guoming Wang · 2017

  About

  The specialisation the paper names for itself. Section 1 opens on differential equations of the form $\dot{\vec{x}} = A\vec{x}$ with $A$ anti-Hermitian, and section 8 says outright that Hamiltonian simulation is a special case of this algorithm with $g = 1$.

  Methods

  The inhomogeneity is absent and the generator is anti-Hermitian, so the solution does not decay and the parameter $g = \max_{t \in [0,T]}\lVert\vec{x}(t)\rVert/\lVert\vec{x}(T)\rVert$ of Eq. (111) — the quantity the success probability and the $O(g)$ amplification rounds are charged against — takes the value one.

  Data

  None found yet.

  Code

  None found yet.

  Results

  The paper uses the specialisation as an argument rather than as an algorithm to run: because Hamiltonian simulation is this case, the no-fast-forwarding theorem implies that the algorithm's dependence on the evolution time $T$ is optimal up to logarithmic factors. Nothing is computed, and no separate construction is written down for it.