Measuring Trotter error and its application to precision-guaranteed Hamiltonian simulations

Abstract

Trotterization is the most common and convenient approximation method for Hamiltonian simulations on digital quantum computers, but estimating its error accurately is computationally difficult for large quantum systems. Here, we develop a method for measuring the Trotter error without ancillary qubits on quantum circuits by combining the $m$th- and $n$th-order ($m<n$) Trotterizations rather than consulting with mathematical error bounds. Using this method, we make Trotterization precision-guaranteed, developing an algorithm named Trotter$(m,n)$, in which the Trotter error at each time step is within an error tolerance $\epsilon$ preset for our purpose. Trotter$(m,n)$ is applicable to both time-independent and -dependent Hamiltonians, and it adaptively chooses almost the largest stepsize $\delta{t}$, which keeps quantum circuits shallowest, within the error tolerance. Benchmarking it in a quantum spin chain, we find the adaptively chosen $\delta{t}$ to be about ten times larger than that inferred from known upper bounds of Trotter errors.

I Introduction

The rapid development of quantum devices in recent years has led researchers to find useful applications with significant quantum advantage [Harrow2017, Cao2019, Cerezo2021]. On top of the eigenvalue problems [Aspuru-Guzik2005, Lloyd1996, Whitfield2011, Higgott2019, Jones2019, Kirby2021], quantum many-body dynamics, or Hamiltonian simulation [Babbush2015, Low2017, Campbell2019, An2021, An2022, Childs2022], is one of the most promising candidates because quantum computers could overcome the exponential complexity that classical computers face [Lloyd1996], enabling us to address intriguing dynamical phenomena like nonequilibrium phases of matter [Swingle2016, Zhang2017, Landsman2019, Joshi2022] and to implement fundamental quantum algorithms like phase estimation [Kitaev1995]. Among several algorithms for the Hamiltonian simulation, Trotterization [Trotter1959, Hatano2005] is and will be used most commonly in the current noisy intermediate-scale quantum (NISQ [Preskill2018]) era and the coming early fault-tolerant quantum computing (FTQC) era because it does not demand additional ancillary qubits or largely controlled quantum gates. Indeed, quantum advantage in Trotterized dynamics simulation has been reported using a 127-qubit NISQ computer only recently [Kim2023].

One major and presumably inevitable issue of Trotterization is the trade-off relation between the simulation accuracy and the circuit depth (see, however, Refs. [Barison2021, Mansuroglu2023] for variational approaches). The $k$-th order Trotterization accompanies an error of $O(\delta{t}^{k+1})$ during a single time step $\delta{t}$, which decreases when $\delta{t}$ is taken shorter. In the meantime, the number of steps to reach a final time increases, meaning a deeper quantum circuit. To suppress the gate depth, it is desirable to choose the largest possible stepsize $\delta{t}$, i.e., the shallowest circuit, within our error tolerance $\epsilon$ preset for our purposes.

However, it is difficult to find the optimal stepsize $\delta{t}$ because the Trotter error is complex in generic many-body systems. According to the previous studies on the Trotter error, its upper bounds [Kivlichan2020, Childs2021] and typical values [Zhao2021] are available. If we choose $\delta{t}$ so that the upper bound is below our tolerance $\epsilon$, the precision is guaranteed, but $\delta{t}$ tends to be too small, as we will see below. On the other hand, if we choose $\delta{t}$ based on the typical values, $\delta{t}$ can be larger, but the precision guarantee is lost. Recently, Zhao et al. [Zhao2022] proposed an approach where $\delta t$ is chosen adaptively in each time step based on the energy expectation value and variance. Yet, the precision guarantee of this method is still elusive, and the applicability is limited to time-independent Hamiltonians .

In this paper, we develop a method to measure the Trotter error on quantum circuits by combining Trotterization formulas at different orders, $m$ and $n(>m)$. Since measured, the estimated error is significantly more accurate than known upper bounds for it and thus allows us to accurately choose the largest possible stepsize $\delta{t}$ so that the error does not exceed our tolerance $\epsilon$. Using this method, we make Trotterization precision-guaranteed, in which almost the largest $\delta{t}$ is adaptively chosen within a preset error tolerance $\epsilon$ (see Fig. 1 for its structure). We name this algorithm Trotter$(m,n)$ in analogy to RFK45 (Runge-Kutta-Fehlberg) for classical simulations, where the fourth- and fifth-order methods are combined. We benchmark Trotter24 in a quantum spin chain under time-independent and -dependent Hamiltonians, finding that the adaptively chosen $\delta{t}$ is about ten times larger than that inferred from the upper bound of Trotter errors.

