API
The stable public API of the GRAPE
consists of the following members:
GRAPE.optimize
as the main function to run an optimizationGRAPE.GrapeResult
as the object returned byGRAPE.optimize
, and accessible in callbacksQuantumControl.optimize
withmethod=GRAPE
, as a higher-level wrapper aroundGRAPE.optimize
with extra featuresGRAPE.Trajectory
as an alias ofQuantumControl.Trajectory
GRAPE.set_default_ad_framework
as an alias ofQuantumControl.set_default_ad_framework
The remaining functions in GRAPE
documented below should not be considered part of the stable API. They are guaranteed to be stable in bugfix (x.y.z
) releases, but may change in feature releases (x.y
).
Note that the GRAPE
package does not export any symbols. All members of the public API must be explicitly imported or used with their fully qualified name.
Index
GRAPE.GrapeResult
GRAPE.GrapeWrk
QuantumControl.Trajectory
GRAPE.evaluate_functional
GRAPE.evaluate_gradient!
GRAPE.gradient
GRAPE.make_grape_print_iters
GRAPE.norm_search
GRAPE.optimize
GRAPE.pulse_update
GRAPE.search_direction
GRAPE.step_width
GRAPE.vec_angle
QuantumControl.optimize
QuantumControl.set_default_ad_framework
Reference
GRAPE.GrapeResult
— TypeResult object returned by GRAPE.optimize
.
Attributes
The attributes of a GrapeResult
object include
iter
: The number of the current iterationJ_T
: The value of the final-time functional in the current iterationJ_T_prev
: The value of the final-time functional in the previous iterationJ_a
: The value of the running cost $J_a$ in the current iteration (excluding $λ_a$)J_a_prev
: The value of $J_a$ in the previous iterationtlist
: The time grid on which the control are discetized.guess_controls
: A vector of the original control fields (each field discretized to the points oftlist
)optimized_controls
: A vector of the optimized control fields in the current iterationsrecords
: A vector of tuples with values returned by acallback
routine passed tooptimize
converged
: A boolean flag on whether the optimization is converged. This may be set totrue
by acheck_convergence
function.message
: A message string to explain the reason for convergence. This may be set by acheck_convergence
function.
All of the above attributes may be referenced in a check_convergence
function passed to QuantumControl.optimize(problem; method=GRAPE)
or GRAPE.optimize
.
GRAPE.GrapeWrk
— TypeGRAPE Workspace.
The workspace is for internal use. However, it is also accessible in a callback
function. The callback may use or modify some of the following attributes:
trajectories
: a copy of the trajectories defining the control problemtlist
: the time grid for the optimizationadjoint_trajectories
: Thetrajectories
with the adjoint generatorkwargs
: The keyword arguments from the call tooptimize
.controls
: A tuple of the original controls (probably functions)pulsevals_guess
: The combined vector of pulse values that are the guess in the current iteration. Initially, the vector is the concatenation of discretizingcontrols
to the midpoints of the time grid.pulsevals
: The combined vector of updated pulse values in the current iteration. All the initialized propagators inside the workspace aliaspulsevals
such that mutatingpulsevals
is directly reflected in the next propagation step.gradient
: The total gradient for the guess in the current iterationgrad_J_T
: The current gradient for the final-time part of the functional. This is from the last evaluation of the gradient, which may be for the optimized pulse (depending on the internal of the optimizer)grad_J_a
: The current gradient for the running cost part of the functional.J_parts
: The two-component vector $[J_T, J_a]$upper_bounds
: Upper bound for everypulsevals
;+Inf
indicates no bound.lower_bounds
: Lower bound for everypulsevals
;-Inf
indicates no bound.fg_count
: A two-element vector containing the number of evaluations of the combined gradient and functional first, and the evaluations of only the functional second.optimizer
: The backend optimizer objectoptimizer_state
: The internal state object of theoptimizer
(nothing
if theoptimizer
has no internal state)result
: The current result objecttau_grads
: The gradients ∂τₖ/ϵₗ(tₙ)fw_storage
: The storage of states for the forward propagation, as a vector of storage contains (one for each trajectory)fw_propagators
: The propagators used for the forward propagationbw_grad_propagators
: The propagators used for the backward propagation ofQuantumGradientGenerators.GradVector
states (gradient_method=:gradgen
only)bw_propagators
: The propagators used for the backward propagation (gradient_method=:taylor
only)use_threads
: Flag indicating whether the propagations are performed in parallel.
In addition, the following methods provide safer (non-mutating) access to information in the workspace
GRAPE.evaluate_functional
— MethodEvaluate the optimization functional encoded in wrk
for the given pulsevals
.
J = evaluate_functional(pulsevals, wrk; storage=nothing, count_call=true)
evaluates the functional defined during the initialization of the GRAPE workspace wrk
, for the given pulse values, using wrk.fw_propagators
. The pulsevals
argument is a vector of Float64
values corresponding to a concatenation of all the controls, discretized to the midpoints of the time grid, cf. GrapeWrk
.
As a side effect, the evaluation sets the following information in wrk
:
wrk.pulsevals
: On output, the values of the givenpulsevals
. Note thatpulsevals
may aliaswrk.pulsevals
, so there is no assumption made onwrk.pulsevals
other than that mutatingwrk.pulsevals
directly affects the propagators inwrk
.wrk.result.f_calls
: Will be incremented by one (only ifcount_call=true
)wrk.fg_count[2]
: Will be incremented by one (only ifcount_call=true
)wrk.result.tau_vals
: For any trajectory that defines atarget_state
, the overlap of the propagated state with that target state.wrk.J_parts
: The parts (J_T
,λₐJ_a
) of the functional
If storage
is given, as a vector of storage containers suitable for propagate
(one for each trajectory), the forward-propagated states will be stored there.
Returns J
as sum(wrk.J_parts)
.
GRAPE.evaluate_gradient!
— MethodEvaluate the gradient $∂J/∂ϵₙₗ$ into G
, together with the functional J
.
J = evaluate_gradient!(G, pulsevals, wrk)
evaluates and returns the optimization functional defined during the initialization of wrk
, for the given pulse values, cf. evaluate_functional
, and write the derivative of the optimization functional with respect to the pulse values into the existing array G
.
The evaluation of the functional uses uses wrk.fw_propagators
. The evaluation of the gradient happens either via a backward propagation of an extended "gradient vector" using wrk.bw_grad_propagators
if wrk
was initialized with gradient_method=:gradgen
. Alternatively, if wrk
was initialized with gradient_method=:taylor
, the backward propagation if for a regular state, using wrk.bw_propagators
, and a Taylor expansion is used for the gradient of the time evolution operator in a single time step.
As a side, effect, evaluating the gradient and functional sets the following information in wrk
:
wrk.pulsevals
: On output, the values of the givenpulsevals
, seeevaluate_functional
.wrk.result.fg_calls
: Will be incremented by onewrk.fg_count[1]
: Will be incremented by onewrk.result.tau_vals
: For any trajectory that defines atarget_state
, the overlap of the propagated state with that target state.wrk.J_parts
: The parts (J_T
,λₐJ_a
) of the functionalwrk.fw_storage
: For each trajectory, the forward-propagated states at each point on the time grid.wrk.chi_states
: The normalized states $|χ(T)⟩$ that we used as the boundary condition for the backward propagation.wrk.chi_states_norm
: The original norm of the states $|χ(T)⟩$, as calculated by $-∂J/∂⟨Ψₖ|$wrk.grad_J_T
: The vector ``∂JT/∂ϵ{nl}, i.e., the gradient only for the final-time part of the functionalwrk.grad_J_a
: The vector $∂J_a/∂ϵ_{nl}$, i.e., the gradient only for the pulse-dependent running cost.
The gradients are wrk.grad_J_T
and wrk.grad_J_a
(weighted by $λ_a$) into are combined into the output G
.
Returns the value of the functional.
GRAPE.gradient
— MethodThe gradient in the current iteration.
g = gradient(wrk; which=:initial)
returns the gradient associated with the guess pulse of the current iteration. Up to quasi-Newton corrections, the negative gradient determines the search_direction
for the pulse_update
.
g = gradient(wrk; which=:final)
returns the gradient associated with the optimized pulse of the current iteration.
GRAPE.make_grape_print_iters
— MethodPrint optimization progress as a table.
print_iters = make_grape_print_iters(; print_iter_info, store_iter_info=[])
generates a print_iters
function that can be passed as callback
to GRAPE.optimize
. It is also used automatically when GRAPE.optimized
is called via QuantumControl.optimize
with print_iters=true
.
The print_iter_info
keyword argument specifies what information should be printed, and defaults to ["iter.", "J_T", "|∇J|", "|Δϵ|", "ΔJ", "FG(F)", "secs"]
. The store_iter_info
similarly specifies what information should be returned from the callback, so that it can be stored in the records
field of the GrapeResult
object.
The available fields for print_iter_info
and store_iter_info
are:
"iter."
: The iteration number"J_T"
: The value of the final-time functional for the dynamics under the optimized pulses"J_a"
: The value of the pulse-dependent running cost for the optimized pulses"λ_a⋅J_a"
: The total contribution ofJ_a
to the full functionalJ
"J"
: The value of the optimization functional for the optimized pulses"ǁ∇J_Tǁ"
: The ℓ²-norm of the current gradient of the final-time functional. Note that this is usually the gradient of the optimize pulse, not the guess pulse."ǁ∇J_aǁ"
: The ℓ²-norm of the the current gradient of the pulse-dependent running cost. For comparison with"ǁ∇J_Tǁ"
."λ_aǁ∇J_aǁ"
: The ℓ²-norm of the the current gradient of the complete pulse-dependent running cost term. For comparison with"ǁ∇J_Tǁ"
."ǁ∇Jǁ"
: The norm of the guess pulse gradient. Note that the guess pulse gradient is not the same the current gradient."ǁΔϵǁ"
: The ℓ²-norm of the pulse update"ǁϵǁ"
: The ℓ²-norm of optimized pulse values"max|Δϵ|"
The maximum value of the pulse update (infinity norm)"max|ϵ|"
: The maximum value of the pulse values (infinity norm)"ǁΔϵǁ/ǁϵǁ"
: The ratio of the pulse update tothe optimized pulse values"∫Δϵ²dt"
: The L²-norm of the pulse update, summed over all pulses. A convergence measure comparable (proportional) to the running cost in Krotov's method"ǁsǁ"
: The norm of the search direction. Should beǁΔϵǁ
scaled by the step withα
."∠°"
: The angle (in degrees) between the negative gradient-∇J
and the search directions
."α"
: The step width as determined by the line search (Δϵ = α⋅s
)"ΔJ_T"
: The change in the final time functional relative to the previous iteration"ΔJ_a"
: The change in the control-dependent running cost relative to the previous iteration"λ_a⋅ΔJ_a"
: The change in the control-dependent running cost term relative to the previous iteration."ΔJ"
: The change in the total optimization functional relative to the previous iteration."FG(F)"
: The number of functional/gradient evaluation (FG), or pure functional (F) evaluations"secs"
: The number of seconds of wallclock time spent on the iteration.
GRAPE.norm_search
— MethodThe norm of the search direction vector in the current iteration.
norm_search(wrk)
returns norm(search_direction(wrk))
.
GRAPE.optimize
— FunctionSolve a quantum control problem using the GRAPE method.
using GRAPE
result = GRAPE.optimize(trajectories, tlist; J_T, kwargs...)
minimizes a functional
\[J(\{ϵ_{nl}\}) = J_T(\{|Ψ_k(T)⟩\}) + λ_a J_a(\{ϵ_{nl}\})\,,\]
via the GRAPE method, where the final time functional $J_T$ depends explicitly on the forward-propagated states $|Ψ_k(T)⟩$, where $|Ψ_k(t)⟩$ is the time evolution of the initial_state
in the $k$th' element of the trajectories
, and the running cost $J_a$ depends explicitly on pulse values $ϵ_{nl}$ of the l'th control discretized on the n'th interval of the time grid tlist
.
It does this by calculating the gradient of the final-time functional
\[\nabla J_T \equiv \frac{\partial J_T}{\partial ϵ_{nl}} = -2 \Re \underbrace{% \underbrace{\bigg\langle χ(T) \bigg\vert \hat{U}^{(k)}_{N_T} \dots \hat{U}^{(k)}_{n+1} \bigg \vert}_{\equiv \bra{\chi(t_n)}\;\text{(bw. prop.)}} \frac{\partial \hat{U}^{(k)}_n}{\partial ϵ_{nl}} }_{\equiv \bra{χ_k^\prime(t_{n-1})}} \underbrace{\bigg \vert \hat{U}^{(k)}_{n-1} \dots \hat{U}^{(k)}_1 \bigg\vert Ψ_k(t=0) \bigg\rangle}_{\equiv |\Psi(t_{n-1})⟩\;\text{(fw. prop.)}}\,,\]
where $\hat{U}^{(k)}_n$ is the time evolution operator for the $n$ the interval, generally assumed to be $\hat{U}^{(k)}_n = \exp[-i \hat{H}_{kn} dt_n]$, where $\hat{H}_{kn}$ is the operator obtained by evaluating trajectories[k].generator
on the $n$'th time interval.
The backward-propagation of $|\chi_k(t)⟩$ has the boundary condition
\[ |\chi_k(T)⟩ \equiv - \frac{\partial J_T}{\partial ⟨\Psi_k(T)|}\,.\]
The final-time gradient $\nabla J_T$ is combined with the gradient for the running costs, and the total gradient is then fed into an optimizer (L-BFGS-B by default) that iteratively changes the values $\{ϵ_{nl}\}$ to minimize $J$.
See Background for details.
Returns a GrapeResult
.
Positional arguments
trajectories
: A vector ofTrajectory
objects. Each trajectory contains aninitial_state
and a dynamicalgenerator
(e.g., time-dependent Hamiltonian). Each trajectory may also contain arbitrary additional attributes liketarget_state
to be used in theJ_T
functionaltlist
: A vector of time grid values.
Required keyword arguments
J_T
: A functionJ_T(Ψ, trajectories)
that evaluates the final time functional from a listΨ
of forward-propagated states andtrajectories
. The functionJ_T
may also take a keyword argumenttau
. If it does, a vector containing the complex overlaps of the target states (target_state
property of each trajectory intrajectories
) with the propagated states will be passed toJ_T
.
Optional keyword arguments
chi
: A functionchi(Ψ, trajectories)
that receives a listΨ
of the forward propagated states and returns a vector of states $|χₖ⟩ = -∂J_T/∂⟨Ψₖ|$. If not given, it will be automatically determined fromJ_T
viaQuantumControl.Functionals.make_chi
with the default parameters. Similarly toJ_T
, ifchi
accepts a keyword argumenttau
, it will be passed a vector of complex overlaps.chi_min_norm=1e-100
: The minimum allowable norm for any $|χₖ(T)⟩$. Smaller norms would mean that the gradient is zero, and will abort the optimization with an error.J_a
: A functionJ_a(pulsevals, tlist)
that evaluates running costs over the pulse values, wherepulsevals
are the vectorized values $ϵ_{nl}$, wheren
are in indices of the time intervals andl
are the indices over the controls, i.e.,[ϵ₁₁, ϵ₂₁, …, ϵ₁₂, ϵ₂₂, …]
(the pulse values for each control are contiguous). If not given, the optimization will not include a running cost.gradient_method=:gradgen
: One of:gradgen
(default) or:taylor
. Withgradient_method=:gradgen
, the gradient is calculated using QuantumGradientGenerators. Withgradient_method=:taylor
, it is evaluated via a Taylor series, see Eq. (20) in Kuprov and Rogers, J. Chem. Phys. 131, 234108 (2009) [22].taylor_grad_max_order=100
: If given withgradient_method=:taylor
, the maximum number of terms in the Taylor series. Iftaylor_grad_check_convergence=true
(default), if the Taylor series does not convergence within the given number of terms, throw an an error. Withtaylor_grad_check_convergence=true
, this is the exact order of the Taylor series.taylor_grad_tolerance=1e-16
: If given withgradient_method=:taylor
andtaylor_grad_check_convergence=true
, stop the Taylor series when the norm of the term falls below the given tolerance. Ignored iftaylor_grad_check_convergence=false
.taylor_grad_check_convergence=true
: If given astrue
(default), check the convergence after each term in the Taylor series an stop as soon as the norm of the term drops below the given number. Iffalse
, stop after exactlytaylor_grad_max_order
terms.lambda_a=1
: A weight for the running costJ_a
.grad_J_a
: A function to calculate the gradient ofJ_a
. If not given, it will be automatically determined. Seemake_grad_J_a
for the required interface.upper_bound
: An upper bound for the value of any optimized control. Time-dependent upper bounds can be specified viapulse_options
.lower_bound
: A lower bound for the value of any optimized control. Time-dependent lower bounds can be specified viapulse_options
.pulse_options
: A dictionary that maps every control (as obtained byget_controls
from thetrajectories
) to a dict with the following possible keys::upper_bounds
: A vector of upper bound values, one for each intervals of the time grid. Values ofInf
indicate an unconstrained upper bound for that time interval, respectively the globalupper_bound
, if given.:lower_bounds
: A vector of lower bound values. Values of-Inf
indicate an unconstrained lower bound for that time interval,
callback
: A function that receives the GRAPE workspace and the iteration number. The function may return a tuple of values which are stored in theGrapeResult
objectresult.records
. The function can also mutate the workspace, in particular the updatedpulsevals
. This may be used, e.g., to apply a spectral filter to the updated pulses or to perform similar manipulations.check_convergence
: A function to check whether convergence has been reached. Receives aGrapeResult
objectresult
, and must return one of the following:- A boolean (
true
if convergence is reached,false
otherwise) - A string with a reason for the convergence, or an empty string if not converged.
- The original
result
object ornothing
, indicating thatresult.converged
andresult.message
may have been modified to indicate convergence
callback
.- A boolean (
prop_method
: The propagation method to use for each trajectory, see below.verbose=false
: Iftrue
, print information during initializationrethrow_exceptions
: By default, any exception ends the optimization, but still returns aGrapeResult
that captures the message associated with the exception. This is to avoid losing results from a long-running optimization when an exception occurs in a later iteration. Ifrethrow_exceptions=true
, instead of capturing the exception, it will be thrown normally.
Experimental keyword arguments
The following keyword arguments may change in non-breaking releases:
x_tol
: Parameter for Optim.jlf_tol
: Parameter for Optim.jlg_tol
: Parameter for Optim.jlshow_trace
: Parameter for Optim.jlextended_trace
: Parameter for Optim.jlshow_every
: Parameter for Optim.jlallow_f_increases
: Parameter for Optim.jloptimizer
: An optional Optim.jl optimizer (Optim.AbstractOptimizer
instance). If not given, an L-BFGS-B optimizer will be used.
Trajectory propagation
GRAPE may involve three types of time propagation, all of which are implemented via the QuantumPropagators
as a numerical backend:
- A forward propagation for every
Trajectory
in thetrajectories
- A backward propagation for every trajectory
- A backward propagation of a gradient generator for every trajectory.
The keyword arguments for each propagation (see propagate
) are determined from any properties of each Trajectory
that have a prop_
prefix, cf. init_prop_trajectory
.
In situations where different parameters are required for the forward and backward propagation, instead of the prop_
prefix, the fw_prop_
and bw_prop_
prefix can be used, respectively. These override any setting with the prop_
prefix. Similarly, properties for the backward propagation of the gradient generators can be set with properties that have a grad_prop_
prefix. These prefixes apply both to the properties of each Trajectory
and the keyword arguments.
Note that the propagation method for each propagation must be specified. In most cases, it is sufficient (and recommended) to pass a global prop_method
keyword argument.
GRAPE.pulse_update
— MethodThe vector of pulse update values for the current iteration.
Δu = pulse_update(wrk)
returns a vector containing the different between the optimized pulse values and the guess pulse values of the current iteration. This should be proportional to search_direction
with the proportionality factor step_width
.
GRAPE.search_direction
— MethodThe search direction used in the current iteration.
s = search_direction(wrk)
returns the vector describing the search direction used in the current iteration. This should be proportional to pulse_update
with the proportionality factor step_width
.
GRAPE.step_width
— MethodThe step width used in the current iteration.
α = step_width(wrk)
returns the scalar α
so that pulse_update(wrk) = α * search_direction(wrk)
, see pulse_update
and search_direction
for the iteration described by the current GrapeWrk
(for the state of wrk
as available in the callback
of the current iteration.
GRAPE.vec_angle
— MethodThe angle between two vectors.
ϕ = vec_angle(v1, v2; unit=:rad)
returns the angle between two vectors in radians (or degrees, with unit=:degree
).
QuantumControl.optimize
— Methodusing GRAPE
result = optimize(problem; method=GRAPE, kwargs...)
optimizes the given QuantumControl.ControlProblem
using the GRAPE (Gradient-Ascent Pulse Engineering) method.
Delegates to
result = GRAPE.optimize(
problem.trajectories, problem.tlist; problem.kwargs..., kwargs...
)
See GRAPE.optimize
for details and supported keyword arguments.
Compared to calling GRAPE.optimize
directly, the QuantumControl.optimize
wrapper adds the following additional keyword arguments:
check=true
: Iftrue
(default), test that all the objects stored in the trajectories implement the required interfaces correctlyprint_iters=true
: Whether to print information after each iteration.print_iter_info=["iter.", "J_T", "|∇J|", "|Δϵ|", "ΔJ", "FG(F)", "secs"]
: Which fields to print ifprint_iters=true
. Seemake_grape_print_iters
store_iter_info=[]
: Which fields to store inresult.records
, given as a list of header labels, seeprint_iter_info
. Seemake_grape_print_iters
These options still allow for the normal callback
argument. With QuantumcControl.optimize
, the callback
can be a tuple of callback functions that will be combined automatically, which GRAPE.optimize
only supports as single callback function.
The GRAPE optimization may also be initiated via QuantumControl.@optimize_or_load
, which additionally adds checkpointing, to ensure that an optimization result is dumped to disk in case of an unexpected shutdown.
QuantumControl.Trajectory
— TypeDescription of a state's time evolution.
Trajectory(
initial_state,
generator;
target_state=nothing,
weight=1.0,
kwargs...
)
describes the time evolution of the initial_state
under a time-dependent dynamical generator
(e.g., a Hamiltonian or Liouvillian).
Trajectories are central to quantum control problems: an optimization functional depends on the result of propagating one or more trajectories. For example, when optimizing for a quantum gate, the optimization considers the trajectories of all logical basis states.
In addition to the initial_state
and generator
, a Trajectory
may include data relevant to the propagation and to evaluating a particular optimization functional. Most functionals have the notion of a "target state" that the initial_state
should evolve towards, which can be given as the target_state
keyword argument. In some functionals, different trajectories enter with different weights [38], which can be given as a weight
keyword argument. Any other keyword arguments are also available to a functional as properties of the Trajectory
.
A Trajectory
can also be instantiated using all keyword arguments.
Properties
All keyword arguments used in the instantiation are available as properties of the Trajectory
. At a minimum, this includes initial_state
, generator
, target_state
, and weight
.
By convention, properties with a prop_
prefix, e.g., prop_method
, will be taken into account when propagating the trajectory. See propagate_trajectory
for details.
QuantumControl.set_default_ad_framework
— FunctionSet the default provider for automatic differentiation.
QuantumControl.set_default_ad_framework(mod; quiet=false)
registers the given module (package) as the default AD framework.
This determines the default setting for the automatic
parameter in the following functions:
QuantumControl.Functionals.make_chi
QuantumControl.Functionals.make_gate_chi
QuantumControl.Functionals.make_grad_J_a
The given mod
must be a supported AD framework, e.g.,
import Zygote
QuantumControl.set_default_ad_framework(Zygote)
Currently, there is built-in support for Zygote
and FiniteDifferences
.
For other packages to be used as the default AD framework, the appropriate methods for make_chi
etc. must be defined.
Unless quiet=true
, calling set_default_ad_framework
will show a message to confirm the setting.
To unset the default AD framework, use
QuantumControl.set_default_ad_framework(nothing)