MCTDH: Output Documentation

Back to Main Menu

Output Documentation


Files are created in a directory specified in the input file. Which files are created depend on the options chosen. The files created by the MCTDH program can be divided into three classes: read-write files, information files, and data files.

The read-write files are used for program section intercommunication:

The information files document what has happened during a calculation:

The data files provide the output needed for analysis of the propagation:

Which data files are created depends on the specification given in the Run-Section of the input file.
In a propagation or relaxation run the following file(s) are always created: check, log, op.log input.
The following files are created by default, but may be omitted: timing, output, speed, stop, update. (update is created for CMF runs only).
The following files may be created on request: auto, autoe, auto1, auto2, ctrace, cross, psi, veigen, expectation, orben, pdensity, gridpop, and steps.

In a diagonalisation run the following file(s) are always created: check and eigval.
The following files may be created on request: timing, lancvec, and eigvec.

In an improved relaxation run the following file is always created: rlx_info.

Various of the read-write and data files have a header section at the start, which contains information defining the data. This header may be composed of what are termed standard definitions, sets of data common to different sections of the program.

The Read-Write Files

The DVR File

The DVR file is written by the gendvr section of the program. It contains the information on the DVR basis, DVR grid, and first and second derivatives. The following is written, unformatted:

HEADER: progver
creation date



DATA: (ort(g),g=1,ortdim)

The Oper File

The oper file is written by the genoper section of the program. It contains the information on the operator, i.e. a description of the Hamiltonian terms, and the operators in the relevant DVR basis required. The following is written, unformatted:

HEADER: progver
operator title

do s=1,nstate






do k=1,ktot

DATA: operator matrices/vectors

The Restart File

Every tout fs the restart file is rewritten. This contains the wavefunction and other relevant information, allowing a calculation to be continued, restarted or recovered. The following is written, unformatted:

HEADER: progver
creation date
integrator information

DATA: time

where progver is the MCTDH program version with which the file was written, and time is the time (in fs) for the wavefunction.

In a diagonalisation run the restart file additionally contains

diagonalisation flag (always true in a diagonalisation run)
second Lanczos vector (the first one is stored in wavefunction, see above)
number of iterations made so far
on-diagonal elements of tridiagonal Lanczos matrix
off-diagonal elements of tridiagonal Lanczos matrix

These data are required to continue the Lanczos iteration with a larger order.

The Stop File

If the keyword stop is given in the Run-Section, the program opens the stop-file and writes the following text to it:

hh:mm  hh:mm:ss  comment 
real    cpu      optional

 Replace "continue" by "stop   comment" or "hhh:mm  hhh:mm:ss  comment"
 The first input (hhh:mm) refers to real-time, the second to cpu-time.
 One of the time entries may be left blank.
 The colon ":" must appear at the correct position, e.g :
032:30             system update
 will stop the run after 08:30 next day (8+24=32)
 NB: The day-count starts after the first read of "hhh:mm",
 e.g. after "continue" has been replaced by "hhh:mm".
 A first line which reads:
        000:56:30  time limit
 will stop the run after the first output later than 56.5 cpu minutes

 When this file (the stop file) is removed, the run will terminate
 after next output. One may use the keywords "tstop" and "tcpu"
 (Run-Section) to set the times "hhh:mm" and "hhh:mm:ss".


In the first line this file contains the word continue. If during the calculation this word is changed to stop, then the calculation will be stopped after the next output of data. If continue is changed to time given in the format hh:mm (followed by ten blanks), then the calculation will be stoped after the first output of data that happens after the indicated time. A hh number that is larger than 24 refers to the next day. The (optional) comment is written to log- and output-files. The string <tout> is overwritten by the actual t-out time.
Similarly, if in column 8-15 one writes the maximal CPU time in the format hh:mm:ss, then the calculation will be stoped after the first output of data that happens after the indicated CPU time has been used. Both times may be given or one of the two time-inputs may be left blank.
The stop file can be edited while the job is running, or the limits for real-time and CPU-time may be supplied via the keywords tstop and tcpu.
If the stop-keyword is given in the input but the stop-file cannot be opened, the program will stop. Hence, alternatively to replacing continue by stop one may just delete the stop-file.
The stop-file is deleted when the program finishes. It may hence serve as a lock-file. A shell script may search for a stop-file and starts a continuation run only when there is no stop-file.

The Information Files

The Log File