Figure 1: Key concept of (fidelity-based) Trotter24. (Left) At each time step $N$, we ask what the optimal stepsize $\delta t_{N}$ for the second-order Trotterization formula $T_{2}(\delta t_{N})$. Here the optimal means the largest with the fidelity error kept less than our tolerance $\epsilon$. (Middle) The fidelity error that $T_{2}(\delta t_{N})$ can be measured using the fourth-order Trotterization $T_{4}(\delta t_{N})$, instead of the exact evolution $U(\delta t_{N})$ when we neglect higher-order corrections in terms of $\delta t_{N}$. (Right) If the measured error is below our tolerance, we accept the value of $\delta t_{N}$ and proceed to the next time step $N+1$. Otherwise, we reject the value and try the same protocol with a smaller value for $\delta t_{N}$ again. In the paper, we also develop an observable-based, rather than fidelity-based, algorithm, establish an efficient scheme for avoiding the possible rejections of $\delta t_{N}$’s, and analyze how errors propagate with time steps.

II Measuring Trotter error

In this section, we present a way to measure a Trotter error on quantum circuits without ancillary qubits when a quantum state is evolved by an $m$-th order Trotterization for a small time step $\delta{t}$. A key idea is making use of a higher-order $n(>m)$-th Trotterization to approximate the unknown exact solution accurately enough.

II.1 Extracting Trotter error using different orders

For simplicity, we first consider a time-independent Hamiltonian $H$ consisting of two parts,

$$H=A+B,$$ (1)

where $A$ and $B$ do not necessarily commute with each other. Generalization to more noncommuting parts is straightforward, and we will generalize the arguments to time-dependent ones later in Sec. LABEL:sec:tdep. We assume that the quantum state at time $t$ is known to be $\ket{\psi(t)}$ and consider evolving it by a small time step $\delta t$,

$$\ket{\psi(t+\delta t)}=U(\delta{t})\ket{\psi(t)}=e^{-iH\delta t}\ket{\psi(t)}.$$ (2)

Trotterization approximately decomposes $e^{-iH\delta{t}}$ into quantum-gate-friendly parts consisting of either $A$ and $B$. We take an $m$-th order Trotterization $T_{m}(\delta{t})$. For example, $T_{m}(\delta{t})$ can be the Lie-Trotter formula $T_{1}(\delta{t})=e^{-iA\delta{t}}e^{-iB\delta{t}}$ for $m=1$ or its symmetrized form $T_{2}(\delta t)\equiv e^{-iA\delta t/2}e^{-iB\delta t}e^{-iA\delta t/2}$ for $m=2$. In general, $T_{m}(\delta{t})$ approximates $U(\delta{t})$ up to the order of $\delta{t}^{m}$ and satisfies

$$T_{m}(\delta t)=e^{-iH\delta t+\Upsilon_{m+1}}$$ (3)

where $\Upsilon_{m+1}=O(\delta t^{m+1})$ is an anti-Hermitian error operator. These relations imply

$$\ket{\psi_{m}(t+\delta t)}\equiv T_{m}(\delta t)\ket{\psi(t)}=\ket{\psi(t+\delta t)}+O(\delta{t}^{m+1}),$$ (4)

meaning that $T_{m}(\delta{t})$ approximates the exact one-step evolution within an error of $O(\delta{t}^{m+1})$.

To quantify the error arising in the one step, we adopt the square root of the infidelity

$${\color[rgb]{0,0,0}\eta_{F}\equiv\sqrt{1-|\braket{\psi(t+\delta t)}{\psi_{m}(t+\delta t)}|^{2}}.}$$ (5)

We can also use other quantities depending on our purposes and make parallel arguments. For example, when we are interested in the expectation value of an observable $O$, we care about the error in it,

$$\eta_{O}\equiv\braket{\psi(t+\delta t)}{O}{\psi(t+\delta t)}-\braket{\psi_{m}(t+\delta t)}{O}{\psi_{m}(t+\delta t)}.$$ (6)

