MCTDH: Direct dynamics

Back to Main Menu

To Howto Menu

Running Direct Dynamics Calculations

In the following a brief description of how to set up and run direct dynamics calculations is given. For more details and options available look in the main manual. An input file will be required, and a directory set up to contain the quantum chemistry files and DB. What is described below is not the only way to set up these calculations, and more control over the integration and other advanced parts of the program can be obtained if desired.

Like all Quantics calculations an input file is required to control what is done. In addition a directory needs to be set up to hold the quantum chemistry database (QC-DB), and relevant files need to be deposited in this according to the calculation to be run.

Three sections are needed in the input file: a RUN-SECTION, an INITIAL-GEOMETRY-SECTION and a DIRDYN-SECTION.

A simple example is shown for the vibration of the O-H bond in ethanol .

The Run-Section

A direct dynamics calculation is specified by the keyword direct in addition to the keyword propagation. Options can be specified such as direct = cartesian or direct = nmodes to specify the coordinates in which the dynamics will be run. The keyword ngwp = N will specify a basis set of N Gaussian wavepackets. Other keywords, such as tfinal, tout and auto are used to specify the calculation in the same way as an MCTDH calculation.

By default, a single-set vMCG calculation is run, with a single set of GWPs describing the nuclear wavepacket. A multi-set calculation can also be made, in which a separate set of functions is used for each state, by setting ngwp = N, multi. The same number of functions is used for each state. Single-set is usually more efficient for direct dynamics as it requires a smaller number of functions, and hence fewer quantum chemistry calculations.

The DirDyn-Section

The DIRDYN-SECTION contains keywords relevant for the direct dynamics. The most basic is the name of the directory containing the QC-DB, defined by data = etoh_dddata in the EtOH example. Also the qcprogram and method keywords must be given to specify the quantum chemistry program to be used and the electronic structure method. This allows the input files to be correctly set up and output files read by the interface between MCTDH and the quantum chemistry program. The final important keyword is subcmd. This gives the command that would be used on the command line to submit the quantum chemistry job. In the example the computer being used submits GAUSSIAN via a script run_g03. The 2 at the end means the command takes 2 arguments - the input file name and directory where the results will go.

The zero point for the energy must be specified by ener0. Other keywords to specify how the DB will be used, e.g. the frequency and method of updating, can also be specified here. Full a full list of options see the main manual.

The Initial-Geometry-Section

The easiest way to set up the initial geometry, as done in the example here, is using an INITIAL-GEOMETRY-SECTION in which the molecule is specified between the keywords cartesian end-cartesian in XYZ format. If the dynamics are run in Cartesians, keywords can then be used to specify, e.g. the width of the GWPs.

If the calculation is to be run in normal modes, the normal mode information needs to be specified between nmode end-nmode. This gives the label, frequency of the modes and any initial displacements. The Cartesian coordinates information is then ignored except the atom types. To set up the transformation matrices between cartesian and mass-frequency scaled normal mode coordinates a file can be read containing the results of a quantum chemistry vibrational frequency calculation. For example if a GAUSSIAN calculation with the keyword freq=hpmodes has been run, this can be put in the QC-DB data directory and the keyword transfile = etohdistfr.log in the DIRDYN-SECTION specifies it should be used. The frequencies used for the coordinates are taken from the Quantics input file rather than the QC output file.

Input for Non-adiabatic Calculations

If the system to be studied has a manifold of electronic states, i.e. a non-adiabatic propagation, more input is required. An example is given for the Butatriene Cation, using normal mode coordinates. In the INITIAL-GEOMETRY-SECTION the keywords nstates = 2 and init_state = 2 specify that the system has 2 electronic states, and the dynamics starts on the upper state. In the DIRDYN-SECTIONi, as CASSCF is being used for the electronic structure (method = cas ), the number of atomic orbital basis functions also needs to be provided by the nbasis = 44. A file with the initial guess orbitals for the CAS may be also be required by the QC programs, e.g. a checkpoint file for GAUSSIAN, or an orbital file for MOLPRO.

The potential surfaces are diabatised and the scheme to be used must be specified using the dd_diab keyword. To use the propagation diabatisation scheme, use dd_diab = global . as in the Butatriene Cation example already described.

