Potfit Input Documentation

General Structure

The input file is divided into sections containing information for the various parts of the program. The sections of the input file may appear in arbitrary order. All input is based on keywords with arguments. The input file ends with the keyword end-input.

Upper and Lower Case Characters

Certain parts of the input are processed in a case insensitive manner whereas others are interpreted in a case sensitive way. In general, all upper keywords are internally converted to lower case characters. Therefore any mixture of upper and lower case characters, e.g. EnD-iNPuT is permitted on input for keywords. However, certain strings like

are interpreted in a case sensitive manner!

Special Characters

All lines containing only blanks and/or hash marks '#' are ignored. All keywords that start with an exclamation mark '!' are ignored. Any text following a hash mark '#' is ignored till the end of the line. Any line which begins with five minus-signs or underscore-signs is ignored.

Keyword Placement

In general, keywords are specified in a free format. They are separated by blanks or semicolons. Hence, the lines

keyword1
keyword2
keyword1 keyword2
keyword1 ; keyword2
keyword1;keyword2

are all valid.

However, keywords specifying a block-structure like run-section must be the only one in a line. The same holds for the keywords in the Primitive-Basis-Section and the Separable-Weight-Section . Otherwise, the order and placement of the keywords within a section is arbitrary.

Arguments to Keywords

If a keyword requires arguments, these are normally indicated by an '=' symbol, and divided by ',' symbols. Brackets '(' ')' as well as spaces can also be used to help readability. e.g.

keyword=arg1,arg2,...
keyword = arg1 , arg2,...
keyword = (arg1 , arg2,...)

are all possible formats.

The arguments are optional in most cases. The equal-sign '=' and/or the comma ',' indicates that the keyword following it is an argument. Each of the three lines given below is a correct example input line

output psi timing
output psi = double timing
output psi = double, natur timing
output psi = (double, natur) timing

since the keyword psi "knows" the arguments double and natur . However,
output psi = double, timing
is an incorrect input line, because timing is not an argument of psi but a keyword of its own.

Note that if a keyword has several arguments the program is likely to accept only one unique sequence.

In some cases there is no '=' following the keyword and there are no commas separating the arguments. In addition, the sequence of the arguments is fixed in these cases. This input format is used in the Primitive-Basis-Section and in the Separable-Weight-Section .

Units

Units are specified in the same manner as in the MCTDH input Units.

Single Precision

Please note that all single precision keywords are now deprecated.



Sections

General Remarks

The section XXX starts with the keyword xxx-section and ends with end-xxx-section.

The various sections are listed below. Those with a STATUS of C are compulsory, O marks optional sections.

XXX Status Description
RUN
O
Whether iteration is started, whether to read DVR information, what files are to be opened, etc.
OPERATOR
C
Which surface to be used, coordinate systems, energy cut-offs, subtraction of one-dimensional potentials, etc.
PRIMITIVE-BASIS
C
Definition of primitive basis.
NATPOT-BASIS
C
Size of natural potential basis, contracted mode, etc.
SEPARABLE-WEIGHT
O
Definition of separable weights.
CORRELATED-WEIGHT
O
Definition of correlated, i.e. non-separable, weights.

Below, tables of keywords are given for some input section.

A first table describes the keywords. The number and type of arguments is specified. The type is S for a character string, R for a real number, and I for an integer. For instance, 'keyword = S,R' indicates that the keyword takes two arguments, the first is a string, the second a real number. The status column indicates whether the keyword is compulsory, C, or optional, O.

A second table defines the default values for optional keywords and arguments.



Run-Section