This file contains information such as type of calculation performed, integrator used, information on DVR-grids, use of explicit uncorrelated Hamiltonian terms, title of operator, specified and minimum possible memory array sizes, which data files are opened, any error messages, and much more "garbage". Since this file is not used by other programs, feel free to add any more information to the log file.

The information provided by the log file can be very helpful, in particular when searching for errors. One should always carefully inspect the log file.

The Op.log File

This file contains information about the used operators. The used description correspond to the internal build of the operators.
For the sake of a better understanding of this file we will clearify its structure with an example op.log file to this related operator input .

The used description correspond to the internal build of the operators.
For the sake of a better understanding of this file we will clearify its structure with an example op.log file to this related operator input . There are two types of operator-terms of each used operator:

Additional information on how to interpret the output of genoper to the op.log file can be found in the code description of genoper.

The Output File

The output contains some standard results and will be saved to file or written to the screen, depending on whether or not the keyword output is given in the Run-Section of the input file. In a propagation or relaxation run every tout fs, the following information is output:

Total: Time (fs), CPU-time (s), A-norm,
Etot (eV), Ecorr (eV), Delta-E (meV)
For each state: population, Ecorr (eV), Etot (eV)
For each state:
For every mode: Natural orbital population
For every mode: Position expectation value <x>, standard deviation Sqrt[<x**2>-<x>**2],
Momentum expectation value <p> (for FFT) or
grid population expectation value <n> (for DVR),
standard deviation of p or n

NB: The expectation values and standard deviations are normalized with the population of the individual electronic state.

The A-norm is the norm of the coefficients, Sqrt[Sum(A**2)]. For each state this is the state population, while the total should remain as 1.0. Ecorr is the correlated Hamiltonian energy, Etot the total energy, and the total Etot should be conserved (not necessarily when mean-fields are held constant). Delta-E is defined as Etot-Etot(t=0). In a multi-packet run, i.e. when npacket > 1, the total energy and the norm of the wavefunction, as given in the "total" part, are averaged over the packets.

At the end of the propagation, a remark is output whether the program terminated normally or has been stopped using the Stop file. Finally, total CPU-time, host, date and time, path of the name-directory and (if specified in the input) the title of the run are printed.

In a diagonalisation run, only the number of Lanczos iterations, the number of eigenvalues found, and the length of the Lanczos vectors is put out, e.g.

Number of Lanczos iterations:      80
Number of eigenvalues found:       60
Length of Lanczos vectors:        128

Again, a remark is given whether the program terminated normally. Finally, total CPU-time, host, date and time, path of the name-directory and (if specified in the input) the title of the run are printed.

The Input File

This file contains a copy of the input data, i.e. of the *.inp and the *.op file.

The Data Files

The Psi File

If the keyword psi is given in the Run-Section, the file name/psi is opened and the wavefunction array psi is saved every tpsi fs. The output is unformatted complex.


DATA: wavefunction

The wavefunction is written in the order in which it is stored, i. e.:

For each state: A-coefficients
For each mode: single particle functions for 1st state
single particle functions for last state

It is written using an implicit loop (in order to save space on disk), i. e.:
write(ipsi) (cmplx(psi(d)), d=1,dgldim)

For a graphical illustration of the structure of the wavefunction arry click here.

Note1: psi may be written in single or double precision or in compact form. See the Run-Section of the input file for further details.

Note2: In the density operator programs the psi array contains the density operator.

Note3: If the wavefunction is in exact format then all elements of the A-vector are set to 1 (there are as many A-vector components as there are electronic states). There is only one (combined) single-particle function, the numerical exact wavefunction.

The Gridpop File

If the keyword gridpop is given in the Run-Section , then every tout fs the population of each grid point and each momentum grid point (or basis representation) for each degree of freedom is output, unformatted real*8, to the gridpop file in the following format:


For each state: (only when gridpop=el is given in the input file. Otherweise the population is summed over all states.)
  For each mode:
    For each grid point: population of the grid point

The Check File

The check file is written to every tout fs. It contains the wavefunction norm (calculated from the coefficients only), the natural orbital population values for each degree of freedom etc. to be read by the Analysis programs. The ouput is unformatted real*8. The first number is the program (i. e. output format) version. Then follows, for each TOutput output time, the data in the order:

HEADER: progver