For the simpler regularized diabatic states, the keyword dd_diab = regularized is specified in the DIRDYN-SECTION and information about two molecular structures: coinfile pointing to an output file with a conical intersection; and dreffile pointing to a file with a reference geometry where diabatic and adiabatic states coincide (often the Franck-Condon geometry). These 2 files are the results from quantum chemistry calculations, and these files must be in the QC-DB directory. An example based on the butatriene cation using CASSCF is given in Here with the data files . The template.dat file includes metastrings to generate the 2 surfaces.

The Template.dat File and the QC-DB

For all calculations a file called template.dat is required in the QC-DB directory. This is a template for the quantum chemistry input files to be generated and used througout the calculations. It should be like the input file that would be required with the information that changes during a propagation, such as the instantaneous geometry, or certain keywords, replaced by metastrings. The template.dat file for the EtOH calculations is for a GAUSSIAN calculation at the RHF/3-21G* level. See below for an explanation of the metastrings used in this file.

The template file for for the butatriene cation is for a CASSCF calculation, again using GAUSSIAN. More keywords are now needed to calculate the nonadiabatic coupling. After a calculation has been run a set of files will be written in the QC-DB directory with the type .db. These are binary file and can be the following (exactly which files are created depends on the calculation):

QC-DB Files
File Description
geo.db The geometries.
pes.db The diabatic potential energies.
gra.db The diabatic potential gradients.
hes.db The diabatic Hessians.
apes.db The adiabatic potential energies as calculated by the QC.
agra.db The adiabatic potential gradients as calculated by the QC.
ahes.db The adiabatic Hessians. Either from the QC or from the Hessian updating.
nact.db The non-adiabatic coupling vectors as calculated by the QC.
nmode.db The structure in normal mode coordinates.
trans.db The adiabatic-diabatic transformation matrices.
If the files have the type .dba they are ascii files that can be read.

Analysis Programs

There are a few programs that are specifically for the anaysis of Direct Dynamics calculations.

gwptraj, gwptraj Plots the trajectory at the centre of the GWPs along a single coordinate as a function of time. Can also plot momenta and phase space trajectories.

ddtraj

Reads the geometries at the centre of a selected GWP as a function of time and writes them to an output file that can be plotted (automatically) as a set of evolving structures by the JMol program.

rddddb

From within the QC-DB directory, reads all the geometries in the QC-DB and plots the structures.

trafodb

Can perform various transformations to the DB data to plot, for example, the non-adiabatic coupling vctor field.

convdb

Converts the QC-DB from binary to ascii (and back) so that the data can be examined. This program can also be used to truncate a QC-DB, or delete singe points.

mergedb

Two QC-DBs from different calculations on the same system are merged into one.

cleandb

A QC-DB can be "cleaned" by removing points that are too close together.

Using showsys need to add -DB. The option -ddraw means that the raw QC data is plotted rather than the model surfaces..

A Typical Protocol

In the following a typical protocol for a non-adiabatic DD-vMCG simulation will be described. It will be assumed that the calculation will be run in mass-frequency scaled normal mode coordinates and propagation diabatisation will be used. First prepare the quantum chemistry by deciding on a level of theory and basis set to be used. This may require some benchmarking and exploratory calculations as you need a balance of accuracy, stability at the range of geometries to be covered, and as cheap as possible. Set up the quantics input file with information as described above for the direct dynamics in the dirdyn-section. The following will need consideration. Run a test calculation, i.e. with the keyword test in the run-section along with propagation and ngwp = 1. This stops after the first energy has been calculated. This may take a while as the Hessians for all states are calculated in addition to energies, gradients and derivative couplings. Check that this is successful. Using this QC-DB as a starting point now build up the calculation. Check the calculation (energy and norm conservation, potential surfaces look correct). If there are problems, look at the structures calculated along the propagation using ddtraj and those in the QC-DB using rddddb. Is the level of theory insufficient? If nothing is obviously wrong, then move on to production runs. The aim is to get a calculation converged in both DB points and GWPs. Convergence wrt the DB will be seen if a new run does not any points to the QC-DB (this is written in the log file). It is very hard to obtain, but near convergence when only one or 2 new structures are added is posisble. Likewise with respect to the GWP basis. The simplest way is to run the calculation one after another, adding more GWPs each time and building up the database. Rediabatise the QC-DB before each calculation to ensure the best description of the surfaces by adding the keyword rediabatise in the dirdyn-section.

Potential Problems

QC Templates