In either case, calculating $\eta_{F}$ or $\eta_{O}$ is difficult because we do not know the exactly evolved state $\ket{\psi(t+\delta{t})}$.

We remark that both $\eta_{F}$ and $\eta_{O}$ are $O(\delta{t}^{m+1})$, but showing $\eta_{F}=O(\delta{t}^{m+1})$ is less obvious. To show this, we note that the leading $O(\delta{t}^{m+1})$ term of $1-\braket{\psi(t+\delta t)}{\psi_{m}(t+\delta t)}$ is pure-imaginary as shown in Appendix LABEL:app:error, and its leading-order contribution of $\eta_{F}$ is given by

$$\eta_{F} =\sqrt{\braket{\psi(t)}{(i\Upsilon_{m+1})^{2}}{\psi(t)}-\braket{\psi(t)}{(i\Upsilon_{m+1})}{\psi(t)}^{2}}$$

$$\qquad+O(\delta{t}^{m+2}),$$ (7)

where $i\Upsilon_{m+1}$ is Hermitian. Equation (II.1) dictates that $\eta_{F}$ is the variance of the “observable” $i\Upsilon_{m+1}$, giving a way to estimate $\eta_{F}$ using $\ket{\psi(t)}$ and the explicit form of $i\Upsilon_{m+1}$. Indeed this is a possible way of measuring $\eta_{F}$, but it requires, for generic many-body Hamiltonians, measuring numerous Hermitian operators involved in $i\Upsilon_{m+1}$ consisting of doubly nested commutators between $A$ and $B$. Hence, in the following subsection, we will also discuss another way to estimate $\eta_{F}$ with less sampling costs.

Our idea of estimating the errors is the following: In calculating $\eta_{F}$ and $\eta_{O}$ in the leading order, we can safely replace the exact $\ket{\psi(t+\delta{t})}$ by a higher-order approximant $\ket{\psi_{n}(t+\delta{t})}$ for $n>m$. Replacing $\ket{\psi(t+\delta{t})}$ by $\ket{\psi_{n}(t+\delta{t})}$ in $\eta_{F}$ and $\eta_{O}$, we obtain the following key analytical results (see Appendix LABEL:app:error for derivation): For the fidelity error,

$$\eta_{F} =\eta_{F}^{(mn)}+O(\delta{t}^{n+1}),$$ (8)

$$\eta_{F}^{(mn)} \equiv\sqrt{1-|\braket{\psi_{n}(t+\delta t)}{\psi_{m}(t+\delta t)}|^{2}},$$ (9)

and, for the observable error,

$$\eta_{O} =\eta_{O}^{(mn)}+O(\delta{t}^{n+1}),$$ (10)

$$\eta_{O}^{(mn)} \equiv\braket{\psi_{n}(t+\delta t)}{O}{\psi_{n}(t+\delta t)}$$

$$\qquad-\braket{\psi_{m}(t+\delta t)}{O}{\psi_{m}(t+\delta t)}.$$ (11)

Given that $\eta_{F}=O(\delta{t}^{m+1})$ and $\eta_{O}=O(\delta{t}^{m+1})$, these results mean that $\eta_{F}^{(mn)}$ ($\eta_{O}^{(mn)}$) coincides with $\eta_{F}$ ($\eta_{O}$) in the leading order since $m<n$. If $n\leq m$, the above equations hold true, but $O(\delta{t}^{n+1})$ contributions are nonnegligible, and $\eta^{(mn)}$ do not give good estimates for $\eta$.

Unlike $\eta_{F}$ and $\eta_{O}$, $\eta_{F}^{(mn)}$ and $\eta_{O}^{(mn)}$ consist of $T_{m}(\delta{t})$ and $T_{n}(\delta{t})$ and are thereby implementable in quantum circuits. In other words, we can estimate the deviation from the exact solution induced by $T_{m}(\delta{t})$ without knowing the solution when supplemented with the fourth-order Trotterization and neglect higher-order corrections. We will discuss in more detail how to measure $\eta_{F}^{(mn)}$ and $\eta_{O}^{(mn)}$ in Sec. II.2.

