What do MM parameters for the quantum region actually do?

From NWChem

Viewed 2047 times, With a total of 2 Posts
Jump to: navigation, search

Clicked A Few Times
Threads 4
Posts 12
Hi everyone...

I'm trying to run a QM/MM simulation (a NEB calculation to be precise) of a simple SN2 reaction:

CH3CH2Br + OH- --> CH3CH2OH + Br-

When I put together a PDB file for this and tried to generate .top and .rst files for it, NWChem gave me back an "Unresolved atom types" error, because the standard AMBER forcefield didn't know what to do with Br or with the OH- ion, so the prepare module didn't know which atom types to assign to these atoms.

For me, this raised two questions. First, how much do the MM parameters for atoms within the QM region actually matter? Since energy calculations for the QM region are all happening within the QM code and electrostatic interactions happen via ESP or something similar, do I just need to make sure I have reasonable VdW parameters for my quantum region atoms? I really don't want to have to parameterize an MD forcefield in order to do a QM/MM calculation!

Secondly, is there a good (i.e. well-commented) example somewhere of how to construct the necessary parameter and topology files for a non-standard molecule? By non-standard I mean that it contains atom types other than those contained in the amber95.par that is distributed with NWChem. Ideally, I'd like to be able to specify atom types and connectivities manually rather than trying to get the prepare module to read my mind. I've looked around in the $NWCHEM_TOP/QA/tests directory and can't find anything that seems to be a clear example of what I'm looking for.

My PDB file (start.pdb) is:

ATOM      1  C1  UNK B   1       0.002  -1.069  -0.182  1.00  0.00            
ATOM      2  H11 UNK B   1       0.569  -1.226   0.733  1.00  0.00            
ATOM      3  H12 UNK B   1       0.607  -0.628  -0.973  1.00  0.00            
ATOM      4  C2  UNK B   1      -1.355  -0.425   0.012  1.00  0.00            
ATOM      5  H21 UNK B   1      -1.948  -0.990   0.740  1.00  0.00            
ATOM      6  H22 UNK B   1      -1.266   0.612   0.378  1.00  0.00            
ATOM      7  H23 UNK B   1      -1.911  -0.403  -0.933  1.00  0.00            
ATOM      8  O   UNK B   1       2.828   3.409   1.030  1.00  0.00            
ATOM      9  HO  UNK B   1       2.771   3.886   0.118  1.00  0.00            
ATOM     10  BR  UNK B   1      -0.300  -3.162  -0.922  1.00  0.00            
END

my NWChem input deck is:

start prep_start
prepare
  system start
  source start.pdb
  new_top new_seq
  new_rst
  modify atom 1:C1 quantum
  modify atom 1:H11 quantum
  modify atom 1:H12 quantum
  modify atom 1:C2 quantum
  modify atom 1:H21 quantum
  modify atom 1:H22 quantum
  modify atom 1:H23 quantum
  modify atom 1:O quantum
  modify atom 1:HO quantum
  modify atom 1:Br quantum
  solvate box 4.0
  update lists
  ignore 
  write start_ref.rst
  write start_ref.pdb
end
task prepare

and my output is:

 Fragment                              UNK
  num name  type   link cntr  grp pgrp    charge     polarizab
    1  C1   CM        0    1    0    1    1   -0.100000    0.000000
    2  H11  HA        0    0    0    1    1    0.050000    0.000000
    3  H12  HA        0    0    0    1    1    0.050000    0.000000
    4  C2   CT        0    0    0    1    1   -0.150000    0.000000
    5  H21  HC        0    0    0    1    1    0.050000    0.000000
    6  H22  HC        0    0    0    1    1    0.050000    0.000000
    7  H23  HC        0    0    0    1    1    0.050000    0.000000
    8  O              0    0    0    1    1    0.000000    0.000000
    9  HO             0    0    0    1    1    0.000000    0.000000
   10  BR             0    0    0    1    1    0.000000    0.000000
                                           ------------
                         total charge    -0.000000
 Connectivity
       1-  2
       1-  3
       1-  4
       4-  5
       4-  6
       4-  7
       8-  9
 Created fragment                      ./UNK.frg_TMP
 Unresolved atom types in fragment UNK
 **********
 *   0: pre_mkfrg failed 9999
 **********

Any ideas at all would be greatly appreciated.

Thanks!

--craig

Clicked A Few Times
Threads 4
Posts 12
OK... so I got around the "Unresolved atom types" error by generating a custom amber.par:
# !!!NOTE: JUNK DATA!!!
This is the AMBER96 user defined parameter file for NWChem 3.2 and ARGOS 7.0
Electrostatic 1-4 scaling factor     0.833333
Relative dielectric constant     1.000000
Parameters epsilon R*
Atoms
BR     79.90000 0.00000E+00 0.00000E+00                            1 1111111111 
      Q      35 0.00000E+00 0.00000E+00
OQ     16.00000 8.80313E-01 1.72100E-01                            1 1111111111 
      Q       8 4.40157E-01 1.72100E-01
HQ      1.00800 0.00000E+00 0.00000E+00                            1 1111111111 
      Q       1 0.00000E+00 0.00000E+00
CQ     12.01000 4.57729E-01 1.90800E-01                            1 1111111111 
      Q       6 2.28864E-01 1.90800E-01
Cross
Bonds
BR    -CQ      0.224041 0.0
CQ    -CQ      0.15260  0.0
CQ    -HQ      0.10900  0.0
OQ    -HQ      0.103079 0.0
Angles
BR   -CQ   -HQ      1.78568 0.0
BR   -CQ   -CQ      1.89522 0.0
Proper dihedrals
Atom types
BR       35  0    0    0  0      0  0      0      0      0
                                 0  0      0      0      0
                                 0  0      0      0      0
OQ        8  0    0    0  0      0  0      0      0      0
                                 0  0      0      0      0
                                 0  0      0      0      0
HQ        1  0    0    0  0      0  0      0      0      0
                                 0  0      0      0      0
                                 0  0      0      0      0
CQ        6  0    0    0  0      0  0      0      0      0
                                 0  0      0      0      0
                                 0  0      0      0      0
End
#

This at least let me start the calculation, although setting the VdW parameters for Br to zero is probably the worst approximation possible. I also set all the bond and angle spring constants to zero and didn't go out of my way to make sure that the prepare module figures out my bonding properly (it has trouble with the C-Br bond, FWIW). Is this a safe thing to do?

--craig

  • Marat Forum:Admin, Forum:Mod, NWChemDeveloper, bureaucrat, sysop
    Profile
    Send PM
Clicked A Few Times
Threads 2
Posts 43
I have added section on how to prepare QM/MM calculations from scratch in http://www.nwchem-sw.org/index.php/QMMM_Appendix#Preparing_QM.2FMM_calculations_from_scrat...

Marat


Forum >> NWChem's corner >> QMMM



Who's here now Members 0 Guests 1 Bots/Crawler 0


AWC's: 2.5.10 MediaWiki - Stand Alone Forum Extension
Forum theme style by: AWC