Optimal Control of Vibrational Cooling

This example illustrates the use of the qdyn library to optimise a pulse for the purpose of cooling an arbitrary initial ensemble of molecules (approximately constrained to a subset of energy levels on the ground electronic surface) to the ground state of the system via optical pumping and consequent decay. Furthermore it is illustrated how to analyse the resulting optimised pulse such as to get useful and detailed information about for example pulse fluence and intensity as well as raw data about the cooling process implemented by that pulse that can be further visualised with tools like MATLAB or Mathematica.

The example is performing an optimisation for the Cs2 molecule using an assembly-line mechanism as the optimisation goal. 11 states, including the total ground state, are considered on the ground electronic surface.

Running the Example

You must compile the example by running make. Then, run the resulting cooling program, which performs the optimisation, by executing

./cooling

Afterwards you should run the analysis program which can be done by executing

./cooling_analysis

Running the MATLAB code FC_plot.m prints out a file called fcf_fancy.eps which shows the Franck-Condon map of the Cs2 molecule.

Executing the Mathematica notebook will generate the files Bar3D.pdf and Bar3D.eps which will show a 3D bar plot of the population in the different channels (i.e. the considered vibrational levels on the electronic ground surface). Furthermore this information is also displayed in two animated 2D bar plots, once for 50 cycles in Simulation_50.gif, and once for 15 cycles in Simulation_15.gif.

If one desires to delete all data except the MATLAB script, the Mathematica notebook, the program itself and the config files one can simply exectute

make clean

Background

All information about the physical background, the philosophy of the cooling mechanisms that are being optimised as well as the specific functional used can be found in the corresponding paper [1].

Details - FORTRAN code

The FORTRAN files for the two programs, as well as the auxiliary globals.f90 module are extensively commented so please refer to those comments for further details. However the role of the config files as well as some important general facts will be summarised in the following.

The potentials that are being used in the optimisation can be found in the Potentials folder. The file gs.dat corresponds to the ground state electronic surface and the file es.dat corresponds to the excited state electronic surface. Both potentials in this example are for the Cs2 atom.

In the run folder r001 which is automatically accessed by the programs one can find two config files. The first file, config, is responsible for the optimisation performed by the cooling program. The second file, config_analysis, is responsible for the analysis performed by the cooling_analysis program.

The optimisation itself considers the first eleven eigenstates on the ground surface including the ground state. Inside the program this is given by the variable m which can be changed freely (for the analysis it should be changed in this program, too). The amount of states used for a cut off when it comes to spanning a unity is given by n. This number cannot be chosen freely since at some points the amount used in the example n=50 is hard-coded. However, if one desires, it is possible to substitute this amount in all occurrences.

In the following the parameters given in the config file will be explained as well as their possible options and corresponding impact on the optimisation.

The optimisation time is around 24 ps with a little over 4000 time steps. A Gaussian is chosen as a guess pulse, centered in the optimisation time interval with a temporal width of around 100 fs. The peak amplitude is chosen with 10^-5 internal units as low but not too low such that the guess already excites relevant transition for the cooling process.

The number of optimisations is set to 100 which will already take a few hours of optimising. The code runs parallelised on up to m cores using the automatic parallelisation features of qdyn if they were enabled when compiling. Furthermore the optimisation is set to complex since an implicit rotating wave approximation is used (note that the ground surface is shifted up in the config file by the asymptote parameter E0), allowing for complex pulse results that become real when combined again with the double-frequency counterpart. Complex optimisation is recommended since it improves convergence behaviour.

The optimisation is performed with krotovpk since the assembly-line functional is convex. While the lock-in functional is non-convex, it is so only slightly such that in general a first-order Krotov optimisation leads to monotonic convergence.

The dipole moment between the two surfaces is constant and set to one. For the proper physical dipole experimental data could be plugged in here freely. A scaling of the dipole would merely lead to a scaling in the resulting electric field.

In the user_reals section the weighting factors of the functional are given. For this test they are all set to one which proved itself to be a fine choice. For the assembly-line functional the assembly_weight is relevant while the symmetry-weight is ignored and vice versa for the lock-in functional. Furthermore in this section the physical asymptotical energy difference between ground state and excited state should be given as the variable asymptote. This is important to calculate the proper Einstein coefficients determining the transition properties due to radiative decay from the excited surface to the ground surface.

In the user_logicals section an option is given to turn the verbose mode on (by default, if not given in the config file, it is set to .false.). The verbose mode prints out information every time the functional is called including the current progress towards the individual goals of the optimisation functional.

Finally, in the user_strings section one can choose between a functional that optimises assembly-line type cooling and a fucntional that optimises symmetrised/lock-in type cooling. The first is triggered by setting functional=assembly, the second is triggered by setting functional=lock.

For the cooling_analysis program the config file should be very similar to the one actually used in the corresponding optimisation from where the to be analysed pulse had been obtained from. Instead of a Gaussian guess pulse for an optimisation the pulse given in the config file should be the pulse given as the output from the optimisation. The weighting factors for the functional are irrelevant since the functional itself is never evaluated here but it is conceivable that additional functionality might require them so they are kept in the config file.

Just like the main program, the analysis program allows for a verbose mode that is deactivated by default. It gives several interesting information for consistency checks as well as additional output that gives information about the system and the action of the optimised pulse on it.

Details - MATLAB script

The MATLAB script uses the file fcf_fancy.dat which is an output file from the cooling_analysis program. It plots the Franck-Condon parabola for the first 15 states on the excited and ground surface. Inside the script the axis description as well as the plot range can readily be changed. However note that by default the output of the Franck-Condon factors is only done for the first 50 levels by default.

Details - Mathematica notebook

The Mathematica notebook uses the gs_to_gs_pathways.dat and the es_to_gs_pathways.dat files to simulate the action of the pulse run through the analysis program on the system for multiple cycles. In the beginning of the notebook the starting distribution of population can be set (by default it is an equidistribution on all levels except the ground state). Afterwards plots and corresponding output files are generated. Note that the notebook 150 cycles by default and is written for 11 states including the total ground state.

Just like the FORTRAN code the Mathematica notebook is annotated so please refer for further details to these comments.

Bibliography

  • [1] D. M. Reich, C. P. Koch Cooling molecular vibrations with shaped laser pulses: Optimal control theory exploiting the timescale separation between coherent excitation and spontaneous emission, http://arxiv.org/abs/1308.0803, 2013.