We emphasize that $\eta_{F}$ and $\eta_{O}$ are the actual Trotter error specific to the current state $\ket{\psi(t)}$. This contrasts the upper-bound arguments on the operator difference $U(\delta{t})-T_{m}(\delta{t})$ [Kivlichan2020, Childs2021]. Such upper bounds apply to arbitrary states and are thus always larger than or equal to the error occurring at a specific state $\ket{\psi(t)}$. The fact that $\eta_{F}$ and $\eta_{O}$ are state-dependent enables us to choose $\delta{t}$ more accurately so that the error is below our tolerance, as we will see in detail below.

We highlight two sets of $(m,n)$ of particular interest in practical use. The first choice is the minimum pair $(m,n)=(1,2)$, for which the Trotterizations $T_{m}$ and $T_{n}$ involve the minimum possible exponentials, i.e., the gate complexity in the quantum circuit. For noisy circuits, such lowest-order Trotterizations are commonly used to avoid gate errors as much as possible. The second choice is the pair of minimum even numbers $(m,n)=(2,4)$, which can be useful when we can use more gates to achieve higher accuracy. Changing $m$ from 1 to 2 increases the number of exponentials in $T_{m}$ from 2 to 3, by which the Trotter error $\eta_{F}$ or $\eta_{O}$ becomes one-order smaller. In this case, using $n=4$ rather than $n=3$ could be beneficial. To see this, let us compare Ruth’s third-order formula [Ruth1983, Hatano2005],

$$T_{3}(\delta t) \equiv e^{-i\frac{7}{24}A\delta t}e^{-i\frac{2}{3}B\delta t}e^{-i\frac{3}{4}A\delta t}e^{i\frac{2}{3}B\delta t}$$

$$\qquad\times e^{i\frac{1}{24}A\delta t}e^{-iB\delta t}$$ (12)

and the fourth-order Forest-Ruth-Suzuki (FRS) formula [Forest1990, Suzuki1990]

$$T_{4}(\delta t) \equiv e^{-i\frac{s}{2}A\delta t}e^{-isB\delta t}e^{-i\frac{1-s}{2}A\delta t}e^{-i(1-2s)B\delta t}$$

$$\qquad\times e^{-i\frac{1-s}{2}A\delta t}e^{-isB\delta t}e^{-i\frac{s}{2}A\delta t},$$ (13)

where $s=(2-2^{1/3})^{-1}$. Here we note that the FRS formula involves only one extra exponential but is one-order more accurate than Ruth’s third-order formula. Since $\eta^{(24)}$ estimates $\eta^{(23)}$ one-order more accurately at the expense of an extra exponential, it can be worth using $n=4$ for $m=2$.

II.2 Quantum circuit implementations

In Sec. II.1 we derived two expressions, Eqs. (II.1) and (8), for estimating the fidelity error $\eta_{F}$ and one (10) for the observable error $\eta_{O}$. We are assuming a quantum advantage regime, where the number of qubits is so large that $\ket{\psi(t)}$ cannot be stored in classical computer memory but is realized on a quantum circuit. In such a regime, all the expressions for $\eta_{F}$ and $\eta_{O}$ cannot be evaluated with classical linear algebraic computations, and we need to develop ways to evaluate them on quantum circuits.

First, let us consider how to measure the leading-order contribution in Eq. (II.1) of the fidelity error $\eta_{F}$, assuming that $\ket{\psi(t)}$ is realized on a quantum circuit accurately enough. For concreteness we focus on the first-order Trotterization $m=1$, for which $T_{1}(\delta{t})=e^{-iA\delta{t}}e^{-iB\delta{t}}$ and $i\Upsilon_{2}=-i[A,B]\delta{t}^{2}$. Since $i\Upsilon_{2}$ and $(i\Upsilon_{2})^{2}$ are Hermitian, one can evaluate the expectation values for them based on samplings (we will discuss the sampling cost below).

Although this method is useful for small quantum systems, it becomes quadratically costly for larger systems. To see this we consider a Hamiltonian $H$ with $A$ and $B$ are linear combinations of $L$ distinct local Pauli strings. Thanks to the commutator scaling [Childs2021], $i\Upsilon_{2}$ consists of $O(L)$ distinct Pauli strings, and its expectation values are evaluated by estimating the expectation values of each Pauli string. Notably, the translation symmetries can reduce the number of combinations of measurements for expectation value estimations, even down to $O(1)$. However, $(i\Upsilon_{2})^{2}$ is no more a commutator and involves $O(L^{2})$ Pauli strings. In most cases symmetries cannot reduce the number of measured operators to be subextensive, and the sampling cost tends to be extensive. This extensive sampling cost also appears in the energy variance estimation in Ref. [Zhao2022, Zhao2023]. Thus Eq. (II.1) tends to be useful only in small quantum systems.