Keyword Status Description
name = S
C
The output files will be written to the directory S. Note: If the -D name option is used (see To Run the Program), the name-string given in the input file is ignored.
title = S
O
The string S is taken as the title of the run and printed to the log and output files. Note: everything that follows the equal sign '=' till the end of the line is taken as title!
ctensor = N
O
The computation of the C-tensor is requested. Eigenvectors and eigenvalues for ALL modes, including the contracted one, will be computed. The resulting coefficients will be written in the binary file 'ctensor'. N is an integer indicating the requested number of SPPs for the contracted mode, which will be identified as 'contr' in the NATPOT-SECTION, as usual.
cutnpot (= natpotname)
O
The previously existing natpot file (natpotname) is read. If no natpotname is indicated the file natpot1, within the name directory, is expected. The number of natural potentials is reduced to the values indicated in the NATPOT-SECTION of the potfit input file. The thus downsized natural potential is written to the file natpot.
rdnpot
O
A previously existing natpot file is read. This is useful to compute the error-measures in particular for natpots downsized with cutnpot. Note that the DVR and grid definitions from the natpot file and the input file must coincide.
ext-natpot
O
Read natural potentials and contracted coefficients from the external files ext-natpot and fdcoeff, respectively, and build a natpot-file. Note that due to the (expected) enormous grid size, the vpot file is (by default) not built and hence no fitting error is computed. If, nevertheless, the fitting error is desired, the keyword (genvpot should be included in the RUN-SECTION.
sp-energy = n1 ... nf
O
Requests potfit to compute the single-point PES energy (exact) corresponding to a given set of DVR indices ni ≤ Ni where Ni is the number of DVR points for the i-th degree of freedom. The indices should be separated by blank spaces. If a natpot file is provided then the corresponding single potfitted energy is also computed.
cutalong = n1 ... nf
O
Requests potfit to compute the PES energy values along a cut defined by the provided DVR indices. One of these indices must be set to zero, and a cut along this DOF will be computed. E.g. for a 3D potential one may give   cutalong = 3 0 5   and a cut along the second DOF will be computed while the values of the remaining coordinates are set to the values of the 3rd or 5th DVR-point, respectively. If a natpot file is provided then the corresponding potfitted PES-cut is also computed. For this, create an empty name-directory and copy the natpot-file to this directory.
niteration = I (,I1)
O
(Can be abbreviated to niter.)
I denotes the number of iteration steps to be performed. I<0 enables the interactive mode where the user is prompted for the number of iterations. If this keyword is not provided, no iteration steps are performed. A Correlated-Weight Section is then ignored. Note that iterations are useless, if there is no Correlated-Weight-Section.
The (optional) second argument, I1, allows to reduce the number of evaluations of the fit error during iteration. These evaluations are quite costly for large grids. The fit error is evaluated for iterations which are integer multiples of I1. The error of the last iteration is always evaluated.
iteration-factor = R
O
R denotes a factor, with which the weights are multiplied when iteration steps are performed. Values of R<1 will be ignored, values of R>2 usually lead to divergence. Values like 1.6 ≤ R < 2 are recommended to speed-up the iteration process.
buffer = I
O
Maximum size (in Mb) of the buffer used when reading the vpot file. Default is 256 Mb. The program is faster when the buffer can load the full potential at once. Compare vpotdim and niobuf which are recorded in the log file.
readdvr (= S)
O
The DVR information will be taken from the DVR file in directory S. Note that the primitive-basis-section of the input file is completely ignored if this keyword is given. A DVR file residing in directory S will NEVER be overwritten or deleted, even if the keywords `deldvr' and/or `gendvr' have been specified. If S is not given, name-directory/dvr is assumed.
gendvr
O
The DVR file will be generated unconditionally. (This is the default)
deldvr
O
The DVR file will be deleted at the end of the calculation.
readvpot (= S)
O
The vpot file S will be read. S will never be overwritten or deleted irrespective of all other keywords. "readvpot" differs from "readdvr": First, S is a FILE, not a directory. Second, the potential parameters written at the head of file S do NOT overwrite the parameters found in the .inp file. Instead, both sets are compared to each other. If they are not identical, the program stops with an error message. The same happens if the vpot cannot be read. If S is not given, name-directory/vpot is assumed..
genvpot
O
The vpot file will be generated unconditionally. A vpot file specified via the "readvpot" keyword, however, will never be overwritten.
delvpot
O
The vpot file will be deleted at the end of the calculation unless it has not been specified via the "readvpot" keyword.
continuation
O
(Can be abbreviated to cont. Equivalent to option -c)
Continuation run. The natpot and vpot files will be read and and additional iterations may be performed. If no iterations are requested, only the fit-error will be computed.
write-only
O
(Can be abbreviated to wo.)
"write only" mode. In this mode only the potential on the product grid is calculated and written to the vpot file. No fitting is done, hence the memory consumption is rather low in this mode. This mode can be greatly sped up by using MPI, see the usempi keyword below. The generated vpot file can be used in subsequent potfit runs by using readvpot.
write-cw
O
The program will write the correlated weights (relevant region) the file name/corrweights and does nothing else (except generating dvr and vpot files if non existing). Note, the weights are written in compact form and in mode order.
overwrite
O
Any files already in name directory may be overwritten. (It is safer not to use overwrite but the option -w).
write_every_fit
O
After every iteration the potential fit will be written to disk. File names are `natpot.iii' where iii is the number of the associated iteration. `natpot.000' contains the initial fit. Useful for fits with long run times and unknown total number of iterations.
output = S
O
The output will be written to the file name/output. The string S may take the values short or long. When long is specified, additional information on grids and weights are printed. short is default, i.e. output is equivalent to output=short .
Note: output is default.
screen
O
The output will be directed to the screen, no output-file is opened. Alternatively to screen one may give the keywords no-output or nooutput.
iteration
O
Some error measures will be written to the file name/iteration.
prodwei
O
If a Separable-Weight-Section has been specified, separable weights will be written to the file name/prodwei.
dvronly
O
The program generates the dvr file and then stops.
no-dcoeff
O
The computation of the contracted coefficients (i.e. D tensor) is avoided.
nofiterr
O
The calculation of errors is disabled. Useful for large systems where the calculation of the error measures is quite time consuming. (See also argument of keyword niteration). It is necessary to set the keyword nofiterr in case of very large systems (number grid points > 2.14d+9) where, due to integer overflow, the vpot file is not computed.
wrallevec
O
The density matrix for ALL modes (Note: including modc) is computed and diagonalized. Then, the subroutine wrevecall writes the resulting natural potentials in the evec (ASCII) and evecb (binary) files.
timing
O
Program timing information will be written to the file name/timing   (default).
no-timing
O
The timing file is not opened.
pestiming = I
O
Special program timing information concerning the calculation of the exact potential and the potential fit will be written to the file name/ptiming (see timer labels vextiming and fittiming ). Ignored if the timing keyword has not been set. If I is > 0 the exact potential and the potential fit are calculated I times over the full grid.
usempi
O
The program runs in a distributed memory parallelised modus using MPI if started with the mpirun command. Currently only the generation of the original potential (vpot) is parallelized, therefore write-only must also be set. The MPI master process is reserved for gathering the data computed by the slave processes and writing the vpot file, hence one must use at least two MPI processes (and at least three to get a speed-up). If vpot-format=2 is used, a more efficient parallelisation strategy is employed -- parts of the potential can be computed out of order, which may be of advantage in heterogenous computing environments.
vpot-format = I
O
Use a different binary format for the vpot file. The default is I=1, which is the traditional format understood by all MCTDH utilities. I=2 is currently only supported in write-only mode. It creates two files, vpot2 which contains the header information, and vpot2.raw which contains the potential data (in DOF order) without any Fortran record markers. The I=2 format is currently only useful for some developers.
transfo
O
Intrinsic coordinate transformation will be used. If given a a transformation-section must be present in the input file.
MGPF related keywords (see note below)
fullcoarse
O
A full representation calculation (i.e. m = Ncoarse, exact calculation) is requested. Natural potentials for each particle are written to the binary files cfullj, where j   (in i2.2 fortran format) is a mode number. Note that the D-tensor is not computed.
printrho
O
The program generates the (binary) rho file including the density matrix elements. This file is needed for MGPF.
onlydcoeff
O
Requests the computation of the contracted coefficients. The D-tensor is written to the binary file fdcoeff. Note that no natural potentials are generated.
finemode = j
O
A single potfit calculation for mode j is performed. The resulting natural potentials are written to the binary file fmodj   (j in i2.2 fortran format). This procedure is devised to provide the natural potentials for mode j on the fine grid while the rest remain on the coarse one.

For more details on the format of the various data files see the Potfit Output Documentation.

A normal run allocates a double precision vector v of full grid size which holds either V'=(V-V1d)*SW or the modified reference potential (MRV) (the latter case occurring during the iterations if a correlated weight is applied). Here, V is the exact potential, V1d denotes an optional 1D potential, and SW is an optional separable weight. Two additional work vectors are allocated, the largest of which has a size reduced by a factor of MAX(GDIM(f)/POTDIM(f)) as compared to v. Here, f runs over all degrees of freedom. The most important steps can be vectorised with this set of vectors. The "raw" potential on the product grid is written to a file "vpot". In a following invocation of potfit in the same directory the program looks for the "vpot" file and, if the latter exists, reads the parameters which precede the potential data. If both sets coincide, the potential is read from the vpot file instead of being recalculated. If the two sets differ, a warning message is issued and the "vpot" file is recreated (just like a dvr file). During a given run, V is loaded several times from the "vpot" file into the vector v.

Note on MGPF calculations: Multigrid Potfit (MGPF) is a variational method by means of which a product form of the potential can be obtained for systems in which the primitive grid size exceeds the integer overflow limit. This is achieved by minimizing the L2 difference between the actual potential and the approximated one. In practice, a MGPF is obtained by successive potfit calculations on two different grids defined by the user: a coarse grid and a fine one The latter is the one eventually used by mctdh). MGPF is actually a Python script which calls potfit several times with different parameters. Here, potfit reads the DVR points from external ASCII files. [The MGPF method will soon be available in the MCTDH package.]



