Reference.Programs.qm_optimal

qm_optimal: Optimal control of bilinear system

This Matlab function solves bilinear control problems as described in the page on the qm_abncd function. The optimal control problem consists of maximizing the following functional

J = J₁ − J₂ − J₃

by means of iterated backward-forward propagations, trying to satisfy the following requirements, see also the corresponding page in the WavePacket central Wiki.

• Maximize output/observable/target at given final time T
  J₁[x] = 2 ℜ ( ⟨c_j|x(T)⟩+ ⟨x_e|D_j|x(T)⟩ ) + ⟨x(T)|D_j|x(T)⟩
  with state vector |x(t)⟩ and where vectors c_j and/or matrices D_j stand for the j-th (linear or quadratic) output/observables and where constant terms ⟨c_j|x_e⟩ and/or ⟨x_e|D_j|x_e⟩ coming from the shift |x⟩ → |x⟩−|x_e⟩ have been omitted.

• Minimize cost/energy associated with input/control field(s) u(t)
  J₂[u] = ∑_k α_k ∫₀^T u_k²(t) / s_k(t)
  with penalty factors α_k > 0 and with prescribed shape function s_k(t).

• Minimize deviation from exact evolution
  J₃[u,x,z] = 2 ℜ ∫₀^T dt ⟨z(t)+x_e | ∂_t−L | x(t)x_e⟩
  with Lagrange multiplier |z(t)⟩.

The total functional J=J₁−J₂−J₃ becomes extremal if the time evolution of state |x(t)⟩ (forward in time) and Lagrange multiplier |z(t)⟩ (backward in time) are given in terms of the matrices A, B, N and C, D, identical to the description given for qm_control see also the corresponding page in the WavePacket central Wiki

∂_t |x(t)+x_e⟩ = L |x(t)+x_e⟩ = ( A+i∑_ku_k(t)N_k ) |x(t)⟩ + i∑_ku_k(t)|b_k⟩, x(t=0)=x₀−|x_e⟩

∂_t |z(t)+x_e⟩ = − L^† |z(t)+x_e⟩ = (−A^†+i∑_ku_k(t)N_k^† ) |z(t)⟩ + i∑_ku_k(t)|b_k⟩, |z(t=T)⟩ = |c_j⟩+D_j|x(T)⟩+D_j|x_e⟩−|x_e⟩

with |b_k⟩ ≡ N_k |x_e⟩ and with A |x_e⟩ = 0. The optimal control field is calculated as

u_k(t)=−s_k(t}ℑ(⟨z(t)+x_e|N_k|x(t)+x_e⟩)/α_k

For a description of numerical techniques to solve the optimal control problem in quantum mechanics, see also the corresponding page in the WavePacket central Wiki.

Source code

The MATLAB function qm_optimal.m can be found here

Note that qm_optimal.m works only with objects of the main classes given in the following table:

class name	description
ket	state vecors, in eigen/energy representation
rho	density matrices, in eigen/energy representation

Input parameters

The following input parameters are required when calling the function qm_propa.m

parameter	description
string1	A text string specifying the numerical method to be employed (ODE solver)
param1	An additional parameter (optional)

• A choice of numerical solvers used for propagations both of objects of class ket or class rho can be found here. Note that when string1 is not specified by the user, the default is to use a Runge-Kutta method of sixth order (for constant stepsize).

File I/O

The required matrices A, B, N and C, D as well as the vectors x_i and x_e etc are read from input data files ket_0.mat or rho_0.mat. Typically they are provided by previously running qm_abncd which generates them for quantum dynamics in either TDSE or LvNE setting, respectively. Alternatively, input data files can come from qm_balance and qm_truncate or from qm_H2model where the dimension of the LvNE problem has been reduced down to lower dimensionality. In that case, there will be specific file name suffixes, see here for more details.

Algorithms

To solve the control problem, we provide a choice simple of self-written integrators based on Runge-Kutta methods (with constant step size)

• MATLAB function for Runge-Kutta method

To solve the optimization problem, we use iterated backward-forward propagations by Rabitz et al., see the WavePacket central Wiki. Of particular interest may be the following :

• MATLAB function for forward propagation
• MATLAB function for backward propagation
• MATLAB function for calculation of optimal fields
• MATLAB function for derivative of optimal fields
• MATLAB function for calculation of optimization functionals
• MATLAB function for increment of optimization functionals
• MATLAB function for calculation of power spectral density
• MATLAB function for plotting optimization functionals
• MATLAB function for plotting power spectral density

Variables

Note that function qm_init should be run previously in order to initialize all variables in use. Of particular interest are the following structures:

• control
• state
• time