Hi Saransh
When Jason Lambert and I added basic OPLSAA support to moltemplate, I did not know how to add automatic atom-typing capability to moltemplate. I assumed that it was a relatively simple procedure to look through the force field file (eg “oplsaa.lt” or “oplsaa.prm”) and choose the correct atom type for every atom in the molecule. This was naive of me. I mostly run coarse-grained simulations.
Since then I feel some measure of guilt. I have created a tool that is incomplete. A tool that offers the illusion of making it easy to set up a simulation, but is actually tool that makes it easy to set up an incorrect simulation. I’m not sure if this is a service or a disservice to the community. (The good news is that, for OPLSAA at least, when you select the wrong atom types, the charge will usually be non-neutral, so you have some warning before you start the simulation.)
Again I’m not an expert on how to run fully atomistic simulations. But I’ll share with you what I know, and make some suggestions:
I used Moltemplate to generate LAMMPS data file using OPLS-AA force-field parameters for deltamethrin (DM) & flubendiamide (FM), both are charge neutral.
However, there is no charge neutrality in the data files generated (~+0.5 for DM & 4.21 for FM). Since, all the atom types present in the molecule are not in OPLS-AA library, I made some approximations while assigning atom-types, which most likely led to this.
please carefully study the OPLS papers to understand the way how (partial) charges are applied. there should be no need to make approximations. OPLS-AA uses an increment system, and thus charges should always add up to zero for neutral molecules or integer charges for charged molecules.
I’ve been confused about this too. It wasn’t clear to me that there are any rules defining how the partial charges were applied. In the original 1996, the partial charges listed in a table according to atom type (See Jorgensen et al 1996, Table3, supplemental.) Now it seems they just use arbitrary PDB templates to assign charges and atom types. (See details below.)
— PDB templates in the official OPLSAA files —
If my understanding is correct, then current version of OPLSAA distributed directly by the Jorgensen group includes a file containing a list of PDB "topology templates". (These old force fields assume everybody is still using PDB files.) For each atom in your PDB file, the atom type string and residue name string in the PDB file is compared with the available atom type and residue type names names in this topology template file. If the topology file contains a definition of the kind of molecule you want to simulate (ie, if there is a matching residue-type string in the topology file), then the atom type and partial charged are read from the corresponding entry in the table. But if the molecule you are interested in does not appear in the PDB topology template file, then you might be out of luck. (The file I referred to above contains mostly amino acids.)
(Perhaps not. I found a “par_opls_ligand.inp” file containing additional atom types and force field parameters for small organic molecules of biological significance, “ligands”. I think it is distributed with CGENFF. Presumably there is another file somewhere containing PDB topology templates for these small ligand molecules which CGENFF uses, but I couldn’t find it. See below for more information on PSFgen and CGENFF.)
— PDB templates in the TINKER .PRM files —
In contrast, as you know, moltemplate currently gets its OPLSAA parameters from the “oplsaa.prm” file which is distributed with TINKER. This version of OPLSAA defines a large number of atom types (906 atom types). Each atom type in this file has a well defined partial charge. (For comparison, the version of OPLSAA distributed by the Jorgensen group is significantly smaller, containing only 81 atom types.) I’m guessing the TINKER developers created a new force field with 906 atom types by scanning through all of the official OPLSAA topology template files, and considering every combination of partial charges in those templates that can be used with each atom type in the original OPLSAA parameter files. At the end of the “oplsaa.prm” and “oplsaa.lt” files, there is a list of templates indicating the atom type (and consequently, its partial charge) for all of the atoms in the common amino acids. (Search for “biotype” at the end of the file.) If you are not simulating one of those, I think you are probably out of luck as well.
— where to find OPLS parameters (and/or templates) —
In general, it’s not so easy:
https://www.researchgate.net/post/Where_can_I_find_parameters_for_all_elements_for_the_OPLS_Force_field
- Try using LigParGen. This is probably the easiest automatic way to build molecules for LAMMPS using the OPLSAA force field, especially if you plan to use moltemplate. It works with a relatively wide range of molecule types. It will assign atom partial charges by performing a CM1A (quantum mechanical) calculation, and it (currently) assigns atom types automatically. There is no way to understand the atom type names in the molecules that it generates. They don’t agree with the atom types distributed by the Jorgensen group, or the Ponder(TINKER) group. However it probably doesn’t matter because all of the force field parameters are included locally in the molecule definition. There’s no need to use an external force-field (such as “oplsaa.lt”. (Incidentally, this is somewhat similar to the approach that ATB uses, although I think the ATB approach sounds a little bit more sophisticated.) Incidentally, I briefly demonstrated how to use LibParGen on slide #191 of this talk.
If you are using LigParGen with MOLTEMPLATE, then you must convert the LAMMPS data file that LigParGen creates into moltemplate format using the “ltemplify.py” program (distributed with moltemplate). Documentation and examples for “ltemplify.py” can be found here and here (and on slide 161 of this talk.)
Unfortunately LigParGen is a “Black Box”. There have been some complaints that it is not working correctly, but the developer is not responding to them and appears to have moved on.
- Try using “EMC”
It’s not open-source, but it is free, and it seems like a reliable and useful tool. Apparently, EMC has a more clever way to infer atom types based on SMILES strings. It might be able to assign atom types to a wider variety of molecule types than PSFgen can. I don’t know how (or if) EMC assigns partial charges to atoms.
- Try using PSFgen. (Documentation located here. Note: PSFgen uses the official OPLSAA files along with a PDB file which you supply. If my understanding is correct, PSFgen creates a PSF file containing the correct OPLSAA atom names and partial charges for the atoms in your PDB file. I’ve never used PSFgen, but my understanding is that you must supply it with the right PDB topology files containing a definition of your molecule. PSFgen just reads your PDB file and looks up the atom type and charge from entries in the PDB topology files that you supply to it. (See above.) Once you have used PSFgen, you must copy the atom partial charges into your moltemplate file and convert the atom type strings into the corresponding numbers used by the TINKER version of OPLSAA. (In the future, I will probably write a short script to do this automatically so that you don’t have to copy this information by hand. In the future, I’d like to have PSFgen and VMD and MOLTEMPLATE working well together. This will happen some time next year.)
Let me know if you find a better way.
I hope this helps.
Andrew
I can try to make better approximations but still there would be some charge left. So what is the alternative to maintain charge neutrality? Is there some way around it rather than deriving force-field parameters?
you have to make certain, that you follow the force field rules.
One more point, on generating data file for a fairly simple molecule (octanol, DMSO), there wasn’t charge neutrality in spite of all atom-types defined in the OPLS library. What could be the reason for this?
you must be incorrectly assigning charges with your force field parameters. you can easily check this manually e.g. for DMSO you have atom types 437 (for alkyl-sulfoxide sulfur), 438 (for alkyl-sulfoxide oxygen), 439 (for sulfoxide methyl carbon) and 85 (for methyl hydrogen). that results in partial charges of 0.13, -0.42, 2x -0.035 and 6x 0.06 = 0.13-0.42-0.07+0.36 = 0
axel.
Sometimes it’s not so hard.
However I have run into examples of molecules (eg, “PAMAM, generation1”) where there were some atoms which just didn’t seem to fit any of the descriptions in the “oplsaa.prm” file (or “oplsaa.lt” file), and the charges did not cancel out either.