Operator-Section

Keyword Status Description
pes = S {pesopts}
C
S is the name of a potential energy surface which must have been encoded in the mctdh package. The name is interpreted in a case sensitive manner. If a surface depends on additional parameters, one can change the encoded default values by adding those parameters via the string '{pesopts}'. Parameter names specified within the braces '{}' are processed in a case sensitive manner. The selected surface determines which parameter names are possible. For a description of the available surfaces see Hamiltonian Documentation -- Available Surfaces. If you want to fit a vpot file that was created by a different program (e.g. adproj or projection) use the "pes = none" option. This means that potfit uses only the vpot file and does not look for a potential energy surface.

Example:
For the lsth surface use 'pes = lsth {tfac = 0.5d0}'. tfac is a mass-dependent parameter which is needed to transform between binding and Jacobian coordinates. Let r1,r2, and r3 denote the binding and rd,rv, and theta the Jacobian coordinates. Then the transformation rules are as follows:

2D:
  r1=rd-tfac*rv
  r2=rv
3D:
  r1=sqrt(rd**2+(tfac*rv)**2
          -2d0*tfac*rd*rv*ctheta)
  r2=rv
  r3=sqrt(rd**2+((1d0-tfac)*rv)**2
          +2d0*(1d0-tfac)*rd*rv*cos(theta))