The template file required to run the quantum chemistry calculations is based on the standard input file for the code to be used. Metastrings are used to allow Quantics to add the relevant information during the propagation, such as the present geometry, or to use different options at differemt times.

Using Gaussian

The above examples use Gaussian09 to run either HF or CASSCF calculations (UHF, ROHF and (TD)DFT calculations are also supported). There are a number of case sensitive metastrings that are used in the template.dat file to modify the calculations. The metastrings will be substituted during the creation of input files as described. Only one metastring can appear on a line. Any other lines will be simply copied. Remember that Gaussian input files are sensitive to blank lines.

Gaussian Metastrings
Type Metastring Description
Metastrings enclosed by < > are substituted by blocks of data and must be put in correct places. <geometry> present geometry will be written here
<mo-coefficients> MO coefficients used as a guess in a CAS calculations written here when guess supplied by $Guess:Guess=Cards$
Metastrings replaced by the relevant values. $chk$ Will be replaced by a checkpoint file name to match the generated input file name. Use with %chk=$chk$
$root$ The number of roots to be calculated. For example in a CAS calculation use the keyword CAS(M,N,nroot=$root$)
$cpus$ Will be replaced by the number of threads available in a parallel job. Use with %nproc=$cpus$.
$bmem: I $ This keyword is slightly different from others of this type in that it will be replaced by the value I multiplied by the number of threads in the calculation ($cpus$). Use with %mem=$bmem:800$mb.
$time$ Replaced by the current simulation time. Used in the context of Ehrenfest simulations.
Metastrings that provide keywords depending on the calculation being run. $Swap:keyword$ Need to swap the order of states to calculate the gradient of the lower state in a CAS calculation. Typically $Swap:IOp(5/97=100,10/97=100)$ is used so that the IOPs are selected when required.
$Freq:keyword$ The keyword to be used to calculate the frequency. For a simple calculation use $Freq:Freq$. For a multi-state CAS calculation $Freq:Freq IOp(5/17=41000200,10/10=700007)$ means that the derivative couplings are also calculated
$Force:keyword$ If the frequency is not being calculated (not required after the first step if using Hessian updating) then the keywords specified by the $Force:$ metastring are used. For example in a multi-state CAS calculation use $Force:Freq IOp(5/17=41000200,10/10=700007,7/10=0)$ to calculate the DCP without diagonalising the Hessian.
$Guess0:keyword$ keyword specifies how the initial guess in a CAS calculation will be made if the DB is empty. Typically "Guess0:Guess=Read" is used and a file start.chk supplied with the initial MOs. This file is copied automatically to the checkpoint file to be used.
$Guess:keyword$ keyword specifies how the initial guess in a CAS calculation will be made if the DB is not empty. Typically "Guess:Guess=Cards" is used and the MO coefficients are taken from the nearest DB point and written in cards format in place of the metastring <mo-coefficients>

It is also possible to use non-standard routes.

Using Molpro

Using Molpro as the electronic structure method, it is possible to use the propagation diabatisation method for CASSCF for more than two electronic states.

A complete example for 4 state water is included in the inputs directory, the input file being h2o-4s-dd.inp and the associated directory for the Database etc. being h2o-4s-dd_dddata. The following explanation uses that example.

In order to perform DD with Molpro using normal modes, it is first necessary to optimise the molecular geometry using CAS with Molpro. It is best to do this using symmetry as otherwise accuracy can be lost. The file produced contains the normal modes which are used in transforming between Cartesians and normal modes in Quantics. The next calculation to perform is another Molpro calculation using the same active space and the optimised geometry, but this time without symmetry, in order to generate a file containing the MO coefficients which will be used as an initial guess for the calculations performed during the DD.

Both the frequency output and orbital files should be copied to the auxiliary directory where refdb.dat, template.dat and the database are kept.

An example template file is given in Molpro template.dat.

In the template file, note that the $Orbread:....$ construct reads in guess orbitals for the calculations. The first calculation requires an MO coefficient file called init-orbital.orb. Subsequent calculations take their orbitals from files created by Quantics (MO coefficients from each Molpro calculation are stored in the DB). The name of the orbital file is filled in at .

If the no_orbread (or no_orbread0) keyword is set in the eindirdyn-section of the input file, the orbitals will not be taken from the .orb file, but generated at each step it is possible to add directives to generate the guess orbitales using the $No_orbread: ...$ metastring. For example, $No_orbread:hf$ would use an initial HF calculation in place of the Molpro default orbital generation.

