From NWChem
Viewed 2047 times, With a total of 2 Posts
|
Clicked A Few Times
Threads 4
Posts 12
|
|
1:18:01 AM PST - Mon, Nov 28th 2011 |
|
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
|
|
7:26:06 PM PST - Mon, Nov 28th 2011 |
|
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
|
|
Clicked A Few Times
Threads 2
Posts 43
|
|
|
|
AWC's:
2.5.10 MediaWiki - Stand Alone Forum Extension
Forum theme style by: AWC