For a homonuclear diatomic molecule tfac=0.5d0, in general tfac=m1/(m1+m2). rv denotes the distance between the two atoms of the diatom, rd the distance between the third atom and the centre of mass of the diatom, and theta the angle between rd and rv.

Note: tfac is not used if binding coordinates are selected for a calculation.

pes = pes-file{PATH}
O
Potfit is requested to read a pes file located at PATH. If PATH is 'name' then the name-directory is understood.
oned_S1 = S2
O
Subtraction of the one-dimensional potential S2 for the degree of freedom with modelabel S1. S1 and S2 are interpreted in a case sensitive manner and must exactly match the name of a mode as specified in the Primitive-Basis-Section and a one-dimensional potential curve encoded in the mctdh package.
order = S/I
O
order defines the ordering of the variables passed to the potential routine. If S=dof, then the ordering is as in the PRIMITIVE-BASIS-SECTION. Note that this is the default, order=dof need not to be given. If S=mode (or S=particle), then the ordering is as in the NATPOT-BASIS-SECTION. Finally, one may give a blank separated list of integers, which define the ordering with respect to the dof-ordering.
E.g.:     order = 3 1 2
passes the third dof as the first argument to the potential routine. The first dof becomes the second and the second dof the third argument. The arguments of order must be a permutation of the numbers 1...ndof. See the log-file to inspect if the order is correct.
NOTE: order DOES NOT work for readsrf.
vcut < R
O
Energy cut-off for exact potential energy surface. All potential energy values greater than R are set to R.
vcut > R
O
Energy cut-off for exact potential energy surface. All potential energy values less than R are set to R.
cspot = R1,R2,R3
O/C
Whether an additional centrifugal potential of the following form is added to the potential energy surface (so far only installed for the 3D LSTH-surface in Jacobian coordinates):
cspot(rd)=[R1(R1+1)-2*R2^2]/[2*R3*rd^2]

