@@ -1122,7 +1122,7 @@ indicates that the terms at the endpoints ``n=0`` and ``n=N_T`` are halved.
11221122
11231123Note that `g_b` is a mandatory keyword argument and must be a function
11241124`g_b(Ψ, trajectory, tlist, n) -> Float64` where `Ψ` is the state ``|Ψ(t)⟩``
1125- at ``t`` corresponding to the the (one-based) time grid point `tlist[n]`, and
1125+ at ``t`` corresponding to the (one-based) time grid point `tlist[n]`, and
11261126`trajectory` is a [`QuantumControl.Trajectory`](@ref) that may hold additional
11271127data in a custom property that is relevant to the calculation of ``g_b``.
11281128
@@ -1131,19 +1131,20 @@ The definition of `J_b` here is compatible with the
11311131"""
11321132function J_b (storage, trajectories, tlist; g_b)
11331133 N = length (trajectories)
1134+ N_T = length (tlist)
11341135 result = 0.0
11351136 for k = 1 : N
11361137 Ψ₁ = get_from_storage (storage[k], 1 )
11371138 dt = tlist[2 ] - tlist[1 ]
11381139 result += g_b (Ψ₁, trajectories[k], tlist, 1 ) * (dt/ 2 )
1139- for n_tl = 2 : length (tlist)
1140- Ψₙ = get_from_storage (storage[k], n_tl )
1141- if n_tl < length (tlist)
1142- dt = 0.5 * (tlist[n_tl + 1 ] - tlist[n_tl - 1 ])
1143- result += g_b (Ψₙ, trajectories[k], tlist, n_tl ) * dt
1140+ for n = 2 : N_T
1141+ Ψₙ = get_from_storage (storage[k], n )
1142+ if n < N_T
1143+ dt = 0.5 * (tlist[n + 1 ] - tlist[n - 1 ])
1144+ result += g_b (Ψₙ, trajectories[k], tlist, n ) * dt
11441145 else
11451146 dt = tlist[end ] - tlist[end - 1 ]
1146- result += g_b (Ψₙ, trajectories[k], tlist, n_tl ) * (dt/ 2 )
1147+ result += g_b (Ψₙ, trajectories[k], tlist, n ) * (dt/ 2 )
11471148 end
11481149 end
11491150 end
0 commit comments