[CPMD-list] how to include water in QM part in a QMMM calculation

Wed May 16 22:43:20 CEST 2007

Dear Axel,

On 5/16/07, Axel Kohlmeyer <akohlmey at cmm.chem.upenn.edu> wrote:
>
> On Wed, 16 May 2007, Yong Zhang wrote:
>
>
> dear yong,
>
> YZ> The system was setup initially using GROMOS. I generated the solute
> YZ> (including NMA and 79 solute water molecules), then solvent was added
> YZ> (solvent water).
> YZ> The whole system was equilibrated in cpmd/gromos using "molecular
> dynamics
> YZ> classical" with solute fixed. This is only a test run, so I suppose it
> is
> YZ> ok.
> YZ>
> YZ> The gromos input was modifed from the one in your example (from
> biocore), in
> YZ> which case amber was used. What is difference between amber and gromos
> form
> YZ> in this QM/MM calcualtion? How can I modify that to a gromos job?
>
> the major difference is the force field. if you set up
> a job with gromos certain sections have 'predefined'
> values in them. e.g.:
> ATOMTYPENAME
> BONDTYPE
> BONDANGLETYPE
> IMPDIHEDRALTYPE
> DIHEDRALTYPE
> LJPARAMETERS
>
> if you prepare a job with amber the contents of
> those are generated from the amber topology files
> from amber7togromos.x

I do not have amber here, but I guess amber7togromos.x converts the amber
TOPOLOGY file into gromos topology, not the input file. (or maybe it does?).
But what is the difference in the gromos input file for amber ff or gromos
ff?

[...]
>
> YZ> I added "flexible water" and the problem gone now. But I still do
> YZ> not understand: what is the "solute water" for in gromos topology?
> YZ> or in other word, under what situation should I use that?
>
> sorry, i cannot follow what you did and i don't really
> think that you've solved your problem (you only avoided the
> symptom, which is worse IMO).
>
> if you look at the gromos input file you have in the SYSTEM
> section NPM and NSM. the first is the number of solute 'molecules'
> (which can consist of more than one molecule) and the solvent.
> the geometry of the solvent is determined by the sections
> SOLVENTATOM and SOLVENTCONSTR in the gromos topology.
> the QM/MM implementation is done in a way, that those
> molecules/atoms _cannot_ be turned into QM atoms. if you
> want them to become quantum, you have to rearrange the
> coordinates in the .crd file, add/change the proper definitions
> to the various sections in the gromos topology and input files.
>
> the FLEXIBLE WATER flag can assist in that conversion, but
> you have to check the resulting topology file carefully
> and adjust the rest of the input.

I realized that solvent atoms can not be treatted in QM section. In that
case, a error message was given clearly demonstrated that. So in this case,
I include some water molecules in the SOLUTE part. I did this by setting
NMOL=80 when builtup the topology using gromos:
============
#number of building blocks, heme, separate molecules, chains
#NRAA   NHEME    NMOL   NRAA0
  0       0        80      0
#
#building block names
#AANM[1..NRAA+NHEME+NMOL]
NMA H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
H2O H2O H2O H2O H2O H2O H2O H2O H2O H2O
==============

More solvent water were added to this "80 molecules" solute. Then in gromos
input for QM/MM run, I have NPM=1 and NSM=429. The NPM=1 here actually
includes one NMA molecules and 79 water molecules which is indicated in the
SUBMOLECULES section in gromos input as
==============
SUBMOLECULES
# NSPM > 0 number of (sub)molecules the solute consists of
# NSP(1.. NSPM) > 0 atom sequence number in topology
# of last atom in submolecule I
#      NSPM  NSP(1.. NSPM)
       80      12 15 18 21 24 27 30 33 36 39 42
45 48 51 54 57 60 63 66 69 72 75 78 81 84
87 90 93 96 99 102 105 108 111 114 117 120
123 126 129 132 135 138 141 144 147 150 153 156
159 162 165 168 171 174 177 180 183 186 189
192 195 198 201 204 207 210 213 216 219 222
225 228 231 234 237 240 243 246 249
END
===============

By doing this way, when added one water into QM part, I did not get the
error message saying that "solvent can not be in QM". But, unfortunately, I
got another error related to the "shake" problem.

Below I attached the two definition of water in gromos, one as solute and
one as solvent. And I used the solute one in my "solute"(the 79) and solvent
on in my "solvent" (429).
But I guess, as you indicated, I do not really need so many complex steps.
I can simply setup a real solute (NMA in my case) and added solvents to it,
then use "flexible water" to convert some solvent into solute and make
necessary changes in the crd file. Then it should work. I will try this.

===================solute one
MTBUILDBLSOLUTE
# building block (residue, nucleotide, etc.)
# RNME
H2O
# number of atoms, number of preceding exclusions
# NMAT,NLIN
    3    0
# preceding exclusions
#ATOM                               MAE MSAE
# atoms
#ATOM ANM  IACM MASS        CGMICGM MAE MSAE
    1 OW      4   16   -0.82000   0   2    2    3
    2 HW1    18    1    0.41000   0   1    3
    3 HW2    18    1    0.41000   1   0
# trailing atoms
#ATOM ANM  IACM MASS        CGMICGM
# bonds
#  NB
    3
#  IB   JB  MCB
    1    2   35
    1    3   35
    2    3   41
# bond angles
# NBA
    0
#  IB   JB   KB  MCB
# improper dihedrals
# NIDA
    0
#  IB   JB   KB   LB  MCB
# dihedrals
# NDA
    0
#  IB   JB   KB   LB  MCB
END
=======================

===================solvent one
MTBUILDBLSOLVENT
#solvent name
#RNMES
H2O
#number of atoms
    3
#atoms
#ATOM ANM   IAC MASS    CG
    1 OW      4   16  -.82
    2 HW1    18    1   .41
    3 HW2    18    1   .41
#constraints
#number
    3
#  IB   JB  LENGTH
    1    2    0.1
    1    3    0.1
    2    3    0.163299
END
=====================

YZ> One more question here, what was changed by "flexible water"? I checked
> the
> YZ> two topology files, but did not realize any difference, like the bond
> YZ> between two water H  atoms are still there.
>
> this is impossible to tell without seeing/trying out
> all the inputs. the documentation in the manual gives
> you some (minimal) hints, but the rest can only be
> found out by reading the code.
>
> it would be _very_ much appreciated when people who actually
> _use_ the cpmd/gromos qm/mm code regularly (i don't) would
> expand that part of the documentation and explain the motivation
> and proper protocols for setting up a QM/MM input. most of
> what i put in there quit a while ago is guessed from what
> the sourcecode does.

I really appreciate all your effort to make QM/MM accessible to more people.
I am sorry I did not help much so far since I am still hesitating at the
very initial steps. But anyway, hope this discussion is somewhat helpful for
others.

Best wishes,
Yong
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://cpmd.org/pipermail/cpmd-list/attachments/20070516/215cd455/attachment-0001.html 