Differences

This shows you the differences between two versions of the page.

--- howto:biochem_qmmm [2017/08/02 07:06] – dvanrompaey
+++ howto:biochem_qmmm [2019/08/22 09:23] – sllabres
@@ Line 93: / Line 93: @@
 </code>
-Now that the QM atoms are defined, CP2K needs to know how to treat these. We will use the semi-empirical AM1 method as our QM method, electrostatically embedded in the MM system: ''[[inp>FORCE_EVAL/QMMM#E_COUPL]]'' COULOMB. This allows for the QM region to be polarized by the MM environment. Mechanical embedding, where the QM region only interacts with the MM region as point charges and no polarization occurs, can be selected through ''E_COUPL NONE''. In DFT calculations, the highly efficient [[doi>10.1021/ct050123f |GEEP]] method can be used as well. For stability reasons, we need to make some slight modifications to our prmtop file. In the TIP3P water model, the hydrogens do not have a separate Lennard-Jones site. The same is true for the amber serine, threonine and tyrosine HO atoms. This can cause unphysical interactions with the QM system: the hydrogen atoms strongly interact with the QM electron clouds, crashing the system. In order to prevent this we will add LJ sites to water hydrogens and serine/threonine/tyrosine hydrogens with the ambertools utility parmed. We will change their LJ parameters to those of a GAFF2 type alcohol using the ''changeLJSingleType'' command. Take care when using this, as it changes the LJ type for all of the atoms that share the LJ type, not just the one you selected. For instance, changing the serine HO also alters all the threonine and tyrosine HO atoms in the system. In this case, this behaviour is useful, but it is strongly advised to use ''printLJTypes'' to get a good idea of what exactly you are modifying. For a more selective manner of doing this, the ''addLJType'' command is useful. This adds an LJ type for your selected atom or group of atoms, rather than modifying all occurrences of a specific LJ type. Confirm that ''@5724'' corresponds to a water hydrogen and that ''@4843'' is a serine HO using ''printDetails @(atomnumber)'' in parmed. If not, use ''printDetails'' to locate the correct atom numbers and modify the changeLJSingleType commands accordingly.
+Now that the QM atoms are defined, CP2K needs to know how to treat these. We will use the semi-empirical AM1 method as our QM method, electrostatically embedded in the MM system: ''[[inp>FORCE_EVAL/QMMM#E_COUPL]]'' COULOMB. This allows for the QM region to be polarized by the MM environment. Mechanical embedding, where the QM region only interacts with the MM region as point charges and no polarization occurs, can be selected through ''E_COUPL NONE''. In DFT calculations, the highly efficient [[doi>10.1021/ct050123f |GEEP]] method can be used as well.
-<code>
+Before running the QM/MM simulation we need to amend our prmtop file, because in the classical forcefield, several hydrogen atom types do not have separate Lennard-Jones parameters or they are set up to 0.0. The AMBER FF atom types are HO (from the TIP3P water model), HG (from serine, threonine and tyrosine residues). This will cause unphysical interactions with the QM region and the simulation will fail due for stability reasons. The hydrogen atoms with these atom types will strongly interact with the QM electron clouds, eventually crashing the system. To prevent this from happening, we have to add the correct LJ parameters to the mentioned atom types using the ambertools utility parmed.
-parmed complex.prmtop
-changeLJSingleType @5724 0.3019 0.047
-changeLJSingleType @4843 0.3019 0.047
+Parmed reads the prmtop file and allows you to modify the parameters defined in it. The parmed command we are going to use is: ''changeLJSingleType'', which changes the LJ type for all of the atoms that share the LJ type, or ''addLJType'', for more selective modifications. We will change their LJ parameters to those of a GAFF2 type alcohol that has ''0.3019 0.047'' values for the LJ type.
+Either if we are using ''changeLJSingleType'' or ''addLJType'', we need to ensure that we are changing the right atom index. To do so, we need to use ''printDetails'' or ''printLJTypes'' that print the information on the selected atom index.
+<code>
+$ parmed complex.prmtop
+changeLJSingleType :WAT@H1 0.3019 0.047
+changeLJSingleType :*@HO 0.3019 0.047
 outparm complex_LJ_mod.prmtop
 quit