where R1 denotes jtot, R2 denotes K and R3 denotes mass_rd. 'cspot' is compulsory if option '-jtot' has been specified on the command line.



Primitive-Basis-Section. Ordering of arguments

The primitive basis is defined in the same way as it is done in the input file of an mctdh run (see the MCTDH Primitive-Basis-Section ). There is, however, one important difference. The order of the degrees of freedom defines the ordering of the arguments passed to the potential routine. To see, how a particular pes is implemented see the file funcsrf.F (just type   mcb funcsrf.F ).  Inspect also the log-file log, which is created when running potfit. There one finds a message like:

BKMP2 PES : 3D model in Jacobian Coordinates
Dissociative Coordinate : rd
Vibrational Coordinate  : rv
Angular Coordinate      : theta
A line like Dissociative Coordinate : theta would hint to a wrong ordering of the degrees of freedom in the Primitive-Basis-Section.
The arguments passed to the potential routine may be reorderd with the aid of the order keyword of the Operator-Section.
NB: In MCTDH the ordering of the degrees of freedom in an HAMILTONIAN-SECTION (i.e. the ordering of the mode-line) must be consistent with the potential definition. There the ordering within the Primitive-Basis-Section is arbitrary.



Natural-Potential-Basis-Section

The input defines for each degree of freedom the modelabel, the number of natural potentials to be used in the product expansion and whether a one-dimensional potential should be subtracted. Exactly one mode must be chosen as contraction mode. No special sequence of the modes is required in this section.

The format for each degree of freedom is

mode_label = I, contr

where I is the number of natural potentials which is not allowed to be greater than the number of primitive grid points in this degree of freedom. The order of the arguments is arbitrary. Exactly one mode must be defined as contracted mode, i.e. must have as argument contr. For the contracted mode, the number of natural potentials, I, is usually omitted. In this case, I is set to the largest possible value (= number of grid points specified for this mode in the Primitive-Basis-Section). More than one mode definition can be written on a line.
Two or several degrees of freedom can be combined into a single mode in the same way as it is done in MCTDH. (See the MCTDH   SPF-Basis-Section). If the resulting natpot file is used in a mctdh run the degrees of freedom combined in potfit must also be combined in the spf-basis-section (in same order!). Assume that there is the following combination scheme, indicated by braces {}:

  MCTDH: {r1,r2,r3} {r4,r5,r6}
  
Then the following combinations are allowed in potfit:
  potfit: {r1,r2,r3} {r4,r5,r6}
  potfit: r1, r2, r3 {r4,r5,r6}
  potfit: {r1,r2,r3} r4, r5, r6
  
Here potfit operates either on a full MCTDH mode or on DOFs.
Moreover, when a potfit does not operate on all DOFs of the system, it may operate only on a part of a mode. E.g:
  potfit: {r2,r3} {r4,r5}
  
But the following potfit combinations are not consistent with the above MCTDH combination scheme:
  potfit: {r1,r2,r3} {r4,r5} r6
  potfit: {r1,r3,r2} r4, r5, r6
  
The first combination schem is not allowed because two potfit terms operate on one MCTDH mode, and the second example shows a wrong ordering of the DOFs within a mode.

Example:

For a 3-mode system, with labels rd, rv, and theta, to define a basis with 10 natural potentials in the second, and 6 in the third mode, which is contracted over the first mode, the natpot-basis-section reads:

natpot-basis-section
rd = contr
rv = 10
theta = 6
end-natpot-basis-section

or

natpot-basis-section
rd = contr    rv = 10    theta = 6
end-natpot-basis-section

To contract the degrees of freedom rd and rv into a single mode:

natpot-basis-section
rd, rv = contr
theta = 6
end-natpot-basis-section

or

natpot-basis-section
rd, rv = 10
theta = contr
end-natpot-basis-section



Separable-Weight-Section

The definition of the separable weights is written one line per degree of freedom. The input is free format with blanks dividing the various parameters. The format for each line is

mode_label   weight_type   para1, para2, ...

mode_label is an alphanumeric string attached to the degree of freedom. This must match one of the labels defined in the Primitive-Basis-Section . The specification of the modes is not subject to any specific order. weight_type is an integer denoting the functional form of the one-dimensional separable weight function. Depending on weight_type a number of parameters must be provided in order to define the separable weights. Where appropriate, energy or length units may be appended to the numerical parameters.