Also note that separate CASSCF calculations are needed to calculate the gradients/Hessians of each state and individual calculations are needed to calculate the non-adiabatic couplings between each pair of states (2-1, 3-1, 3-2 etc.)

The Quantics input file for Molpro is much the same as for Gaussian. It can be viewed at h2o-4s-dd.inp.

Of note in the input file are the keywords:

use_dd_sym: this turns on the use of symmetry when using DD in normal modes. It relies on the term symbols for each mode to tell whether the PES along that direction should be symmetric. In order for this to work properly a full shell of GWPs should be used so that GWPs are evenly distributed along the modes (the log file tells you whether the shell is full). The minimal number of GWPs for an N mode system is 2N+1 i.e. 1 in the centre, 2 along each mode either side of the centre. N.B> This only works when the centre of the initial wavepacket is at 0 along all symmetric potentials. In this case the wavepacket is initially centred at 0.0 along the 3B2 mode, which has a plane of symmetry at 0.0 in the PES.

grid X Y Z: In the initial-geometry-section, this set the length of the grid (from Y to Z) and number of points (X) along each mode for analysis purposes. The default is X=101, Y=-10, Z=10.

default: This forces the program to use the integrator accuracy defined rather than performing the automatic loosening of a factor of 10.

gwp_eoms=precond: This causes the program to use a preconditioning method to generate the EOMs for the GWPs (to solve iCL=Y).

qcprogram=molpro: Alerts Quantics to look for Molpro inputs/outputs.

dd_diab=global: This command causes the program to use the propagation diabatisation method to generate diabatic PESs from the Molpro data. This must use the hess_upd keyword (see below).

nbasis=X: This gives the total number of basis functions (contracted) used by the electronic structure calculations. Needed to set up array sizes properly in the program.

subcmd=...: Gives the location (on your machine) of the script to run Molpro.

hess_upd: The program uses the Powell update method to create approximate Hessians of the PESs at each DB point after an initial Hessian has been created at the first DB geometry. The propagation diabatisation only works with this. It can also be used with other methods, and can save a huge amount of time.

dbsave: Causes the DB to be saved in RAM during the course of the dynamics. This saves a lot of time. The DB is still written to file at each new point.

