$QM Command for SIGMA

$QM command for SIGMA
****************************

BOUND [PSEUDO | LINK] <newline>
SMEAR [ON/YES | OFF/NO] [FACTorsmear_factor] <newline>
QMTEMPlate<newline>
QM_Template_File_Name<newline>
QMJOBNAME <newline>
QM_Job_File_Name<newline>
BASISFILE<newline>
QM-basis-filename <newline>
CHITEMPlate<newline>
CHI_Template_File_Name<newline>
CHIJOBname<newline>
CHI_Job_File_Name<newline>
KAPPATEMPlate<newline>
KAPPA_Template_File_Name<newline>
KAPPAJOBname<newline>
QM_Job_File_Name<newline>
ECPFILE <newline>
ECP_File_Name<newline>

First of all, this is an extended version of QM command. It allows user to choose either the SCCDFTB program naturally built in Sigma or Gaussian 03 programs through a "system" call. If the QM calculation is very expensive, the user may wish to use shared-memory version of Gaussian or Linda. For both the shared-memory and Linda version of Gaussian job, one may specify number of processors in the "Qm_Template_File"; however, for the Linda job, one has to specify the "nodelist" before starting Sigma simulation job. Gaussian web site has an instruction of how to setup Linda job.

For a reference to the use of SCC-DFTB in Sigma, please see the following reference: Hao Hu, Marcus Elstner and Jan Hermans, Comparison of a QM/MM force field and molecular mechanics force fields in simulations of alanine and glycine "dipeptides", Proteins: Structure, Function and Genetics, 50: 451-463 (2003)

The keywords "PROG" and "METHOD" are straightforward.

