Adding constraint terms to the cost functional for Krotov's method

[HuelsMat]

Using Krotov's method, I'd like to add constraint terms $\int_0^T dt \left( g_a[{ \epsilon_k(t) } ] \right)$ and $\int_0^T dt \left( g_b[{ |\psi_j(t) \rangle } ] \right)$ to the cost functional J. Thereby, $g_a$ is supposed to penalize pulse forms that violate amplitude or frequency constraints. $g_b$ should penalize a population transfer to undesired states during the time evolution.

I was thinking about defining a cost functional J_T(Ψ, trajectories), where I compute the time evolution with the information given in trajectories and then use it to compute my constraint terms. However, this would be time consuming as most probably the time evolution would be calculated twice in each iteration.

Is there another way to implement these constraint terms?

When running an optimization with Krotov's method, we can track the value of $\int g_a (t) dt$. Here I assume that g_a quantifies the distance of the current pulse ϵ⁽ⁱ⁺¹⁾ to the pulse of the previous iteration ϵ⁽ⁱ⁾.
Is it possible, to directly define a custom g_a(ϵ⁽ⁱ⁺¹⁾) that takes the current pulse ϵ⁽ⁱ⁺¹⁾ as input?

[goerz]

Using Krotov's method, I'd like to add constraint terms $\int_0^T dt \left( g_a[{ \epsilon_k(t) } ] \right)$

Krotov's method does not support arbitrary running costs. Fundamentally, you should consider it a method suitable for unconstrained optimization. Some forms of g_a are possible, see, e.g., Eq. (3.63) in my thesis. This feature is not currently implemented, and would be very low priority (I have not ever found it to be useful).

Also, penalties do not enforce constraints. To ensure amplitude constraints with Krotov's method, "pulse parameterization" should be used, see Pulse Parametrization for Krotov's Method. However, since this makes the control non-linear, it is important to take small update steps (large λₐ), slowing down the conversion significantly. Non-linear parameterization can also introduce local traps in the optimization landscape, so you have to try different guess pulses.

Frequency constraints are mathematically possible, see J. P. Palao, D. Reich, and C. Koch, "Steering the optimization pathway in the control landscape using constraints", Phys. Rev. A 88, 053409 (2013). However, this significantly complicates the implementation, and there are no plans to implement this. In practice, a better way to impose spectral constraints is to apply a filter (via a callback function) to the optimized pulses in each iteration. Again, larger values of λₐ are recommended in this case to ensure that the filter does not have to modify the pulses too much.

State-dependent running costs g_b are also mathematically possible, but again complicate the method (although not quite as much as spectral constraints), see the paper linked above. Specifically, they require an inhomogenous backward propagation. I'm not planning to implement this in the near future, but I'm open to PRs (it would be a very big project, though). Some state-dependent running costs can be mapped onto the final-time functional J_T, and that's usually the way to go. For example, in order to avoid population in a given level, remove or decay population in that level during the propagation. The resulting state will have a norm < 1, which will increase the value of most J_T functionals. Simply by minimizing J_T, Krotov's method will effectively search for solutions that avoid populating that level. There is an example for that in the Python implementation of Krotov.

Either type of running cost is significantly easier to realize with GRAPE, see M. H. Goerz, S. C. Carrasco, and V. S. Malinovsky, Quantum optimal control via semi-automatic differentiation, Quantum 6, 871 (2022). That's also why implementing any of these for Krotov is very low priority. The state-dependening running cost described in the paper hasn't been implemented in GRAPE.jl yet. It's on my to-do list, but won't happen until I actually need the feature myself (again, I'm open to PRs, and that's a much smaller project than adding the same capability to Krotov).

I was thinking about defining a cost functional J_T(Ψ, trajectories), where I compute the time evolution with the information given in trajectories and then use it to compute my constraint terms.

No, that completely breaks the derivation of Krotov's method. Remember that the method needs to know ∂J/∂⟨Ψ|, and calculating that derivative will not be possible with the approach you're describing.

I assume that g_a quantifies the distance of the current pulse ϵ⁽ⁱ⁺¹⁾ to the pulse of the previous iteration ϵ⁽ⁱ⁾.

Yes, that is correct, and is (up to the variations discussed above), a core component of Krotov's method

Is it possible, to directly define a custom g_a(ϵ⁽ⁱ⁺¹⁾) that takes the current pulse ϵ⁽ⁱ⁺¹⁾ as input?

No, since that would give an entirely different update equation (or, in most cases, not provide an equation at all). In any case, the update has to be hand-coded for a particular choice of g_a.

So the bottom line is that this implementation of Krotov's method does not support any kind of running costs. It is unlikely to gain such support in the near future (unless someone else really motivated wanted to implement it), since GRAPE is more effective at dealing with running costs, and there are simpler alternative approaches within Krotov's method to achieve the desired outcome.

[HuelsMat]

Thanks for the detailed answer and references.
That has clarified a lot!