GitHub

Composing Templates

Problem templates in Piccolo.jl are designed to be composable. You can chain them together to build sophisticated optimization pipelines that combine multiple capabilities.

Composition Overview

The templates form a hierarchy:

Base Problems (create from trajectory):
├── SmoothPulseProblem
├── BangBangPulseProblem
└── SplinePulseProblem
        │
        ▼
Wrapper Problems (wrap existing problem):
├── SamplingProblem (adds robustness)
└── MinimumTimeProblem (adds time optimization)

Any wrapper can wrap another wrapper, enabling combinations like:

  • MinimumTimeProblem(SamplingProblem(SmoothPulseProblem(...)))
  • MinimumTimeProblem(SplinePulseProblem(...))

Full Pipeline Example

Here's a complete pipeline: base optimization → robust optimization → time-optimal control.

using Piccolo

Step 1: Setup

# Define nominal system
H_drift = 0.5 * PAULIS[:Z]
H_drives = [PAULIS[:X], PAULIS[:Y]]
sys_nominal = QuantumSystem(H_drift, H_drives, [1.0, 1.0])

# Create initial trajectory
T, N = 20.0, 100
times = collect(range(0, T, length = N))
pulse = ZeroOrderPulse(0.1 * randn(2, N), times)
qtraj = UnitaryTrajectory(sys_nominal, pulse, GATES[:X])
UnitaryTrajectory{ZeroOrderPulse{DataInterpolations.ConstantInterpolation{Matrix{Float64}, Vector{Float64}, Vector{Union{}}, Float64}}, SciMLBase.ODESolution{ComplexF64, 3, Vector{Matrix{ComplexF64}}, Nothing, Nothing, Vector{Float64}, Vector{Vector{Matrix{ComplexF64}}}, Nothing, SciMLBase.ODEProblem{Matrix{ComplexF64}, Tuple{Float64, Float64}, true, SciMLBase.NullParameters, SciMLBase.ODEFunction{true, SciMLBase.AutoSpecialize, SciMLOperators.MatrixOperator{ComplexF64, Matrix{ComplexF64}, SciMLOperators.FilterKwargs{Nothing, Val{()}}, SciMLOperators.FilterKwargs{Piccolo.Quantum.Rollouts.var"#update!#_construct_operator##2"{QuantumSystem{Piccolo.Quantum.QuantumSystems.var"#53#54"{SparseArrays.SparseMatrixCSC{ComplexF64, Int64}, Vector{SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, Int64}, Piccolo.Quantum.QuantumSystems.var"#55#56"{Vector{SparseArrays.SparseMatrixCSC{Float64, Int64}}, Int64, SparseArrays.SparseMatrixCSC{Float64, Int64}}, @NamedTuple{}, Vector{DriftTerm{SparseArrays.SparseMatrixCSC{ComplexF64, Int64}, Returns{Float64}, Returns{Float64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}}, Val{()}}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Piccolo.Quantum.Rollouts.PiccoloRolloutSystem{Union{Int64, AbstractVector{Int64}, CartesianIndex, CartesianIndices}}, Nothing, Nothing}, Base.Pairs{Symbol, Vector{Float64}, Nothing, @NamedTuple{tstops::Vector{Float64}}}, SciMLBase.StandardODEProblem}, OrdinaryDiffEqLinear.MagnusAdapt4, OrdinaryDiffEqCore.InterpolationData{SciMLBase.ODEFunction{true, SciMLBase.AutoSpecialize, SciMLOperators.MatrixOperator{ComplexF64, Matrix{ComplexF64}, SciMLOperators.FilterKwargs{Nothing, Val{()}}, SciMLOperators.FilterKwargs{Piccolo.Quantum.Rollouts.var"#update!#_construct_operator##2"{QuantumSystem{Piccolo.Quantum.QuantumSystems.var"#53#54"{SparseArrays.SparseMatrixCSC{ComplexF64, Int64}, Vector{SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, Int64}, Piccolo.Quantum.QuantumSystems.var"#55#56"{Vector{SparseArrays.SparseMatrixCSC{Float64, Int64}}, Int64, SparseArrays.SparseMatrixCSC{Float64, Int64}}, @NamedTuple{}, Vector{DriftTerm{SparseArrays.SparseMatrixCSC{ComplexF64, Int64}, Returns{Float64}, Returns{Float64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}}, Val{()}}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Piccolo.Quantum.Rollouts.PiccoloRolloutSystem{Union{Int64, AbstractVector{Int64}, CartesianIndex, CartesianIndices}}, Nothing, Nothing}, Vector{Matrix{ComplexF64}}, Vector{Float64}, Vector{Vector{Matrix{ComplexF64}}}, Nothing, OrdinaryDiffEqLinear.MagnusAdapt4Cache{Matrix{ComplexF64}, Matrix{ComplexF64}, Matrix{ComplexF64}, Matrix{ComplexF64}, Nothing}, Nothing}, SciMLBase.DEStats, Nothing, Nothing, Nothing, Nothing}}(QuantumSystem: levels = 2, n_drives = 2, ZeroOrderPulse(Number of drives = 2, T = 20.0), ComplexF64[1.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 1.0 + 0.0im], ComplexF64[0.0 + 0.0im 1.0 + 0.0im; 1.0 + 0.0im 0.0 + 0.0im], ComplexF64[1.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 1.0 + 0.0im;;; 0.9949908641865791 - 0.09983297298085714im 0.0013991200978624883 - 0.004959854231740058im; -0.001399120097862489 - 0.00495985423174006im 0.9949908641865789 + 0.09983297298085711im;;; 0.9798867132626343 - 0.1987377280285841im 0.018018704700881 - 0.0008190911371980083im; -0.018018704700881003 - 0.00081909113719801im 0.9798867132626345 + 0.19873772802858414im;;; … ;;; -0.8525249543266193 + 0.34934637043930183im 0.022308781537524667 + 0.388150272417244im; -0.022308781537525 + 0.38815027241724437im -0.8525249543266186 - 0.34934637043930167im;;; -0.8078207059964332 + 0.4243998944186706im 0.07518886305747637 + 0.402066003851305im; -0.07518886305747677 + 0.40206600385130553im -0.8078207059964327 - 0.42439989441867065im;;; -0.7687780993192775 + 0.4990666877845331im 0.12790691320550931 + 0.3788832230450973im; -0.12790691320550981 + 0.3788832230450978im -0.7687780993192772 - 0.4990666877845332im])

Step 2: Base Problem (with free time enabled)

qcp_base = SmoothPulseProblem(
    qtraj,
    N;
    Q = 100.0,
    R = 1e-2,
    Δt_bounds = (0.05, 0.5),  ## Required for MinimumTimeProblem
)
solve!(qcp_base; max_iter = 100)

fidelity(qcp_base)
0.989797705664722
sum(get_timesteps(get_trajectory(qcp_base)))
21.966699602348307

Step 3: Add Robustness

# Create perturbed systems (±5% drift variation)
sys_high = QuantumSystem(1.05 * H_drift, H_drives, [1.0, 1.0])
sys_low = QuantumSystem(0.95 * H_drift, H_drives, [1.0, 1.0])

qcp_robust = SamplingProblem(qcp_base, [sys_nominal, sys_high, sys_low]; Q = 100.0)
solve!(qcp_robust; max_iter = 100)

fidelity(qcp_robust)
3-element Vector{Float64}:
 0.9999794169785005
 0.9999917626760731
 0.9999924228552589

Step 4: Minimize Time

qcp_mintime = MinimumTimeProblem(qcp_robust; final_fidelity = 0.95, D = 100.0)
solve!(qcp_mintime; max_iter = 100)

fidelity(qcp_mintime)
3-element Vector{Float64}:
 0.9997574757773747
 0.9973204700506145
 0.9953739450195523
sum(get_timesteps(get_trajectory(qcp_mintime)))
13.125869584312529

Common Composition Patterns

Pattern 1: Robust Gate

Optimize for parameter uncertainty without time constraints.

qcp_base = SmoothPulseProblem(qtraj, N; Q=100.0)
solve!(qcp_base)

qcp_robust = SamplingProblem(qcp_base, systems)
solve!(qcp_robust)

Pattern 2: Fast Gate

Minimize time without robustness requirements.

qcp_base = SmoothPulseProblem(qtraj, N; Q=100.0, Δt_bounds=(0.01, 0.5))
solve!(qcp_base)

qcp_fast = MinimumTimeProblem(qcp_base; final_fidelity=0.99)
solve!(qcp_fast)

Pattern 3: Fast + Robust Gate

The full pipeline for production-quality gates.

qcp_base = SmoothPulseProblem(qtraj, N; Q=100.0, Δt_bounds=(0.01, 0.5))
solve!(qcp_base)

qcp_robust = SamplingProblem(qcp_base, systems)
solve!(qcp_robust)

qcp_final = MinimumTimeProblem(qcp_robust; final_fidelity=0.95)
solve!(qcp_final)

Pattern 4: Spline Warm-Start Pipeline

Start with smooth problem, refine with splines.

# Initial optimization with piecewise constant
qcp_smooth = SmoothPulseProblem(qtraj_smooth, N; Q=100.0)
solve!(qcp_smooth)

# Extract optimized pulse and convert to spline
optimized_pulse = get_pulse(qcp_smooth.qtraj)
spline_pulse = CubicSplinePulse(optimized_pulse)
qtraj_spline = UnitaryTrajectory(sys, spline_pulse, U_goal)

# Refine with spline problem
qcp_spline = SplinePulseProblem(qtraj_spline; Q=100.0)
solve!(qcp_spline; max_iter=50)  # Quick refinement

Iteration and Refinement

You can iteratively refine solutions:

# First pass: coarse optimization
qcp = SmoothPulseProblem(qtraj, 50; Q=10.0)
solve!(qcp; max_iter=50)

# Second pass: increase resolution
# ... resample to higher N ...

# Third pass: tighten tolerances
solve!(qcp; max_iter=100, tol=1e-8)

Accessing Results Through the Chain

Each wrapper preserves access to the underlying trajectory:

# Access trajectory at any level
traj = get_trajectory(qcp_mintime)
sys = get_system(qcp_mintime)

# Fidelity evaluation uses the innermost trajectory type
fid = fidelity(qcp_mintime)

# Get optimized pulse
pulse = get_pulse(qcp_mintime.qtraj)

Order Matters

The order of composition affects the optimization:

MinimumTimeProblem(SamplingProblem(base)):

  • First achieves robustness, then minimizes time
  • Time minimization respects the robust solution
  • Generally preferred for production gates

SamplingProblem(MinimumTimeProblem(base)):

  • First minimizes time, then adds robustness
  • May require re-solving if time-optimal solution isn't robust
  • Less common, but useful for exploring trade-offs

Tips for Complex Pipelines

1. Solve Each Stage

Always solve! after each composition step:

qcp_base = SmoothPulseProblem(qtraj, N)
solve!(qcp_base)  # Important!

qcp_robust = SamplingProblem(qcp_base, systems)
solve!(qcp_robust)  # Important!

qcp_mintime = MinimumTimeProblem(qcp_robust)
solve!(qcp_mintime)  # Important!

2. Monitor Progress

Check fidelity at each stage to catch issues early:

for (name, prob) in [("Base", qcp_base), ("Robust", qcp_robust), ("MinTime", qcp_mintime)]
    println("$name: Fidelity = $(fidelity(prob))")
end
Base: Fidelity = 0.9997574757773747
Robust: Fidelity = [0.9997574757773747, 0.9973204700506145, 0.9953739450195523]
MinTime: Fidelity = [0.9997574757773747, 0.9973204700506145, 0.9953739450195523]

3. Adjust Parameters at Each Stage

Different stages may need different settings:

# Base: prioritize finding a solution
qcp_base = SmoothPulseProblem(qtraj, N; Q=100.0, R=1e-2)
solve!(qcp_base; max_iter=200)

# Robust: may need more iterations
qcp_robust = SamplingProblem(qcp_base, systems; Q=100.0)
solve!(qcp_robust; max_iter=300)

# MinTime: typically faster since starting from good solution
qcp_mintime = MinimumTimeProblem(qcp_robust; final_fidelity=0.95, D=100.0)
solve!(qcp_mintime; max_iter=100)

See Also

This page was generated using Literate.jl.