MOLSCAT User's Manual -
The parameters input in
&BASIS specify the
collision type, the quantum numbers and energies of the levels
to be used in the basis set, and the dynamical approximations
(if any) to be used in constructing the coupled equations. The
input parameters common to all collision types are as follows:
specifies the collision type, as described in
Section 1.1. It is constructed
ITYP, which specifies the nature of the
colliding particles, and
IADD, which specifies the
dynamical approximations to be made (if any). The values allowed
ITYP are as follows:
Approximate methods are specified as follows:
ITYPE = ITYP- for full close-coupling
ITYPE = ITYP+10- for effective potential method
ITYPE = ITYP+20- for coupled states approximation
ITYPE = ITYP+30- for decoupled l-dominant approximation
ITYPE = ITYP+100- for infinite order sudden approximation
NLEVEL .gt. 0, it generally indicates the number
of levels to be included in the basis set, as specified by the
JLEVEL array. If
NLEVEL = 0, the
quantum numbers for the levels to be included in the basis set
will be calculated internally. For exceptions and more details
see discussion for individual collision types below.
is an integer array specifying the quantum numbers of the levels
to be included in the basis set. The
is structured differently for different values of
ITYP as described below.
JLEVEL is a
one-dimensional array of dimension 4000, although it is
conceptually two-dimensional for
ITYP .gt. 1.
is an array of
NLEVEL molecular level energies,
corresponding to the levels in the
JLEVEL array. If
all the elements of
ELEVEL are 0.0, or
NLEVEL is 0, the energies are calculated using
molecular constants as described below, and
is not used.
ELEVEL is currently dimensioned at
EUNITS, EUNITC -
EUNITS is an integer specifying the energy units
ELEVEL and for other molecular energy constants
&BASIS (for example,
WEXE, etc., discussed below).
EUNITC is a
character*4 variable which may be used for the same purpose.
Both are interpreted in the same way as the
EUNITC parameters of the
&INPUT NAMELIST block (see
Section 2.4) including the fact
that a nonblank value for EUNITC will override any value input
EUNITS. Although these
input variables have the same names as values input via the
&INPUT data set, they are distinct; a value
example, has no effect on
(default 0) if nonzero specifies a scratch data set unit used to
hold the coupling array. This can provide a great savings in
required memory at the cost of a great deal of I/O. This option
is currently implemented only for
ITYP = 6 and
ITYP = 9.
The parameters used to control the automatic generation of basis
sets and energy levels from quantum numbers and molecular
constants will be described separately for the different values
ITYP below. The methods of specifying the
molecular levels to be included are independent of any dynamical
approximations employed, so that the information given for each
ITYP below is applicable to
ITYPE = ITYP,
ITYP+10, ITYP+20 and
ITYP+30. For IOS cases
ITYPE = ITYP+100) the required input is generally
the same, except that rotational energies are not required for
The input data are used to generate a
((JLEV(LEV,I),LEV=1,NLEV),I=1,NQN) table which
contains a description of the quantum numbers and other
necessary data for each 'level'. The number of columns,
NQN, for each
ITYPE and the meaning of
each column are described below for the different collision
types. Note that the last column,
is always a pointer to the index of level
the cross section array (
For atom - linear rigid rotor systems, the
array is usually generated from input parameters
JSTEP, which have the obvious
meanings. For special purposes,
NLEVEL may be set
greater than zero and a list of
NLEVEL j values
supplied in the array
JLEVEL is calculated from
JSTEP, the energy levels are
generated from the input parameters
ALPHAE (vibrational dependence of
DE (centrifugal distortion
constant). The rotational constant actually used is, of course,
BE-0.5*ALPHAE, and this value may be input
ALPHAE = 0.0 (the
default) if preferred.
NLEVEL .gt. 0, the energy levels may be
specified in the input array
ELEVEL is not supplied, they are generated from
BE etc. as above.
NQN = 2 and
JLEV(LEV,1) = j, the
rotational quantum number.
For atom - vibrating diatom scattering, the
array must contain a list of (j,v) pairs specifying the
rotational and vibrational quantum numbers of the levels to be
included in the basis set. There is no option for generating
this list automatically.
If the level energies are not supplied in
they are generated from the input parameters
(rotational constant) and
ALPHAE (alpha(e)), using
the standard formulae.
DE, the centrifugal
distortion constant is assumed to be independent of vibrational
level here. Note that energies calculated from
etc. are set with the zero of energy at v = j = 0.
For IOS calculations the vibrational quantum numbers are
selected from the
JLEVEL array, and any rotational
quantum numbers supplied are ignored. If vibrational energies
are supplied in
ELEVEL, the first one encountered
for each vibrational state is kept. If all
are zero, energies are generated from
WEXE. If the input
NLEVEL is negative,
a user-supplied routine
GET102 is called to set
NQN = 3;
JLEV(LEV,1) = j, the
rotational quantum number, and
JLEV(LEV,2) = v, the
vibrational quantum number.
For linear rigid rotor - linear rigid rotor scattering, the
JLEVEL array must contain a list of (j1,j2) pairs
describing the rotational quantum numbers of the two molecules.
NLEVEL = 0 on input,
J1MIN, J1MAX, J1STEP, J2MIN, J2MAX
J2STEP, which have the obvious meanings for
molecules 1 and 2 respectively.
If the two molecules are identical, the basis functions
corresponding to (j1,j2) and (j2,j1) are indistinguishable. In
this case, the input variable
IDENT (default 0)
should be set to 1. Only distinguishable pairs (i.e.
.gt. J2) should be supplied. It is also necessary to
specify the statistical weights
to be applied to the different symmetry combinations when calculating cross
sections. The statistical weights for anti-symmetric and symmetric
combinations of (j1,j2) and (j2,j1) may either be specified explicitly in the
WT array as WT(1) and WT(2) respectively, or (if both WT(1) and WT(2) are
zero) may be calculated from the single nuclear spin input in SPNUC. Note
that MOLSCAT will not perform a scattering calculation for any symmetry type
with a statistical weight of zero.
Energy values may be specified in
ELEVEL or (if all
ELEVEL are zero) will be computed from
and the values for
I=2 are zero, the program sets
them equal to the values for
NQN = 4;
JLEV(LEV,1) = j1, the
rotational quantum number of the first
molecule, JLEV(LEV,2) = j2, the rotational quantum number of the second
molecule, and JLEV(LEV,3) = j12, the vector sum of j1 and j2. N.b.
JLEV(LEV,4), the index of the cross section, depends only on j1 and j2,
not on j12.
The basis functions here are combinations of asymmetric top
rotor wavefunctions and linear rotor wavefunctions. In general,
input of the asymmetric top functions follows the capabilities
ITYP = 6 which are
described more fully below. If rotor constants (
C(1)) are specified, the asymmetric
top wavefunctions will be calculated, subject to limitations in
JMIN, JMAX, JSTEP (synonyms are
ISYM. Otherwise, wavefunctions will be input from
NLEVEL .gt. 0, only
this many levels will be input, else input to end-of-file. If
functions are generated from rotor constants and if a valid
IASYMU unit number is also specified, the functions
will be written to this unit in a format suitable for future input.
The asymmetric top wavefunctions are combined with linear rotor
functions specified by
J2MIN, J2MAX, J2STEP;
energies of the linear rotor are obtained from the rotation
constant which must be supplied as
also used as for
ITYP = 1. If
0 is specified, it is used to limit basis functions to
energies (asymmetric top plus linear rotor) less than
In order to allow for finer control of the basis set, specifying
.lt. 0 will try to obtain from
IASYMU the basis functions specified by
JLEVEL(3*i-2) = J1,
JLEVEL(3*i-1) = ITAU,
JLEVEL(3*i) = J2), i = 1,
NLEVEL). The asymmetric top functions for all requested
J1, ITAU must
be available on unit
IASYMU. In this case it is also possible to specify the
energies (asymmetric top plus linear rotor) in the
ELEVEL array; otherwise
energies are taken from
IASYMU for the asymmetric top functions and computed
BE(2) for the linear rotor.
NQN = 8; JLEV(LEV,1) = j12, the vector sum of asymmetric top and linear rotor angular momenta, JLEV(LEV,2) = j2, JLEV(LEV,3) = j1, JLEV(LEV,4) = ITAU, JLEV(LEV,5) = PRTY (a parity designation discussed under ITYP = 6), JLEV(LEV,6) and JLEV(LEV,7) are pointers to the asymmetric top wavefunctions.
Rotor basis sets may be taken either from
JLEVEL input or from quantum number restrictions.
In the former case,
NLEVEL should be greater than 0
and rotor levels should be specified in
JLEVEL(3*i-2) = j,
JLEVEL(3*i-1) = k,
Note that parity is indicated by an even or an odd integer for
even and odd combinations of (+k) and (-k); for k=0, only even
parity is allowed.
Alternatively, the basis set may be specified as including levels for j from
JMIN to JMAX:
JSTEP is not used. The additional input variable KSET (default
0) may be used to limit the k values included: if KSET is negative, only
levels with k=|J2MAX| are included; if KSET is zero or positive, only levels
with k less than or equal to KSET are included. If ALL k levels are required,
KSET should be set to at least JMAX (999 recommended).
The corresponding energy levels may either be input in
ELEVEL, or they may be calculated from the zeroth
order near-symmetric top equation using rotational constants
supplied in the input variables
A, B and
C. Note that these must
correspond, respectively, to moments of inertia about the x, y,
and z axes used in the description of the potential (see below);
this might not correspond with the descending order standardly
adopted in spectroscopy. In particular, the program requires
that the z axis be the symmetry axis of the top (or the axis of
near symmetry) and that the x-z plane be a reflection plane. For
a symmetric top, then,
the unique constant.
NQN = 4; JLEV(LEV,1) = j, the rotational quantum number, JLEV(LEV,2) = k, the projection of j on the symmetry axis, and JLEV(LEV,3) is a parity designation.
Asymmetric top wavefunctions are expanded in terms of symmetric top functions. The necessary functions may either be calculated internally or supplied as data.
If all three rotation constants (A, B and C) are specified, the program
constructs the asymmetric top Hamiltonian and diagonalises it. Note that
(as for ITYP = 5) A, B and C MUST correspond to the x, y and z axes used in
the potential expansion, and are NOT necessarily in descending order of
magnitude. In this case centrifugal distortion constants DJ, DJK, and DK
may also be used; they are input as ROTI(7), ROTI(8), and ROTI(9),
respectively. Rotational functions will be calculated for J = JMIN, JMAX,
JSTEP and kept only if they have energy below EMAX (if EMAX.gt.0 has been
specified). Functions will also be kept only if they meet the symmetry
restrictions imposed by ISYM, which is described below. If a valid unit
number is specified for
IASYMU the resulting wavefunctions will be written
in the format described below for input from
If rotation constants are not input, the program (subroutine
SET6) will try to obtain the wavefunctions by
reading data from unit
IASYMU (specified in
&BASIS, default = 5, the standard input unit).
Each level is described on
ELVL (2I5,F15.10) where
J is the rotational
ITAU is an index that is not
actually used, and
ELVL is the rotor energy (in
&BASIS EUNITS has specified
otherwise). For each level,
coefficients are required (corresponding to k=-J,-J+1,...,+J) and
these must follow the J,ITAU,ELVL card in format (6F12.8) on (NK+5)/6
subsequent cards. Coefficients need not be normalized, but they will be
checked for valid symmetries (by subroutine IPASYM). Note, again, that the
coefficients must correspond to symmetric top functions in the coordinate
system of the potential energy (cf. ITYP = 5 above); in particular, the
rotational constants used to diagonalize the rotor hamiltonian must
correspond to the x, y, and z-axes in which the potential is described.
IASYMU is the standard input unit, data cards
should follow the
&BASIS ... &END set and
&POTL set; in this case a positive
NLEVEL must also be set in
&BASIS, to specify the number of expected rotor
IASYMU is an auxiliary unit, it is possible to
read levels until an end-of-file condition occurs, and
NLEVEL .le. 0 will cause
SET6 to do this. If
NLEVEL = 0, input
functions will be skipped if they do not fall in the range
JMIN, JMAX, JSTEP (if
JMAX .gt. 0 is
specified) or if they are above energy
EMAX .gt. 0 is specified). It is possible to
choose specific levels from
IASYMU by setting
NLEVEL to a negative value. In
that case, input functions will be skipped unless JI and ITAU match values in
NLEVEL). In this case indexing order in the
SIGMA (cross section) array will correspond to the ordering in
Note that for IOS calculations (ITYPE = 106) rotor wavefunctions are used
only to calculate state-to-state cross sections from the "generalized IOS"
cross sections. For this case input of rotor wavefunctions is currently
limited; in particular, generation of rotational wavefunctions from the
rotational constants is not supported and they must be explicitly supplied as
IASYMU if state-to-state cross sections are desired.
Asymmetric top functions are characterised by two symmetries:
Internally, the program assigns one of four symmetry types to
each asymmetric top function (in subroutine
If the monomer functions are calculated from rotational
constants (but not if they are input explicitly from
IASYMU), the symmetry types that are actually
included in the basis set may be restricted using the input
ISYM. The selection criteria that are
possible in MOLSCAT version 11 and later are much
more general than in earlier versions.
In the current version,
ISYM is interpreted
bitwise; a particular asymmetric rotor function is included if
it passes all the appropriate tests.
.gt.3 are excluded
ISYM, start with 0 (including all
functions) and add 2**n for each class of functions to be excluded.
To decide which states are coupled, consider the symmetries (l,m) present in the expansion of the intermolecular potential:
If either or both of these symmetries is present, it is most efficient to do separate calculations for the (two or four) different sets.
For the special cases of
J=0 and coupled states
K=0, there are additional restrictions on the
functions that are coupled: j+j'+l must be even, and epsilon
must be conserved. Do not be misled by these into believing that
the symmetries hold in more general cases.
For asymmetric rotors, all functions are in principle of
degeneracy 1. However, near-degeneracies can occur, and the
program interprets these as degenerate if the levels (from
diagonalisation) are within a tolerance coded in subroutine
SET6C (usually 1.D-8). To be safe, do not set any
high-bit flags for asymmetric rotors.
The flags that test degeneracy are intended for use with spherical tops, and allow A, E and F levels to be included selectively.
The levels to be included may also be restricted using the
EMAX if desired.
NQN = 6; JLEV(LEV,1) = j, the rotational quantum number, JLEV(LEV,2) is an index to distinguish different asymmetric top levels with the same J-value (e.g., the 'tau' label; this value is actually not used by the program), JLEV(LEV,3) is the PRTY value discussed above, JLEV(LEV,4) and JLEV(LEV,5) are pointers to coefficients of the asymmetric top wavefunctions expanded in terms of symmetric top functions.
The basis set for ITYP = 7 is specified in the same way as that for ITYP = 2.
NQN = 3; JLEV(LEV,1) = j, the rotational quantum number, and JLEV(LEV,2) = v, the vibrational quantum number.
For atom - corrugated surface scattering, the
JLEVEL array must contain a list of (g1,g2) pairs
specifying the surface reciprocal lattice vectors to be included
in the basis set. If
NLEVEL = 0 on input, these
will be calculated from
J1MAX, J2MAX and
EMAX; g1 loops from
+J1MAX, and g2 loops from
+J2MAX. However, not all these functions will be
included in the basis set, as described below.
The dimensions and symmetry of the surface unit cell are specified in the input array ROTI. ROTI(1) and ROTI(2) are the lengths of the real space lattice vectors in Angstroms (irrespective of the units of length used elsewhere). ROTI(3) is the lattice angle in degrees.
To understand MOLSCAT's methods of generating basis sets for surface
scattering, it is important to realise that the operation is carried out in
two steps. At the time that the
&BASIS NAMELIST block is read, MOLSCAT sets
up a master list of basis functions which MAY be included in subsequent
scattering calculations. This takes place before entering the loops over
angles and energies, so must be angle- and energy-independent. For each
basis function G = (g1,g2) covered by the loops over g1 and g2 described
above, MOLSCAT checks that the parallel kinetic energy
(0.5*hbar**2/URED)*|G|**2 is less than the input parameter EMAX, and
includes the basis function in its master list only if it passes this check.
Note that the test involving EMAX is bypassed if
NLEVEL .gt. 0 and the
array is read in explicitly.
Subsequently, for particular values of
ENERGY, MOLSCAT calculates the
parallel component of the incident momentum
selects from the master list those basis functions for which
(0.5*hbar**2/URED)*|K+G|**2 is less than the input
EMAXK. This occurs even if the
JLEVEL array is specified explicitly in
Thus, only those basis functions which satisfy both the
EMAXK criteria will
ultimately be included in the basis set. Since the
EMAXK criterion is usually the most sensible
EMAX serves principally to keep an
automatically generated master list within manageable bounds;
the current dimensioning of MOLSCAT does not allow
the master list to contain more than 200 functions.
For the common case of surface scattering with the incident beam approaching along a symmetry direction, MOLSCAT automatically constructs the appropriate symmetrised linear combinations of basis functions. This takes place as part of the second step described above, since it is not until that stage that the program knows the incident angle.
For normal incidence scattering from square or hexagonal lattices, enormous
savings in computer time are possible if the full symmetry is accounted for.
Special versions of the routine
SURBAS which take advantage of the symmetry
are available from JMH.
To facilitate writing user-supplied code which handles different collision types than those described above, skeleton routines are provided, with numerous comment cards, which provide a template for the required code.
JHALF (default 1) is introduced to facilitate
half-odd-integer JTOT values in BASE9 codes.
The following variables are specifically relevant to decoupling approximations.
(default -1). By default, CS/DLD calculations
are performed for all values of the body-fixed
projection number m up to the largest diatom j
value in the basis set. However, if cross sections
are required only between the lowest few levels,
this is unnecessary. If
JZCSMX is set .gt. -1 on
input, m will be limited by
JZCSMX instead of JMAX.
Cross sections involving levels with j .gt.
are then not valid.
(default 0). This is a historical remnant and
values other than the default are not recommended.
In CS calculations it sets the orbital quantum number
in each channel to L(I) = IABS(JTOT + JZCSFL*J(I)),
with allowed values of JZCSFL = -1, 0, 1).
(default 0). For spectroscopic applications using
the helicity decoupling (coupled states) approximation,
it is useful to be able to include the diagonal part of
the Coriolis coupling exactly. This is achieved if IBOUND
= 1 on input. Note, however, that the asymptotic matching
still uses Bessel functions for integer l, so that the
convergence with respect to RMAX may be slightly irregular
when this option is used. This option is supported only
ITYPE = 21, 22, 23 and 27.
Another effect of
IBOUND = 1 is to suppress coupled
states calculations for projection quantum numbers greater than
an integer array of dimension 3. It specifies the number of
(Gauss) integration points to use for orientation dependence of
the scattering matrix.
IOSNGP(1) is the number of
points for theta. In
ITYPE = 105,/code> and 106,
IOSNGP(2) is the number of points for phi; in
ITYPE = 103, the 3 values are for theta(1),
If values are not given (defaults are 0) the program will try to
choose the minimum number needed for requested values of
LMAX, MMAX (see
Section 2.1) and/or input basis
set rotor levels.
(default 0) controls the type of numerical quadrature on phi for
ITYPE = 103, 105 and
106. The default
requests Gauss-Mehler (recommended) while
0 requests Gauss-Legendre.
Forward to Section 4
Back to Section 2
Back to the Table of Contents