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:
ITYPE

specifies the collision type, as described in
Section 1.1. It is constructed
from ITYP
, which specifies the nature of the
colliding particles, and IADD
, which specifies the
dynamical approximations to be made (if any). The values allowed
for ITYP
are as follows:
Approximate methods are specified as follows:
ITYPE = ITYP
 for full closecouplingITYPE = ITYP+10
 for effective potential methodITYPE = ITYP+20
 for coupled states approximationITYPE = ITYP+30
 for decoupled ldominant
approximationITYPE = ITYP+100
 for infinite order sudden
approximationNLEVEL

If 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.
JLEVEL

is an integer array specifying the quantum numbers of the levels
to be included in the basis set. The JLEVEL
array
is structured differently for different values of
ITYP
as described below. JLEVEL
is a
onedimensional array of dimension 4000, although it is
conceptually twodimensional for ITYP .gt. 1
.

is an array of ELEVEL
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 ELEVEL
is not used. ELEVEL
is currently dimensioned at
1000.
EUNITS, EUNITC

EUNITS is an integer specifying the energy units
for ELEVEL
and for other molecular energy constants
input in &BASIS
(for example, BE, WE,
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 EUNITS
and 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
for EUNITS
. Although these &BASIS
input variables have the same names as values input via the
&INPUT
data set, they are distinct; a value
input in &INPUT
EUNITC
, for
example, has no effect on ELEVEL
values.
IVLU

(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
of 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
IOS calculations.
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, JLEV(LEV,NQN)
,
is always a pointer to the index of level LEV
in
the cross section array (SIGMA
).
For atom  linear rigid rotor systems, the JLEVEL
array is usually generated from input parameters JMIN,
JMAX
and 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
.
If JLEVEL
is calculated from JMIN,
JMAX
and JSTEP
, the energy levels are
generated from the input parameters BE
(rotational
constant), ALPHAE
(vibrational dependence of
BE
), and DE
(centrifugal distortion
constant). The rotational constant actually used is, of course,
simply BE0.5*ALPHAE
, and this value may be input
directly in BE
with ALPHAE = 0.0
(the
default) if preferred.
If NLEVEL .gt. 0
, the energy levels may be
specified in the input array ELEVEL
. If
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 JLEVEL
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 ELEVEL
,
they are generated from the input parameters WE
(omega(e)), WEXE
(omega(e)x(e)), BE
(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 WE
,
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 ELEVEL
are zero, energies are generated from WE
and
WEXE
. If the input NLEVEL
is negative,
a usersupplied routine GET102
is called to set
NLEVEL
, JLEVEL
and ELEVEL
values.
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.
If NLEVEL = 0
on input, JLEVEL
is
generated from J1MIN, J1MAX, J1STEP, J2MIN, J2MAX
and 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. J1
.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 antisymmetric 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 BE(I),
ALPHAE(I)
and DE(I)
. If IDENT=1
and the values for I=2
are zero, the program sets
them equal to the values for I=1
.
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
for ITYP = 6
which are
described more fully below. If rotor constants (A(1),
B(1)
and C(1)
) are specified, the asymmetric
top wavefunctions will be calculated, subject to limitations in
JMIN, JMAX, JSTEP
(synonyms are J1MIN
,
J1MAX and J1STEP
), EMAX
and
ISYM
. Otherwise, wavefunctions will be input from
unit IASYMU
; if NLEVEL .gt. 0
, only
this many levels will be input, else input to endoffile. 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 BE(2)
. If
specified, ALPHAE(2)
and DE(2)
are
also used as for ITYP = 1
. If EMAX .gt.
0
is specified, it is used to limit basis functions to
energies (asymmetric top plus linear rotor) less than
EMAX
.
In order to allow for finer control of the basis set, specifying NLEVEL
.lt. 0
will try to obtain from IASYMU
the basis functions specified by
(JLEVEL(3*i2) = J1
, JLEVEL(3*i1) = ITAU
, JLEVEL(3*i) = J2)
, i = 1,
abs(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
from 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 NLEVEL
,
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*i2)
= j, JLEVEL(3*i1)
= k,
and JLEVEL
(3*i)=parity, i=i,NLEVEL
.
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 nearsymmetric 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 xz plane be a reflection plane. For
a symmetric top, then, A=B
and C
is
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 IASYMU
.
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 IASYMU
by J, ITAU,
ELVL (2I5,F15.10)
where J
is the rotational
quantum number, ITAU
is an index that is not
actually used, and ELVL
is the rotor energy (in
1/cm unless &BASIS EUNITS
has specified
otherwise). For each level, NK=2*J+1
expansion
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 zaxes in which the potential is described.
If IASYMU
is the standard input unit, data cards
should follow the &BASIS ... &END
set and
precede the &POTL
set; in this case a positive
value of NLEVEL
must also be set in
&BASIS
, to specify the number of expected rotor
functions.
If IASYMU
is an auxiliary unit, it is possible to
read levels until an endoffile condition occurs, and
specifying 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
(if
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
JLEVEL
(2*I1),JLEVEL
(2*I),I=1,abs(NLEVEL
). In this case indexing order in the
SIGMA (cross section) array will correspond to the ordering in JLEVEL
.
Note that for IOS calculations (ITYPE = 106) rotor wavefunctions are used
only to calculate statetostate 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
data on IASYMU
if statetostate 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 IPASYM
):
PRTY  even/odd k  epsilon 

0  even  +1 
1  even  1 
2  odd  +1 
3  odd  1 
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
variable 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
To construct 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
with 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, neardegeneracies can occur, and the
program interprets these as degenerate if the levels (from
diagonalisation) are within a tolerance coded in subroutine
SET6C
(usually 1.D8). To be safe, do not set any
highbit 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
variable 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 Jvalue (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
to
+J1MAX
, and g2 loops from J2MAX
to
+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 energyindependent. 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 JLEVEL
array is read in explicitly.
Subsequently, for particular values of THETA, PHI
and ENERGY
, MOLSCAT calculates the
parallel component of the incident momentum K
, and
selects from the master list those basis functions for which
(0.5*hbar**2/URED)*K+G**2
is less than the input
parameter EMAXK
. This occurs even if the
JLEVEL
array is specified explicitly in
&BASIS
data.
Thus, only those basis functions which satisfy both the
EMAX
and EMAXK
criteria will
ultimately be included in the basis set. Since the
EMAXK
criterion is usually the most sensible
physically, 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 usersupplied 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

JHALF (default 1) is introduced to facilitate
halfoddinteger JTOT values in BASE9 codes.
The following variables are specifically relevant to decoupling approximations.
JZCSMX

(default 1). By default, CS/DLD calculations
are performed for all values of the bodyfixed
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. JZCSMX
are then not valid.
JZCSFL

(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).
IBOUND

(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
for ITYPE
= 21, 22, 23 and 27.
Another effect of IBOUND = 1
is to suppress coupled
states calculations for projection quantum numbers greater than
JTOT
.
IOSNGP

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),
theta(2), delta(phi).
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 &INPUT
,
Section 2.1) and/or input basis
set rotor levels.
IPHIFL

(default 0) controls the type of numerical quadrature on phi for
ITYPE = 103, 105
and 106
. The default
requests GaussMehler (recommended) while IPHIFL .ne.
0
requests GaussLegendre.
Forward to
Section 4
Back to
Section 2
Back to the
Table of Contents