In contrast, $\eta_{F}^{(mn)}$ in Eq. (8) avoids measuring numerous distinct Pauli strings. To show this we write $\ket{\psi(t)}$ as $\ket{\psi(t)}=U_{T}(t)\ket{\psi(0)}=U_{T}(t)U_{\mathrm{prep}}\ket{\bm{0}}$. Here $\ket{\bm{0}}$ denotes the initialized state, $U_{\mathrm{prep}}$ is an initial state preparation unitary, and $U_{T}(t)$ is some Trotterized unitary propagation from time 0 to $t$. We assume $U_{\mathrm{prep}}$ and $U_{T}(t)$ have appropriate circuit realizations. With these notations $\eta_{F}^{(mn)}$ can be obtained by

$$\eta_{F}^{(mn)} =\sqrt{1-p_{\bm{0}}},$$ (14)

$$p_{\bm{0}} \equiv|\braket{\bm{0}}{U_{\mathrm{prep}}^{\dagger}U_{T}^{\dagger}(t)T_{n}^{\dagger}(\delta{t})T_{m}(\delta{t})U_{T}(t)}{\bm{0}}|^{2}.$$ (15)

Our assumptions tell us that $U_{\mathrm{prep}}^{\dagger}U_{T}^{\dagger}(t)T_{n}^{\dagger}(\delta{t})T_{m}(\delta{t})U_{T}(t)$ has a circuit realization, so $p_{\bm{0}}$ can be interpreted by the probability of finding $\ket{\bm{0}}$ after $\ket{\bm{0}}$ being evolved by the circuit (see also the middle panel of Fig. 1). Note that $p_{\bm{0}}$ is a kind of the Lochmidt echo [Wisniacki2012]. A challenge is that $p_{\bm{0}}$ becomes exponentially small when the system size $L$ increases, so this procedure makes practical sense only in small systems. This issue is not technical but intrinsic in the fidelity error because it is the most stringent quantifier among many-body wave function errors.

The observable error $\eta_{O}$ is more useful for extensive quantum systems. Typically we are interested in local Pauli strings or their linear combinations of $O(L)$. The expectation value of each Pauli string can be estimated by measuring the wave function $\ket{\psi(t)}$ in its eigenbasis. Unlike fidelity, the expectation value does not decay exponentially with the system size. Also, symmetries, such as the translation one, can reduce the number of measured strings to possibly $O(1)$. When one is only interested in local observables, using $\eta_{O}$ is the cheapest option to guarantee precision.

In each method, the Trotter error estimation is based on sampling, and the statistical error can be a bottleneck to achieve high accuracy. The statistical error in general scales as $\propto\mathcal{N}^{-1/2}$, with $\mathcal{N}$ denoting the number of samples (or shots). For the error estimation to be successful, this error should be smaller enough than the error, $\eta_{F}$ or $\eta_{O}$, so $\delta{t}^{m+1}\gtrsim\mathcal{N}^{-1/2}$ must hold true. When we use $\eta^{(mn)}$ as an estimator, the estimation error is determined by the larger of $\sim\mathcal{N}^{-1/2}$ and $\sim\delta{t}^{n+1}$ (see Eqs. (8) and (10)). So, one cannot make the error estimation infinitely accurate by increasing $n$ when the available number of shots $\mathcal{N}$ is limited finitely. Rather, $\delta{t}^{n+1}\sim\mathcal{N}^{-1/2}$ would hold for a reasonably chosen set of $\delta{t}$, $n$, and $\mathcal{N}$. In Section LABEL:sec:finiteShots, we will benchmark our way to estimate $\eta_{O}$ for $(m,n)=(2,4)$ in an example spin system and show that it works with $\mathcal{N}\sim 10^{5}(10^{7})$ for $\epsilon_{O}=10^{-2}(10^{-3})$.

III Precision-guaranteed Trotterization