rediabatise: This is commented out here as it is a post-processing option. Along with the genoper keyword in the run-section (don't use propagation), this takes a pre-created DB of adiabatic points and performs the propagation diabatisation again, starting a the central point and moving outwards, diabatising in order of distance from the centre (updating the gradients and Hessians as well). This can smooth surfaces out particularly if a GWP has left the central region then come back again, creating DB points as it goes. The diabatisation is approximate, so shorter paths are to be preferred between points, this allows that.

HAMILTONIAN-SECTION_system: See section below as well. This appends a CAP operator to the system Hamiltonian (see also the LABELS-SECTION) allowing the wavepacket to be annihilated as dissociation occurs.

Using Molcas

It is also possible to use Molcas as the quantum chemistry program. As with Molpro, the propagation diabatisation method for CASSCF allows the study of dynamics with more than two electronic states involved. HF can be used for adiabatic dynamics.

The Quantics input file for Molcas is the same as for Molpro or Gaussian (adjusting the qcprogram keyword in the DIRDYN-SECTION). An example can be viewed at tol.inp. In this example, the Molcas electronic structure calculations are submitted via the script molcasrun. The log files generated by Molcas need to have the .out extension.

A template example file is given in Molcas template.dat, where the properties of all electronic states are computed in a single Molcas calculation. Molcas has a very versatile input format, and users are referred to the Molcas documentation for details. Case sensitive XML type metastrings are used by quantics to modify Molcas input. The following metastrings are available:

Molcas Metastrings
Type Metastring Description
<tag/> : standalone metastring (note the / after tag). They are replaced by data. <geometry/> Present geometry will be written here in xyz format. Note that the number of atoms needs to be written in the template itself.
<cpus/> Replaced by the number of threads available in a parallel calculation. Use for example with:
>>EXPORT OMP_NUM_THREADS=<cpus/>
or
>>EXPORT MOLCAS_NPROCS=<cpus/>
<orbfile/> Replaced with the file name of an orbital file with the same stem as the current input file. Usually in combination with the <Orbread>...</Orbread> metastring (see below).
<tag> body </tab> : open-close metastrings. Replaced by the body depending on the type of calculation required. <Force>...</Force> Used when gradients of the potential energy surfaces are needed. Usually encloses calls to the ALASKA module and its options.
<Freq>...</Freq> Used when second derivatives of the potential energy surfaces are required. Molcas does not do analytic frequencies for state-averaged CASSCF calculations. Usually encloses calls MCKINLEY module and its options, but due to numerical frequencies distortions it may enclose extra SEWARD and RASSCF calls as well.
<NACT>...</NACT> Used in calculations requiring non-adiabatic derivative coupling to be calculated. Usually enclosing an ALASKA call.
<Orbread0>...</Orbread0> Used to specify the path to a orbitals file with the initial guess orbitals in a CASSCF calculation if the database is empty.
<Orbread>...</Orbread> Used to specify the orbitals file build by quantics from orbital coefficients in the database, when a database is not empty. Usually used in combination with the <orbfile/> standalone metastring (see above).

Using QChem

It is also possible to use QChem for TDDFT electronic structure calculations. Propagation diabatisation is the recommended method to diabatize the surfaces.

The Quantics input file using QChem is very similar to the ones for Molpro or Gaussian. An example can be found at water.inp. The script that is executed to launch electronic structure calculations is defined with the keyword subcmd, like for all electronic structure programs. In this example, the QChem electronic structure calculations are launched via the script run_qchem, which has to be adapted for your system and QChem installation. Quantics calls the script with two additional arguments, the name of the inputfile and the directory where the inputfile is located. If you don't have a running QChem version,you can download a trial version from the qchem website.

You may want to create a folder, where the database is stored. In this example we called it "water-dd_dat". You have to put the refdb.dat file and the template file for the QChem calculations in this directory. A template example file for QChem is given in template.dat. Finally, for a calculation in normal modes you have to provide the reference file with the Hessian at the initial point. An example for that can be found in water.out. At the moment, QChem only provides analytical hessians for very few functionals including B3LYP. Any other functional can be used, but the initial calculation of the hessians will take longer. The calculations for hessians and gradients for different states have to be put in separate calculations. They should be put all together in one template file separated by @@@, as done in the example template file. For Quantics to be able to produce correct input-files a comment section has to be added before each hessian calculation, containing 'freq' or 'FREQ'. Calculations following such a comment won't be launched except at the reference point. QChem input files are in principle case insensitive, however for Quantics to read in the output correctly the 'cis_state_deriv' keyword has to be either written as 'cis_state_deriv' or as 'CIS_STATE_DERIV'. Comment sections have to be marked as '$comment' or '$COMMENT' and molecular geometries should start with '$molecule' and end with '$end'. The atomic coordinates will be added by Quantics. For the first geometry the charge and spin multiplicity has to be given.

Now you are ready to start your DD-vMCG calculation using QChem's TDDFT module for the electronic structure calculations. Just use the example files given here. Copy "water.inp" and "run_qchem" to your current directory, make a directory called "water-dd_dat" and copy the "template.dat", "refdb.dat" and the "water.out" file there. Type "quantics -mnd water.inp" to start the calculation. (Type: "quantics -w water.inp" if the directory water-dd already exists, e.g. from a previous calculation. The "-w" option will overwrite all outputfiles in the output directory.) The only file you might have to adapt is "run_qchem" such that it fits to your set up of QChem. (Using it as it is, it is assumed that you can start QChem by typing "qchem" in your shell and 4 cores are used for the QChem calculation.) After the simulation you can compare your results with the output file of the reference calculation, which is output.

Using External Program

It is possible to use any available external program related to the molecule under investigation as the quantum chemistry program. The energy, gradients and Hessians must be generated by the chosen external program at each point. In the case where the Hessians cannot be computed by the external program, a method to calculate the Hessian matrix when gradients are known was developed and can be used for any external program. The source code can be viewed at hessians, which has to be adapted for your external program.

The Quantics input file is the same as for Molpro or Gaussian. An example input file is given in phenol.inp.

Of note in the input file (DIRDYN-SECTION) are the keywords:

qcprogram=external: Alerts Quantics to look for external program's inputs/outputs.

method=ext: This command causes the program to use the method from the external program to generate diabatic PESs.

subcmd=...: Gives the location (on your machine) of the script to run the external program.

The script that is executed to launch the quantum chemistry calculations is defined with the keyword subcmd, like for all electronic structure programs. Before this script is called, Quantics creates a file that uses the name and the location of the Direct Dynamics directory (ddname_in_out) that contains the geometry of the new point in Cartesian coordinates that needs to be calculated. The script should make use of the get_command_argument intrinsic subroutine to read the geometry and then to save in the same file the energy, gradients and hessians so Quantics can then read the former variables. Note that the script should save the results in the file in the desired way so Quantics can read it. In the following example, the source code of the external program is givel in run_phenol.f90, which has to be adapted for your system.

A folder where the database is written must be created as in any other quantum chemistry program. This directory, phenol_data in this example, should contain the refdb.dat, the template.dat and the phenol_freq.log that sets up the transformation matrices between cartesian and normal mode coordinates. The quantum chemistry vibrational frequency calculation can be computed with any desired method. In this example Gaussian was used.

Use of Symmetry

If the keyword db_point_group = S is given in the dirdyn-section then the QC-DB will be generated to fulfill the symmetry for the point group given by the label S. Each time a new point is added to the DB, its symmetry replicas will also be created and added. This not only produces better potentials, but also saves QC calculations. Not many point groups have been added, and new ones must be added to the lib/utilities/symoplib.f90 library.

Partitioning GWPs

For large calculations, better convergence can be achieved by "partitioning" the GWPs to make them lower dimensionality. For example, in calculation of ethene, rather than using a single set of 12D GWPs, two sets of 6D GWPs could also be used. This then uses a G-MCTDH formalism rather than vMCG, with a multi-configurational wavefunction set up using these basis functions. The points generated and used in the Direct Dynamics are then given by the centres of the configurations. Thus if 10 functions are used for each of the 6D GWPs, this gives 100 configurations, and 100 "trajectories" in phase space. These trajectories are no longer independent. See Coonjobeharry et al Phil. Trans. Roy. Soc. (21) for an example on cyclohexadiene. If the keyword partition = I is added to the RUN-SECTION, then I sets of GWPs will be used, with the modes evenly distributed between then. The number of GWPs is still controlled by the ngwp, but this is now the total number of configurations. However, this allocates the degrees of freedom to the basis functions according to order of listing and this may not be a very efficient way. The question is, which degrees of freedom should be combined together in the different sets? As with MCTDH, the best results (fewest GWPs needed) will be obtained if the dynamically important modes, and modes strongly coupled, are together in one set of coordinates. If normal modes are being used, this can be worked out using an analysis of the potential surfaces to see which modes are initially excited and where the couplings are. The Vctrans program can be used for this by setting up a linear vibronic coupling (LVC) model Hamiltonian. and ordering the coupling parameters.

The log file contains the information on the key modes. This can be interpreted to see which are the dynamically most important modes using the common LVC nomenclature. From the ordering of the kappa, lambda and gamma strengths, take the top few from each set into one GWP set, with a preference for the tuning modes, then coupling modes, then gamma couplings. Give priority to the initial excited state and states closer in energy to that one Thus a suitable partitioning can be achieved by adding the following section to define the GWPs in the Quantics input file. 
SBASIS-SECTION
5A,6A,7A,8A,9A,10A = 10
1A,2A,3A,4A,11A,12A = 10
end-sbasis-section

Setting up Operators for Analysis

As no operator file is used, if an operator is required for analysis, e.g. to calculate an expectation value, this can be added to the input file. Simply add the operator in a section HAMILTONIAN-SECTION_oper, where oper is the name of the operator. If you want to append an operator e.g. a CAP to the system Hamiltonian then oper is set to be system. Both a system operator and an extra operator use all the MCTDH operator inout and a PARAMETER-SECTION and LABELS-SECTION can also be used.

Adding A Light Field

In direct dynamics the system Hamiltonian is generated automatically. To add extra terms, such as a light field, add the terms in the input file in a section HAMILTONIAN-SECTION_system. This uses all the MCTDH operator inout and a PARAMETER-SECTION and LABELS-SECTION can also be used. The operators qchemdipx, qchemdipy, qchemdipz and qchemdip are available to use the components or absolute value of the dipole moment as calculated from the electronic structure calculations. The keyword dipole must be added in the DIRDYN-SECTION to generate the database of dipole moments. An example input for a local control calculation is here