DATA (for every propagation step):  (spop(s),ecorr(s),etot(s),s=1,nstate)

   do s=1,nstate
      do m=1,nmode

   norm, trace                                    ( only for densities,
                                                    since Ver. 83 Rev. 5 )
   do s=1,nstate
      do m=1,nmode
         do n=1,nspfdof(m)
            q(i,s),dq(i,s),expn(i,s),dexpn(i,s)   (since Ver. 82 Rev. 2)
For each state: norm, correlated energy, total energy
For each state:
For each mode: natural orbital populations

The Auto Files

If the keyword auto is given in the Run-Section, every tout fs the auto-correlation function

c(2t) = <psi*(t)|psi(t)> with t=n*tout

is output. This is used to calculate the absorption spectrum. The argument error opens additionally the autoe file which contains an auto-correlation function computed with the least important natural orbital omitted.

If the parameters order1 or order2 are given in the input file, the first and second order auto-correlation functions

c(n)(2t) = <psi*(t)|Hn|psi(t)> (n=1,2)

are written to the files auto1 and auto2, respectively. These files are employed in the filter-diagonalisation method.

In a multi-packet run, i.e. when npacket > 1, the files auto, auto1, and auto2 contain cross- rather than auto-correlation functions. The cross-correlation functions are defined as

cab(n)(2t) = <psia*(t)|Hn|psib(t)>,

where a and b enumerate the wavepackets. These matrices are stored columnwise in the corresponding files, e.g. c11,c21,c12,c22 if there are two packets, with one matrix element per line.

The output of the files is formatted and in the format

#  time[fs]     Re(autocorrel)     Im(autocorrel)     Abs(autocorrel)
   0.00000000   1.00000000000000   0.00000000000000   1.00000000000000
   1.00000000  -0.20997305604067  -0.94892462752682   0.97187788995839
   2.00000000  -0.81853614184398   0.36544890012746   0.89641191095902
   3.00000000   0.44558991083218   0.65579825196272   0.79285668056263
   4.00000000   0.50003471497158  -0.46061533315106   0.67985380877846
   5.00000000  -0.43329600407928  -0.37165953561295   0.57085570642949

Note that for auto=twice the autocorrelation functions are written every 0.5 TOutput fs, yielding c(t) instead of c(2t). See the Run-Section of the input file for further details.

The Cross File

If the keyword cross is given in the Run-Section, every tout fs the cross-correlation function

c(t) = <psiref|psi(t)> with t=n*tout

is output to the cross file. See the Run-Section of the input file for further details.

The output of the file is formatted and in the format

#    time[fs]         Re(cross-corr)     Im(cross-corr)     Abs(cross-corr)
    0.000000000       0.01022653944480   0.00060353257912   0.01024433309640
     1.00000000      -0.00797138392581  -0.01637809685553   0.01821496687622
     2.00000000      -0.01698410940811   0.02026027317266   0.02643744770997
     3.00000000       0.03269164125237   0.01153504157016   0.03466699571349
     4.00000000       0.00089263280612  -0.04341396275159   0.04342313847620

The Ctrace Files

If the keyword ctrace is given in the Run-Section, every tout fs the function

c(t) = Tr( rho(S1,S2) ) with t=n*tout

is output, where S1,S2 are given states. This is used to calculate the absorption spectrum. The output of the file is formatted as the auto files.

The Expectation File

If the keyword expect=S is given in the Run-Section, every tout fs the expectation value of the operator S is written to the expectation file in ASCII format. If multiple operators are given, for every operator there are two columns (real and imaginary part) in the expectation file. Note, output is the norm of the WF (not norm**2) and the expectation value(s), i.e. <psi|O|psi>/<psi|psi>. The information in this file is equivalent to the information produced with the analyse program expect .

The Pdensity File

If the keyword pdensity is given in the Run-Section, every tout fs the one-particle density is written to the pdensity file. If no arguments are given with the pdensity keyword, one-particle densities for all modes are output. If pdensity=I is specified, only the one-particle density for the (contracted) mode of number I is output. One may specify up to four mode numbers. This may be done in order to keep the pdensity file small. The pdensity file is read by the analyse program pexpect. Note that the pdensity file contains the one-particle density with respect to the (contracted) modes. The gridpop file, for comparison, contains the diagonal part of the density for all degrees of freedom.
The data is written to the pdensity file in the following form:

         write(ipden) progver, log2int(lgpel)
         write(ipden) nmode,nstate,((dim(m,s),m=1,nmode),s=1,nstate)
         call wrrhoinfo(ipden)
         write(ipden) s,m, s,m, (avec(k),k=1,dim(m,s))
         write(ipden) (psi(k),k=1,subdim(m)*dim(m,s))