In the previous section, we developed the methods to evaluate Trotter errors on quantum circuits. Given that the Trotter error is known, one can make Trotterization precision-guaranteed: In each time step, one can make sure that the Trotter error is within a preset accuracy target $\epsilon$. In this section we develop such an algorithm consisting of error measurements and step size optimization, as illustrated in Fig. 1 (The figure is for the fidelity version, but it works, in parallel, for the observable version.).

Our algorithm uses either $\eta_{F}^{(mn)}$ and $\eta_{O}^{(mn)}$ as the Trotter error estimator, and we name it Trotter$(m,n)$. For concreteness we set $(m,n)=(2,4)$ and describe Trotter24 since the generalization to other $(m,n)$ is straightforward (we may use Trotter$(m,n)$ and Trotter$mn$ interchangeably). In current NISQ devices, $(m,n)=(1,2)$ could be more realistic. Since the argument goes in parallel, we first focus on the fidelity error and will address the observable error later in this section.

Our overall task is to simulate the time evolution according to the Hamiltonian $H$ from the initial time $t_{\mathrm{ini}}$ to the final time $t_{\mathrm{fin}}$, starting from an initial state $\ket{\psi_{0}}$. We set an error tolerance $\epsilon$ for the fidelity error in each time step. Initially, we have no a priori information about the appropriate time step, so take a reasonably small trial stepsize $\delta{t}_{0}$, say, $\delta{t}_{0}=0.1J^{-1}$ with $J$ being a typical energy scale of $H$. One could also choose $\delta{t}_{0}$ so small that $T_{2}(\delta{t}_{0})$ never gives larger error than our tolerance $\epsilon$, as guaranteed by a mathematical bound (see Eq. (LABEL:eq:dt_bound) and Appendix LABEL:sec:errorbound for detail). One can also use this $\delta{t}_{0}$ to know the upper bound, in advance, for the required quantum resources for the calculation, although they tend to be too pessimistic, as we will see below.

For the trial $\delta{t}_{0}$ taken in either way, we implement $T_{2}(\delta{t}_{0})$ and $T_{4}(\delta{t}_{0})$ and calculate $\eta_{F}^{(24)}$ using a quantum circuit. Basically, we aim the stepsize to be so small that

$$\eta_{F}^{(24)}<\epsilon.$$ (16)

If this is true, we accept our trial $\delta{t}_{0}$ and evolve our state as $\ket{\psi_{2}}=T_{2}(\delta{t}_{0})\ket{\psi_{0}}$. If $\eta_{F}^{(24)}\geq\epsilon$ instead, our trial $\delta{t}_{0}$ is too large and we need a smaller $\delta{t}_{0}^{\prime}$. In choosing $\delta{t}_{0}^{\prime}$ appropriately, we invoke the leading-order scaling relation $\eta_{F}^{(24)}\approx\alpha\delta{t}_{0}^{3}$ for some unknown $\alpha$ independent of $\delta{t}_{0}$. We can use this relation to estimate $\alpha$ by $\alpha\approx\eta_{F}^{(24)}/\delta{t}_{0}^{3}$ since we measured $\eta_{F}^{(24)}$. For $\delta{t}_{0}^{\prime}$, we expect $\eta_{F}^{(24)}{}^{\prime}\approx\alpha(\delta{t}_{0}^{\prime})^{3}\approx\eta_{F}^{(24)}(\delta{t}_{0}^{\prime}/\delta{t}_{0})^{3}$, which we wish is smaller than $\epsilon$. Thus, the condition $\eta_{F}^{(24)}{}^{\prime}<\epsilon$ leads to $\delta{t}_{0}^{\prime}\approx\delta{t}_{0}(\epsilon/\eta_{F}^{(24)})^{1/3}$ as an optimal choice within our error tolerance. For a safety margin, we introduce a constant $C$ ($0<C<1$) and set $\delta{t}_{0}^{\prime}=C\delta{t}_{0}(\epsilon/\eta_{F}^{(24)})^{1/3}$ as an updated trial $\delta{t}_{0}$. We repeat this update procedure until $\eta_{F}^{(24)}$ gets smaller than $\epsilon$ and accept the latest $\delta{t}_{0}$ to evolve our state as $\ket{\psi_{2}}=T_{2}(\delta{t}_{0})\ket{\psi_{0}}$.