The parameters to be input depend on the weight_type as follows:

Weight Type Parameters and Functional Dependence
0   w=1
1   w=(q/(para1+q))**4
2 para1=name w=(1/(1+((wpot(q)-wpot_min)/para2)))**para3
3   w=exp(-para2*(q-para1)**2)
4   w=cos(0.5*(q-para1))
5   w=0.5+1/Pi*atan((q-para1)/para2)
6 obsolete  
7   w=0.5-1/Pi*atan((q-para1)* tan((0.5-para3)*Pi)/para2)
8   w=1-(1-para1)*sin(q)**2
9   w=cos(q-para1)**2+para2

For weight_type=2, name represents the name of the 1D potential used for the weight, i.e. for wpot. This potential must be chosen from the same collection as the one specified via string S2 with the keyword oned in the Operator-Section. wpot_min represents the minimum of wpot on the grid points. The minimun is calculated by the program.



Correlated-Weight-Section

In this section one can specify correlated, i.e. non-separable, weights. The correlated weight is either zero or one, and it is thus sufficient to specify the relevant region, i.e. the region where the weight is one. The relevant region is defined through inequalities, limiting either the coordinates directly or through the energy they represent. The weights can also be given in an ascii file through the keyword readcw. The file contains 0's and 1's in a single column. The order is the same as the corresponding potential energies read with pes=readsrf{} in the operator section. (See Hamiltonian Documentation -- Available Surfaces). Note that if readcw is given, all other statements of the Correlated-Weight-Section will be ignored.

Example:

v < -1.0,eV             # Energy bound defining relevant regions.
                        # Gridpoints associated with energies > -1.0eV have zero weights.
v > -5.0,eV             # Meaning similar as above.
rd < 6.54               # Coordinate cut-offs defining relevant regions.
rv < 3.74,au            # Instead of "<" the symbol ">" is also valid with opposite meaning.
theta < 0.9,sin         # Only gridpoints with sin(theta) < 0.9 and
theta > -0.8,cos        # cos(theta) > -0.8 are relevant for the fit.
readcw = [path/to/file] # if readcw is given, all other lines of this section will be ignored.
                        #
  or
readcw = [path/to/file],no-map # If the weights on file are unpacked and in mode order.
                               # One must set "no-map" to inhibit transformation from DOF to mode order.
                               # Note that the weights which are are generated by POTFIT
                               # are packed and in mode order. Unpacked weights
                               # are usually generated by an external program.

All entries are optional.



transformation-section

In this section one can specify the polyspherical coordinates of the system in terms of vectors and their order. Moreover, it is possible to indicate frozen coordinates and change the order of atoms and their cartesian coordinates. Implementing new surfaces using the transformation is similar to implementing new surfaces except for using trasrf.F instead of usersrf.F.

Keyword Status Discription
order = n1 ... nf O This keyword enables changing the order of Polyspherical coordinates and works exactly similar to the order transformation of TNUM.
cart = n O This keyword enables using cartesian coordinates for one vector. See Zundell and H3O2- for examples.
fro = n R O This keyword allows defining the frozen coordinates of the system, the nth dof would be kept constant at value R.
cord = n1 ... nN O This keyword allows to sort the cartesian coordinates of the atoms according to the PES.

Only the vector definition part is not optional.

First, the vectors of the system are defined. The example is provided for a seven dimensional model of the bihydroxide anion:

Example:

transformation-section
1  1  0.5               X   X  # The first vector, the first dummy atom is on the middle of this vector
                               # and the start and final points of the vector are two new dummy atoms.
2  2  0.0592649997001   O   H  # The second vector, the second dummy atom is on the center of mass
                               # of this vector and the start and final points of the vector are Oxygen and Hydrogen atoms, respectively.
3  3  0.0592649997001   O   H
4  1  0.0               X   H  # The fourth vector, the first dummy atom (middle point of the first vector)
                               # is the starting point of the fourth vector and the final point is the Hydrogen atom.
order = 1 7 8 9 3 5 6
cart = 4
cord = 5 1 2 3 4
fro = 2 1.85
fro = 4 1.85
end-transformation-section

  
The trasrf.F file for H3O2- is given for an example of a surface implemented using trasrf.F.



Example Input Files