The following keywords are effective for both QM programs:
SUBSystem specifies a set of atoms that will be treated by QM. The set starts from "atm1", and ends at "atm2" (both inclusive). You may use multiple "SUBS" keyword to specify different QM regions. (Since the specification is dependent on the order of atoms in Sigma, you may wish to use "$tree" command to examine the atomic orders before you start QM calculation, unless you are an experienced Sigma user/developer)
USEQMCHARGEspecifies that we will use QM charges in MD simulation. Considering the multiple-time-step algorithm employed in Sigma, one may use averaged QM charges instead of instantaneous QM charges of one particular conformation. The number of steps to average QM charges can be specified by "STEP". A "Qm_Charge_Step" must always be given.
USEPOINTCHARGEis a special keyword, and shall make no difference for SCCDFTB method. In Gaussian job, the QM-MM electrostatic forces acting on MM atoms can be calculated by two means: we can compute the electrostatic potential from QM atoms at certain (MM) positions, and hence the force on MM atoms; or, we can compute the force by using ESP point charges for QM atoms. The former is most accurate (for example, it obeys Newton's third law); however, in some circumstances such as the user doesn't want to do QM calculation every MD step (QM is too expensive), then he may wish to use the ESP charges to compute electrostatic interaction.
USEPOLARCHARGEis a special keyword for using polarized QM charges. In QM/MM calculations, we sometimes decompose the full QM/MM electrostatic interaction into two parts: Eele(QM/MM) corresponding to the electrostatic interaction between QM and MM part which is often approximated by the point charge interaction between MM charges and QM ESP charges, and E'(QM) which represents the intrinsic energy of (polarized) QM subsystem. The Eele(QM/MM) term is calculated as the sum of Columbic interactions between MM point charges and QM ESP charges when USEPOINTCHARGEis specified, which is really an approximation. The accuracy of Eele(QM/MM) term can be improved to allow the QM ESP charges to be polarized by the change of MM point charges which is achieved in this USEPOLARCHARGEkeyword. When USEPOLARCHARGEis specified, one also needs to specifieCHITEMPlate, CHIJOBname,KAPPATEMPlate, andKAPPAJOBname for calculating two kernels: CHI and KAPPA. Therefore, USEPOINTCHARGEand USEPOLARCHARGE are exclusive.
RECORDCHARge: write out atomic charges every so many steps. (IfFILE is specified, the list will be in a separate file, whose filename will be based on job file name. Else the list will be part of the standard output from the sigma program.)
QMFREQis a special keyword and shall NOT be used in normal simulations. The user must understand the consequence of this command. This keyword specifies that we will only do one QM (force and energy) calculation for every "Qm_Cal_Step" MD steps. Therefore, it only makes sense when the QM subsystems are frozen during those "Qm_Cal_Step" MD steps of simulation period. This command then shall only be used with combination of "$neb" command.
FULLMMis a special keyword and shall NOT be used in normal simulations. By default, SIGMA generates the list of MM atoms based on the cutoff specified in the job file (by "$md" or "$setup"); if "FULLMM" is activated, then the charges of all MM atoms will be used in the QM calculations as background charges.

The following keywords are mostly related to Gaussian programs:
BOUND specifies how to treat the QM-MM crossing bond. Right now we have only one choice as "pseudo-bond" approach developed by Zhang et al.. When one uses pseudo bond approach, one has also to specify "BASISFILE" which stores information of pseudo-bond parameters and will append to "QMJOBNAME".
SMEAR specifies how to smear the MM charges that are close to QM atoms. For those MM atoms that have covalent interaction with QM atoms, i.e., 1-2, 1-3, and 1-4 interactions, their charges shall be smeared (by a smear_factor) in QM claulcations, otherwise, the strong electrostatic interaction from near MM point charges will make extra difficulties for QM SCF iterations.
QMTEMPLATE specifies a file that will be used as template for building a Gaussian input file. Usually one will put everything of a normal Gaussian input file before the coordinate part into that file, the Sigma program will append coordinates after that. (Remember, you have to calculate and specify charges and multiplicity of your QM system there!!!)
QMJOBNAME specifies the name of Gaussian input file that will be created by Sigma. This keyword is important when one has several Gaussian jobs concurrently running in the same directory which one needs to specify different names for each job.
BASISFILE specifies a file which stores information about basis set and pseudo-bond parameters. Its contents will be appended into Gaussian input file. So one usually won't need this if the QM system is isolated.
ECPFILE specifies a file which stores information about effective core potentials (ECP). Its contents will be appended into Gaussian input file. So one usually won't need this if the the basis set is already specified in the QMTEMPLATE file.

The following keywords are mostly related to SCCDFTB program:
BASISFILE:The file QM-basis-filename is a helper file that describes the atom types and gives the filenames where the basis sets are given. The file contains:
n = number of atom types <newline>
n* [name = 1 or 2-char chemical name of the atom type <newline>]
n integers = number of shells for each atom type <newline>
n*n [basis-set-filename = description of interaction between atom types i and j <newline>]

The order in which the files containing the tables are given must correspond to the order of the atom types. E.g., if the atom types are given in the order O, C and H, then the files must be given in the order
table for OO
table for OC
table for OH
table for CO
table for CC
table for CH
table for HO
table for HC
table for HH

A library of basis-set-files is provided in a directory $SIGMA/QMDATA/SLKO. E.g. the file co.slk is for interaction of C with O. (It has a different content than the file oc.slk.).
A basis-set-filename can be given as (i) absolute path (begins with /) (ii) path relative to the location of QM-basis-file. (iii) absolute path beginning with a sigma macro (begins with %).
SUBSystem: Select the range of atoms (atm1 through atm2) that is to be treated with QM
SCFTOL:Specify tolerance for the SCF calculation
SCFDEBUG: Control debug messages.
CHARGE: PSF = charges of MM atoms are found in the psf file
VDW: ON,DAMPed = include damped intramolecular Lennard-Jones 1/r6 term, i.e. only at distance exceeding the van der Waals distance; FULL = include full intramolecular LJ6 term.
DETOL: Specify the shift applied to the atomic coordinates to estimate the gradient (as E/x). Default value is .001.

Example:

$QM
methodscctb
subsystem 1 ha1 3 ha3 # alanine dipeptide
scftol 2.e-9
basisfile # for O, N, C and H
%(SIGMA)/QMDATA/basis_onch.dat
chargepsf
vdw on
scfdebug off
recordcharge 8
end

RETURN TO TABLE OF COMMANDS