Next, we move on to the second step, using a time step $\delta{t}_{1}$. In choosing this, we again use the latest $\eta_{F}^{(24)}$ obtained at the end of the previous time step. Since $\ket{\psi_{2}}\approx\ket{\psi_{0}}$, we can expect the error scaling coefficient $\alpha$ to be almost the same in the present and previous steps. Therefore, like in the updated trials within the previous time step, we have $\delta{t}_{1}=C\delta{t}_{0}(\epsilon/\eta_{F}^{(24)})^{1/3}$ as a good candidate for the optimal stepsize in the present time step. We note that $\eta_{F}^{(24)}$ here is what was measured in the previous step, and we have not made any measurements in the present step yet. Using this $\delta{t}_{1}$ as a trial stepsize, we implement $T_{2}(\delta{t}_{1})$ and $T_{4}(\delta{t}_{1})$ and calculate $\eta_{F}^{(24)}$ using a quantum circuit. Depending on whether $\eta_{F}^{(24)}$ is less or greater than $\epsilon$, we accept or update $\delta{t}_{1}$ like in the previous step.

The following iteration is straightforward and repeated until the accumulated evolution time $t_{\mathrm{ini}}+\delta{t}_{0}+\delta{t}_{1}+\dots$ exceeds the final time $t_{\mathrm{fin}}$. We summarize a pseudocode for the algorithm in Algorithm III.

Fidelity-based Trotter24

Input: Initial and final times, $t_{\mathrm{ini}}$ and $t_{\mathrm{fin}}$, an initial state $\ket{\psi_{0}}$, a Hamiltonian $H=A+B$, an error tolerance $\epsilon$, an initial stepsize $\delta t_{0}$, a safety constant $C$ ($0<C<1$), an oracle function $\mathrm{FIDELITY}(\ket{\phi},\ket{\psi})$ that calculates $|\braket{\phi}{\psi}|^{2}$.

Output: An ordered list of unitaries $U_{\mathrm{list}}$ that approximates $e^{-iH(t_{\mathrm{fin}}-t_{\mathrm{ini}})}$ within the error tolerance for each time step.

1:$t\leftarrow t_{\mathrm{ini}}$

2:$\delta t\leftarrow\delta t_{0}$

3:$U_{\mathrm{list}}=\{\}$ (empty list)

4:while $t+\delta t<t_{\mathrm{fin}}$ do

5:     $\ket{\psi(t)}\leftarrow\prod_{k}(U_{\mathrm{list}})_{k}\ket{\psi_{0}}$

6:     do

7:         $T_{2}(\delta t)\leftarrow e^{-iA\delta t/2}e^{-iB\delta t}e^{-iA\delta t/2}$

8:         $T_{4}(\delta t)\leftarrow e^{-i\frac{s}{2}A\delta t}$ $e^{-isB\delta t}$ $e^{-i\frac{1-s}{2}A\delta t}$ $e^{-i(1-2s)B\delta t}$ $\times e^{-i\frac{1-s}{2}A\delta t}$ $e^{-isB\delta t}$ $e^{-i\frac{s}{2}A\delta t}$

9:         $\eta\leftarrow 1-\mathrm{FIDELITY}(T_{4}(\delta t)\ket{\psi(t)},T_{2}(\delta t)\ket{\psi(t)})$

10:         $\delta t\leftarrow C\cdot(\epsilon/\eta)^{1/3}\delta t$

11:     while $\eta>\epsilon$

12:     Prepend $T_{2}(\delta t)$ to the ordered list $U_{\mathrm{list}}$

13:     $t\leftarrow t+\delta t$ return $U_{\mathrm{list}}$

Let us make a parallel argument for the observable error $\eta_{O}$ instead of the fidelity error $\eta_{F}$. At each time step, we measure $\eta_{O}^{(24)}$ and judge if the condition

$$|\eta_{O}^{(24)}|<\epsilon_{O}\|O\|$$ (17)

is met. This is an analog of Eq. (16), and we introduced the operator norm $\|O\|$ as a reference scale and put the subscript $O$ on the tolerance as $\epsilon_{O}$ to avoid confusion. The iteration scheme is parallel to the fidelity case.We summarize a pseudocode for the observable-based algorithm in Algorithm III.

Observable-based Trotter24