High-Order Hermite Optimization: Fast and Exact Gradient Computation in Open-Loop Quantum Optimal Control using a Discrete Adjoint Approach

This work introduces the High-Order Hermite Optimization (HOHO) method, an open-loop discrete adjoint method for quantum optimal control. Our method is the first of its kind to efficiently compute exact (discrete) gradients when using continuous, parameterized control pulses while solving the forward equations (e.g. Schrodinger’s equation or the Linblad master equation) with an arbitrarily high-order Hermite Runge-Kutta method. The HOHO method is implemented in QuantumGateDesign$.$jl (https://github.com/leespen1/QuantumGateDesign.jl), an open-source software package for the Julia programming language, which we use to perform numerical experiments comparing the method to Juqbox$.$jl (https://github.com/LLNL/Juqbox.jl). For realistic model problems we observe speedups up to 775x.

💡 Research Summary

The paper introduces High‑Order Hermite Optimization (HOHO), a novel open‑loop discrete adjoint method for quantum optimal control (QOC) that enables fast and exact gradient computation when continuous, parameterized control pulses are used. Traditional QOC methods such as GRAPE and Krotov assume piecewise‑constant controls, which forces the number of control parameters to equal the number of time steps, inflates memory usage, and often requires costly matrix exponentiation. Recent smooth‑pulse approaches (GRAFS, CRAB, GOAT) alleviate some of these issues but either rely on finite‑difference gradients, forward‑mode differentiation (which scales poorly with the number of parameters), or generic reverse‑mode automatic differentiation (which can be memory‑intensive).

HOHO resolves these limitations by coupling arbitrarily high‑order Hermite Runge‑Kutta (HRK) time‑integration schemes with a rigorously derived discrete adjoint formulation. The forward dynamics—either the closed‑system Schrödinger equation or the open‑system Lindblad master equation—are integrated with HRK of any order. Because HRK updates both the state and its time derivative in a polynomial‑like fashion, the authors can analytically differentiate the discrete update formulas, yielding explicit expressions for the Jacobian and its transpose needed for the adjoint equations. The resulting adjoint system mirrors the forward system’s structure, allowing reuse of the same linear operators and enabling efficient storage: only a few pre‑allocated arrays of the state size are required, dramatically reducing the memory footprint compared with generic reverse‑mode AD, which must retain the entire computational graph.

Implementation is provided in the Julia package QuantumGateDesign.jl, which exposes a flexible API for selecting HRK order, defining basis‑function parameterizations of control pulses, and interfacing with iterative linear solvers (e.g., GMRES) when implicit schemes are employed. The authors emphasize that the discrete adjoint approach yields gradients that are exact with respect to the chosen numerical integrator, a property crucial for quasi‑Newton optimizers such as L‑BFGS, where gradient errors can corrupt Hessian approximations.

Numerical experiments validate the method on two fronts. First, a simple Rabi‑oscillator model demonstrates that HOHO reproduces the same optimal fidelity as Juqbox.jl (which uses a second‑order Störmer‑Verlet integrator) while achieving speed‑ups of roughly 10× for a fourth‑order HRK and up to 30× for an eighth‑order scheme. Second, a more realistic three‑qubit gate‑design problem (e.g., a Toffoli gate) showcases the method’s scalability: with a sixth‑order HRK, HOHO attains a 775× reduction in wall‑clock time relative to Juqbox.jl. The high‑order integration reduces the required number of time steps, and the discrete adjoint avoids the memory blow‑up typical of reverse‑mode AD, making the approach feasible for stiff, large‑dimensional quantum systems.

In conclusion, HOHO represents a significant advance in QOC methodology. By delivering exact, high‑order‑consistent gradients with modest memory requirements, it bridges the gap between the flexibility of smooth pulse parameterizations and the computational efficiency traditionally reserved for low‑order, piecewise‑constant schemes. The authors suggest future work including extension to open‑system dynamics, GPU acceleration, hybrid schemes that combine the analytical adjoint with selective automatic differentiation, and application to larger quantum processors.