Computational Biophysics & Bioengineering Section
Department of Physics, Tougaloo College, MS
|gmx-3.3.1_qmmm-1.3.2||Bug-fixed (related to MM-layering) version of gmx-3.3.1_qmmm-1.3.1||Last Updated: August 1, 2008|
|CPMD-3.11.1||Link to the official CPMD source code (free) for QM/MM with Gromacs-3.3.1||Version Released: March 2006|
|QM/MM Examples||Examples on water-dimer, ethane, BH4-TIP4P, Imidazole-CO, and Heme-Imodazole-CO||Posted: July 19, 2006|
Brief introduction to the Gromacs-CPMD QMMM Interface
The interface can be used to perform Energy Minimization, Molecular Dynamics Simulation, Proton Transfer Reaction, Biased Proton Transfer Reaction with Gromacs Umbrella Sampling Potentials, Transition State Search, etc. Further information is given at the end of the Installation Guide section
NOTE: The interface is meant for hybridization with all atom force-fields. For Gromos united atom ff, though it is conjectured that it can be used with some manipulation to the QM atom definitions (see united-ff item below), it has never been pursued.
How to Configure, make and install the above Gromacs versions so as to enable the interface with CPMD
The QMMM compatible Gromacs package provided above and the CPMD package need to be installed separately in your system. The configuration of gmx-3.3.1_qmmm-1.3.2 is like that in Gromacs but needs the qmmm flag: '--with-qmmm-cpmd'. The configure flag '--with-qmmm-cpmd' activates CPMD related codes in Gromacs. To install, entering your gmx-3.3.1_qmmm-1.3.2 folder use the following:
configure --with-qmmm-cpmd [you can use additional configure flags as you need] make make install
How to Compile CPMD
Follow usual compilation steps as in CPMD. Make sure that you DO NOT USE any QMMM option because the qmmm option in CPMD is for interfacing with GROMOS.
How to start a QMMM run
Using 'make_ndx' create index.ndx file and add your [QM] group there (see supplied examples). Note: QM atoms id's should be in accordance with the .gro file of Gromacs.
Make sure your working directory, where you have your .pdb or .gro or .top, contains the files:
- 'index.ndx' (as mentioned above)
- 'CPMD_inp.tmpl', (supplied; may need to be edited for new CPMD atom/Pseudo-Potential types)
- 'runcpmd' (supplied; edit to suit your CPMD run) and
- all the CPMD pseudopotential files which CPMD will look for during runtime.
Now edit the mdp file for grompp instructions for CPMD (see rgmx script from example).
Finally, to remove the QM atoms from the SHAKE list of Gromacs (which constraints the atoms to some fixed bond lengths), set the BONDTYPE to 5 in the topology files. This can be done by opening the appropriate file and changing the third column in the bonds section to 5 for all atoms which are part of the QM system.
You are now all set for qmmm. Execute the rgmx script; grompp will setup the qmmm environment and mdrun will invoke CPMD once cpmd executable cpmd.x is accessible from your current working directory.
How to choose link atoms
Once you choose your [QM], link atoms are automatically set by the interface code. Right now only cut accross C-C bond is supported by the code. If you need to cut accross other bond type and add additional link options, you may contact us.
How Gromacs accesses CPMD and QMMM hybrdizaion works
First, grompp sets the QMMM environment through the mdp file parameters (see in the rgmx script). When mdrun is invoked, the qmmm interface call creates CPMD input file CPMD_inp.run, an empty file QMSTART and invokes cpmd through the command in the 'runcpmd' file (which should be present in your working directory. If you are using a parallel version for CPMD, runcpmd file must contain parallel run commands (see supplied runcpmd file). Once CPMD finishes its job, it writes QM output in qmoutput.out and creates an empty file GMXCONTINUE and sleeps. Gromacs continues with GMXCONTINUE. The process repeats.
How the MM info goes to QM
Via the &EXTELEC section of 'CPMD_inp.run'. The first part of 'CPMD_inp.run' contains the QM information including the link atoms, if any. This file and the &EXTELEC section is created automatically during 'mdrun' and is updated for the coordinates at each step of simulation/energy-minimization.
How to provide parameters and instructions to control the CPMD run
The job is divided into two parts. Some of the basic parameters like QM simulation box-size, planewave-cutoff, etc are fed through the mdp file of gromacs (please see the 'rgmx' scrpit supplied in the examples). For other standard CPMD specifications and the type of Pseudo Potentials, Total Charge, Spin, Multiplicity, etc, are fed through a template file called CPMD_inp.tmpl (please see in the worked out examples).
How to use Gromos united atom force field in Gromacs-CPMD QMMM
It is conjectured that the missing non-polar hydrogen for the QM may be addressed by: i) duplicating the residue in the rtp file with a new residue type (which will be treated by QM) where the non-polar hydrogen atoms are replaced with a new polar hydrogen type. ii) and by setting the LJ coefficients of this new hydrogen type to zero. If someone interested to test pl. contact us.
Further information about this interface
The interface uses a s-wave partial wave expansion for MM atom point charges to avoid non-physical polarization of QM by MM. Further development with a s3p type expansion of MM atom point charges that would include charge-charge, charge-dipole, and charge-quadrupole interactions for H_QM/MM will soon be submitted to CPMD and is expected to be available accordingly.
The interface divides the MM atoms into three layers with respect to each QM atom where the inner layer contains all the residues within the chosen inner_layer_radius and all the POLAR residues and IONS disregarding their distances from QM atoms. Thus, a cutoff of QMMM interactions with a inner_layer radius of (10-15 Angstroms) should provide a reasonable polarization effect on the QM by the MM considering that the polarizing effect of the non-polar residues in the intermediate and outer layers will not be appreciable as the Charge(QM)-Dipole(MM) interaction (1/r^2) decays faster than Coulomb (1/r). However, if you try to include the effects of the intermediate and outer layers through some updating frequencies (other than one; which means including full QMMM interaction) please check the conservation of energy whether it is satisfactory for your system.
In its current form, gmx-3.3.1_qmmm-1.3.2 contains Gromacs interface to other QM codes as well (developed by Dr. Gerrit Groenhof). To access the CPMD, it is important that one uses the flag '--with-qmmm-cpmd' while configure (as mentioned above). Details of QM/MM procedures are given in the manual inside qmmm_project and is also available in the manual of CPMD-3.11.1 distribution.
Gromacs Interface: The version gmx-3.3.1_qmmm-1.3.1 had a typo in the function set_MM_layers() which was causing the MM layers to be decided based on the first QM atom instead of all the QM atoms. This bug reported by Marius Retegan has been fixed in the version gmx-3.3.1_qmmm-1.3.2 (last updated on August 1, 2008).
CPMD: If you are interested to consider the effect of intermediate and outer MM layers on QM, you need the following fixes in CPMD-3.11.1.
- In line 784 of egointer.F, change the index of QINTM(2,I) to (3,I)
- In line 808 of egointer.F, change the index of QOUTR(2,I) to (3,I)
- In line 249 of eextern.F, change the index of QINTM(2,I) to (3,I)
- In line 290 of eextern.F, change the index of QOUTR(2,I) to (3,I)
A likely cause for this is you are attempting to run CPMD on more than one processor, but you are using the serial compiled version of cpmd.x. If you have CPMD parallel version compiled, check the path to this cpmd.x, and if not, recompile CPMD with mpi and run again.
CPMD starts but stops with: "CONTROL| INPUT FILE NOT FOUND ...PROGRAM STOPS IN SUBROUTINE CONTROL 999"
Code Development and Maintainance: The code was developed by P. K. Biswas while at Gogonea group with financial suppoort from DoE: grant# DE-FG02-03ER15462. The QMMM code is now maintained by the Computational Biophysics & Bioengineering group at Tougaloo College, headed by Dr. P. K. Biswas, Assistant Professor of Physics.
Developmental Support: Sridhar Kalubandi, Kiran Tatapudi, and Robert Fiske; Gogonea group, Cleveland State University, Ohio.
Improvement & Maintainance Support:
- Dr. Axel Kohlmeyer, Center for Molecular Modeling, University of Pennsylvania, USA.
- Mr. Christian Seifert, Department of Biophysics, University of Bochum, Germany.
- Dr. Gerrit Groenhof, Max Planck Institute for Biophysical Chemistry, Germany.
- Mr. Marius Retegan, Departement de Chimie Moleculaire, France.
- Dr. Samantha Jenkins, University West, Trollhättan, Sweden.
Webpage Development: Aritra Biswas