Here avec contains the natural populations and psi the natural orbitals. wrrhoinfo writes the DVR- and the grid-definitions (dvrdef and grddef)to file.
Use the analyse routine prhosub to transform the data on the pdensity file to coordinate representation, i.e. to rho(x,x').

The Steps File

If the keyword steps is given in the Run-Section, the steps file is opened. The steps file contains information on the step sizes and the corresponding errors made by the integrators. Note: the information supplied by the steps file may be very useful, but the steps file may become very long! The output is formatted.

Time in fs: 1.00.

In the following line the integrator is output. One is also informed whether the complete MCTDH equation is integrated or only the coefficients or one (and which) mode. Typically this line has one of the following forms:

Complete equation: Adams-Bashforth-Moulton
Mode 1: Bulirsch-Stoer
A: Lanczos

The rest of the output depends on the integrator used.

For the Adams-Bashforth-Moulton (ABM) predictor-corrector it looks like:

Complete equation: Adams-Bashforth-Moulton
    Step         Order      Stepsize in fs          Error
     68            6        0.733652E-002       0.109774E-009
     69            6        0.124267E-001       0.961385E-010
     70            6        0.210485E-001       0.147705E-008
     71            6        0.337200E-001       0.555701E-008
     73            6        0.307112E-001       0.309803E-007
     74            6        0.318530E-001       0.369706E-007
     75            6        0.322134E-001       0.691726E-007

Step counts the integration steps. Note that a failed step is counted but not put out. Order is the integration order of the corresponding step, Stepsize the size of the step in femtoseconds, and Error the error made for this step.

For the Bulirsch-Stoer (BS) extrapolation method the output has the form:

Mode 1: Bulirsch-Stoer
  Iteration  Small steps    Stepsize in fs           Error
      1            2        0.500000E+000       0.000000E+000
      2            4        0.500000E+000       0.301975E-004
      3            6        0.500000E+000       0.161627E-005
      4            8        0.500000E+000       0.593421E-007

Iteration counts the columns in the extrapolation tableau. The integration order equals 2*iteration+1. Small steps is the number of steps made to integrate over a large step the size of which is shown in the Stepsize column. Error denotes the error. For the first iteration the error is always zero and has no meaning. The Stepsize will not change unless the extrapolation does not converge within the allowed number of iterations (given as order in the input file).

For the Runge-Kutta integrators of fixed order 5 and 8 (RK5/RK8), the output looks like this:

Mode 2: Runge-Kutta order 5
   Step        Accepted?    Stepsize [fs]      Error (scaled)       Time [fs]
      1            0        0.900933E-01       0.131160E+01             0.000
      2            1        0.774299E-01       0.602994E+00             0.000
      3            1        0.525410E-01       0.618417E-01             0.077
      4            1        0.743764E-01       0.301910E+00             0.130
      5            1        0.456527E-01       0.206736E-01             0.204

Step counts the integrations steps for the current integration interval. Rejected steps are also put out; Accepted? indicates whether the step was accepted (1) or not (0). Stepsize shows the size of the integration step in femtoseconds. Error (scaled) gives a measure for the error made in this step; if it's above 1.0, the step was rejected, otherwise accepted. Time shows the starting time for the current integration step in femtoseconds. (It should increase by Stepsize when a step is accepted, or stay constant when the step is rejected.)

The short iterative Lanczos or Lanczos-Arnoldi (SIL) integrator can only be applied to the MCTDH-coefficients (when the mean-fields are kept constant), or in a numerically exact calculation, since the Lanczos scheme requires a linear differential equation. The output looks like:

A: Lanczos
   Step          Order      Stepsize in fs           Error
      4            1        0.250000E+000       0.593782E-001
      4            2        0.250000E+000       0.329024E-002
      4            3        0.250000E+000       0.193616E-003
      4            4        0.250000E+000       0.952046E-005
      4            5        0.250000E+000       0.418954E-006
      4            6        0.250000E+000       0.190529E-007

Step counts the Lanczos steps, Order is the order of integration (i. e. all Krylov-vectors up to H^order|psi> are computed). The rest should be clear.

The Davidson integrator (used for improved relaxation only) produces an output similar to the Lanczos one. However, the Stepsize column shows the squared overlap beta rather than the step-size (except for the very first line).

The Update File

If the keyword update is given in the Run-Section , and a constant mean-fields calculation with adaptive update-interval has been selected, the update file is opened. It contains information on the update steps and the corresponding errors. The output is formatted and looks typically like:

Update-step  Updatetime   Error of A      Error of phi   Time in fs
     1       0.40000000   0.139600E-008   0.167591E-007    0.0000
     2       0.40000000   0.329630E-008   0.440428E-007    0.4000
     3       0.80000000   0.406385E-006   0.444941E-006    0.8000
     4       0.80000000   0.130705E-005   0.100519E-005    1.6000
     5                    0.723068E-004
     5       0.40000000   0.251409E-006   0.754364E-006    2.4000
     6       0.40000000   0.122548E-005   0.338185E-006    2.8000
     7       0.20000000   0.511988E-005   0.235645E-005    3.0000

Each line contains the number of the current update step, the update-time, the error of the coefficients and the single-particle functions, and the time (in fs). For a failed step (like step 5), only the number of the step and the error of either coefficients or single-particle functions (whichever caused the repetition of the step) is output.

The error estimate of the MCTDH-coefficients requires a backward propagation. If the Lanczos integrator did not converge within the maximum order but needed more than one step, the Krylov space does not belong to correct psi anymore. In that case, the back-propagation which uses the "old" Krylov vectors would lead to wrong results, so the update-interval is reduced and the step is repeated. For such a failed update step, only the number of the step is given in the update file (as in step 7).

The Timing File

If the keyword timing is given in the Run-Section, then various timing statistics are written to this file. For each listed subroutine the number of times the routine is called and the total CPU time (user and system) spent in the routine and in all subroutines below it are written, as well as the user time per function call.
Note: not all entries occur during every calculation!

Name Routine-Description
ABM The ABM integration scheme
BS Bulirsch-Stoer integrator
DAV Davidson diagonaliser (improved relaxation)
cross Calculation of cross-correlation
Density Calculatin of density matrices
FBRint Computation of FBR-integrals in a spherical harmonics basis
Funk0 Calculates the mean-field and inverse density matrices and the time derivative of the MCTDH-coefficients without computing the time derivative of the single-particle functions
Funk1 Computes the time derivative of the single-particle functions without recomputing the density and mean-field matrices
Funk2 Computes - after having propagated the single-particle functions - the matrix elements needed for the propagation of the A-coefficients
Funk3 Performs an evaluation of the right hand side of the MCTDH working equation for the MCTDH-coefficients. Used when only the A-coefficients are propagated.
Funka MCTDH module that calls the routines that compute the right hand side of the MCTDH working equations
Funkex Performs an evaluation of the right hand side of the MCTDH working equation for the MCTDH-coefficients. Used when a numerically exact wavepacket is propagated.
Funkphi N Computes the time derivative of the single-particle functions of mode N with mean-field and inverse density matrices held constant. The projector is applied after the uncorrelated Hamiltonian terms have been added
Funkr Calculates the derivative wrt time of the MCTDH wavefunction
genDVR The build of the DVR
genOPER The build of the PES and the operator
initWF The initialization of the wave-function
Lanczos Diagonalizing a real symmetric or complex Hermitian matrix using the Lanczos algorithm
MFields Calculate mean-fields and - if required - dtpsi (adim)
opWF The operation on the inital wavepacket by a defined operator
orben Computation of the averaged mean-field and the energy orbitals
Output The system requirments for the various output (e.g. expectation values)
SIL The short iterative Lanczos (SIL) integration routine
Startup The whole startup and initialization procedure
RK5 Runge-Kutta integrator of order 5
RK8 Runge-Kutta integrator of order 8
trafo Transformation to natural or energy orbitals (CMF)
phihphi Correlated Hamiltonian matrix element computation (ptiming)
calcha Action of the correlated Hamiltonian on the MCTDH-coefficients (ptiming)
summf Sum of the mean-field matrices times the single-particle potential operators (ptiming)
hlochphi1m Product of the mean-field matrices and h|phi> (ptiming)
Funkphi Computes the time derivative of the single-particle functions (ptiming)
dsteqr Computes all eigenvalues eigenvectors of a symmetric tridiagonal matrix (ptiming)

The different statistics give very rich informations for the optimization of a MCTDH calculation. Important are the statistics for the Bulirsch-Stoer (BS) integrator and the numerical effort for the mean-field (MFields) calculations. Click here for an example of a timing file.

The Speed File

If the keyword speed is given in the Run-Section, the speed file is opened. The speed file contains the CPU-time spent between two consecutive outputs and thus measures the speed of the propagation. The data is written every tout fs to the file (except for t= tinit). The output is formatted and consists of four columns:

#    Time[fs]    CPU-time[s]     Delta-CPU      Real time
#        0.00            0.6                    Jan 26 21:40:13 
         0.50            1.0         0.430      Jan 26 21:40:13 
         1.00            1.8         0.800      Jan 26 21:40:14 

The entries are the propagation time, the total CPU-time, the CPU-time since the last output, and the real date and time. In the example, the program needed 0.6 CPU-seconds for setting up the Hamiltonian, the primitive and single-particle function basis, and the initial wavefunction. The propagation of the first 0.5 fs required 0.43 CPU-seconds, that of the second 0.8 CPU-seconds.

The Veigen File

If the keyword veigen is given in the Run-Section and the keywords build and eigenf are chosen in the INIT_WF-SECTION, the eigenvalues and eigenfunctions generated when building the initial wavepacket are written to the file veigen_modelabel. The output is formatted and contains the following information:
The first two lines give some information on the Hamiltonian which was diagonalised and the DOF. Then follows a list of all eigenvalues, and then all eigenvectors are printed (first column: coordinate, second column: value).

The Eigval File

In a diagonalisation run, the computed eigenvalues, intensities, and error estimates of the eigenvalues are written to the eigval file. The file has the form

# Number   Eigenvalue [eV]     Intensity         Eigenerror [eV]
       1   1.0000000000000027  5.645794565E-007  6.764435784E-012
       2   2.9999999999999978  4.608811890E-006  6.764435784E-012
       3   5.0000000000000968  2.083853133E-005  6.764435784E-012
       4   7.0000000000055600  6.858159072E-005  1.919798403E-010
       5   9.0000000002631015  1.829218512E-004  8.624632946E-009
       6  11.0000000085855483  4.184441941E-004  2.654579026E-007
       7  13.0000002009788531  8.498156632E-004  5.821914155E-006
       8  15.0000034612080793  1.567951643E-003  9.312889649E-005
       9  17.0000444999398681  2.671605744E-003  1.099214779E-003
      10  19.0004294936165792  4.256557063E-003  9.576380504E-003
      11  21.0031004573831979  6.413115556E-003  6.078459278E-002
      12  23.0165505793464291  9.264673560E-003  2.738499830E-001

Note that if the keyword energy-not-ev was specified in the input, the eigenvalues and corresponding errors are given in atomic units.

The Lanczvec File

In a diagonalisation run the Lanczos vectors are written (if desired) to the lanczvec file. The file is in binary format. Each Lanczos vector is written as double precision complex numbers using an implicit loop, i.e.:

write(ilancz) (lanczvec(i), i=1,lanczdim)

The length of the Lanczos vectors is given in the output and log files. The number of Lanczos vectors equals the number of iterations specified in the input file, which can also be found in the output and log files.

The Eigvec File

In a diagonalisation run those eigenvectors of the tridiagonal Lanczos matrix that correspond to correct (i.e. not spurious) eigenvalues are written (if desired) to the eigvec file. The file is in binary format. Each eigenvector is written as double precision real numbers using an implicit loop, i.e.:

write(ieigvec) (eigvec(i), i=1,eigdim)

The length of the eigenvectors equals the number of Lanczos iterations, which is given in the output and log files. The number of eigenvectors equals the number of eigenvalues found, which is also given in the output and log files.

The Orben File

The orbital energies and their weights are written at the end of each output interval to the file orben, if the keyword orben is set in the RUN-SECTION. The orbital energies are defined as the eigenvalues of an averaged 1D Hamiltonian, Hav , in the basis of the single-particle functions. The trace of the mean-field operators is taken as the averaged 1D Hamiltonian. The orben file contains the eigenvalues and the weights of the corresponding "energy orbitals". The weights allow to characterise a wavefunction by their single-particle quantum numbers. This is in particular useful in a relaxation run. The average energy is defined as the sum of the eigen-energies times weights. This energy differs slightly from the energy printed in the output, because here one computes the expectation value of Hav , not of H.

The Rlx_info File

The rlx_info file contains information on an improved relaxation run. This file contains slightly different information, depending whether the SIL or the DAV "integrator" is used to diagonalise H in the set of the SPFs. Because DAV is much more efficient, we discuss here only the DAV output. (This output is produced by running $MCTDH_DIR/inputs/hono.dav.inp)

At the beginning and at the end of the rlx_info file the program prints the eigenvalues (in eV) of the very first and the very last diagonalisation. In between these two blocks it prints some information on the orbital relaxation and the diagonalisation steps. This information is much easier to read (because it gets sorted) by running the script rdrlx . rdrlx produces an output like:

 Lost SPFs :   dim - Sum[square of orbital overlap matrixi elements]
    time   s    m=1     m=2     m=3     m=4     m=5    m=5     m=7     m=8
    1.000  1  0.0055  0.3229  0.4793  0.7916                                 
    3.714  1  0.0355  0.3411  1.6491  1.6834                                 
   11.857  1  0.0343  0.3406  2.8734  1.8022                                 
   20.000  1  0.0535  0.3405  3.1553  1.8627                                 
   40.000  1  0.0899  0.3409  2.9795  1.9431                                 
   60.000  1  0.0921  0.3406  3.6452  2.1011                                 
After Davidson   T =   0.000      0.9866963606 eV,  -17.0485200565 meV
After Davidson   T =   1.000      0.9866478851 eV,    0.0017504888 meV
After Davidson   T =   3.714      0.9866355764 eV,    0.0000698917 meV
After Davidson   T =  11.857      0.9866353736 eV,    0.0004596297 meV
After Davidson   T =  20.000      0.9866351837 eV,    0.0011511134 meV
After Davidson   T =  40.000      0.9866347004 eV,    0.0004938465 meV
After Davidson   T =  60.000      0.9866354600 eV,    0.0031944302 meV
After Relaxation T =   1.000      0.9866461346 eV,   -0.0502260015 meV
After Relaxation T =   3.714      0.9866355065 eV,   -0.0123785839 meV
After Relaxation T =  11.857      0.9866349140 eV,   -0.0006624254 meV
After Relaxation T =  20.000      0.9866340326 eV,   -0.0013410604 meV
After Relaxation T =  40.000      0.9866342065 eV,   -0.0009771635 meV
After Relaxation T =  60.000      0.9866322656 eV,   -0.0024347370 meV
   time  order   q   beta*1000     Energy[ev]         ovl*1000    Delta-E  
   -20.0    0    0   0.000E+00    1.003 744 880 671  1000.0000   0.000E+00 
   0.000  208   80   923.06088    0.986 696 360 614   923.0609  -1.705E-02 
   1.000  217   84   998.67071    0.986 647 885 112   922.1870  -4.848E-05 
   3.714  207   83   999.83633    0.986 635 576 422   922.3830  -1.231E-05 
  11.857  220   86   999.91423    0.986 635 373 623   922.2775  -2.028E-07 
  20.000  220   86   999.85515    0.986 635 183 684   922.3543  -1.899E-07 
  40.000  211   84   999.88363    0.986 634 700 376   922.3439  -4.833E-07 
  60.000  227   88   999.72612    0.986 635 462 900   922.2617   7.625E-07 

Let us discuss the last of the four blocks first. The first line, which is given some (artificial) negative time, gives the expectation value of the initial WF. The other data in this line is meaningless. The program then diagonalises the Hamiltonian in the set of initial SPFs. The result is shown in the second line for T=0. (By the way, the energy unit may be set through the keyword rlxunit in the Run-Section). Hence, 187 Davidson iterations were needed for convergence, the seeked eigenvalue is the 74-th one (counting from 0) and the overlap-squared of the new A-vector with the A-vector of the previous diagonalisation (in this case the starting vector) is beta=0.923108 . The column ovl displays the squared overlap with the initial WF. This column is filled with zeros, if the keyword lock is not given as an extension to relaxation. The last column shows the energy change over one iteration step.

The very first diagonalisation lowers the energy by 17 meV as shown by the first line of the second block (After Davidson T = 0.0). Then the orbitals are relaxed by imaginary time propagation over 1 fs. This leads to another lowering of the energy by 0.05 meV (After Relaxation T = 1.0). The diagonalisation within this new set of SPFs increases the energy by 0.00267 meV. And so on. One can show that orbital relaxation will always lower the energy. (If there should be a small energy increase, this is due to insufficient accuracy of the integrator. It can be cured by increasing the integrator accuracy parameter and by lowering the regularisation constant (eps_inv) of the density matrix.) The diagonalisation, however, may lead to a positive or negative energy change.

The first block provides information on how much the orbitals have changed. For each mode dim - tr{P(t) P(0)} is displayed, where dim denotes the number of SPFs, P(t) the projector onto the SPFs of time t, and P(0) the projector onto the initial SPFs. This number may be interpreted as the number of orbitals removed from the initial set of orbitals and replaced by new ones. If this numer is large, the choice of the set of initial SPFs was poor.

Note that the script rdrlax may be called with an option to print only one of the above four blocks. See rdrlax -h.

If the complex version, CDAV, is run, there will be two lines for each time in the last block. The additional lines give the imaginary parts of the energies.

Some numbers may be augmented with a star, *, to indicate that convergence may be lost. A star is added to order when the maximal number of iterations is reached, to beta when it is lower than 2/3 and to ovl when it is lower than 1/2. If beta is very close to one, 1-beta is printed in E-format, rather than beta*1000 in F-format.

An additional block of information may be printed before the last block. This happens, when the Davidson iterates on an other vector than the one with the largest beta. Such an additional block looks like:

  time   order   q  q1    beta   beta1    ovl     ovl1    E           E1
   0.000   82   33  34    1.52  925.96    0.00 1000.00  0.984355  0.987941
   0.000  112   45  47    0.00  925.96    0.00 1000.00  0.984041  0.987941

This block was created because the keyword full was given as second argument to relaxation. In this case the Davidson does not stop after the seeked state is converged, but continues to iterate on the state below. Hence at time=0 (full changes the behaviour of the Davidson only for the very first diagonalisation, i.e. only for T=0) and after order=82 iterations the seeked state q1=34 has an overlap beta1 and an energy E1. The iteration is then continued with the q=33 state with (not necessarily converged) energy E. This state converges after order=112 iterations. The seeked state is now the q1=47 one, but the iteration continues with the q=45-th state, i.e. two below the seeked one. The process continues till the maximal dimension (given in the Integrator-Section) is reached. All eigen-energies of the first diagonalisation are printed to the rlx_info file.
A similar situation may occur, when the keyword lock is given as argument to relaxation. In this case the Davidson iterates (similar as for follow) on the state of maximal overlap with the previous A-vector (i.e. maximal beta). When this state converges it is checked if it has the maximal overlap with the initial WF. If not, q is adjusted and information similar to above is printed to the rlx_info file.

The Spwneval File

If the keyword spawn is given in the  Run-Section, the spwneval file is opened. The spwneval file contains the eigenvalues derived from the optimal generation of SPFs and arranged in decreasing order. Every spawning event the corresponding time is recorded:

  Time in fs:      0.00000
and a summary of the methodological details is given:
  Optimal SPFs (Lanczos):    mode =  1 , state =  1 , maxit =    35
The rest of the output consists of five columns:

  spf        eval*tau^2              eval[au]       eval(i)/eval(1) thr
    1    0.000003885350     0.000003678208304     1.000000000000000   T
    2    0.000003212026     0.000003040781551     0.826701834116186   T
    3    0.000001834031     0.000001736252144     0.472037470718575   T
    4    0.000000063548     0.000000060160451     0.016355912029162   T
    5    0.000000004771     0.000000004516625     0.001227941604539   T
    6    0.000000000176     0.000000000166813     0.000045351823654   F

The entries are the SPF number, the corresponding eigenvalue times tau square (adimensional), the corresponding eigenvalue in atomic units, the corresponding eigenvalue divided by the largest one and a logical string indicating if the corresponding SPF was added (T) or not (F) during the enlargement of the wavefunction on-the-fly.

The Nspf File

If the keyword spawn is given in the  Run-Section, the nspf file is opened. The nspf file contains the number of SPFs as a function of time.

The Standard Definition Data Sets

The following sets of data comprise standard definitons for the calculation. There are three separate sets: for the DVR, the grid (i.e. after any combinations are made), and for the wavefunction.

The dvrdef information set



The grddef information set

dentype                                (since Ver. 81 Rev. 11)
do m=1,nmode


The psidef information set

dentype                                (since Ver. 81 Rev. 11)
if (dentype.eq.1) then
   do s=1,nstate
      do t=1,nstate
   do s=1,nstate

psitype                                (since Ver. 81 Rev. 4)

citype                                 (since Ver. 82 Rev. 0
if (citype.qt.0) then                   -------- " -------- 
   cicut, (ciweight(m),m=1,nmode)       -------- " -------- 
endif                                   -------- " -------- )