alanwilter / acpype

OFFICIAL: AnteChamber PYthon Parser interfacE
https://alanwilter.github.io/acpype/
GNU General Public License v3.0
201 stars 46 forks source link

Adding atoms on top of existing atoms #25

Closed mquevill closed 3 years ago

mquevill commented 3 years ago

I have had this issue with a few molecules that I have tried passing into ACPYPE. Even a small test case generates these errors. The

Generate input file with obabel -:"[Li+]" -o mol2 -O li.mol2 --gen3D. The contents are:

@<TRIPOS>MOLECULE
*****
 1 0 0 0 0
SMALL
GASTEIGER

@<TRIPOS>ATOM
      1 LI          1.0009    0.0999    0.0039 Li      1  UNL1        0.0000
@<TRIPOS>UNITY_ATOM_ATTR
1 1
charge 1
@<TRIPOS>BOND

When I try to run ACPYPE via acpype -i li.mol2 -n 1 (since +1 charge), it gives the following output:

=========================================================================================
| ACPYPE: AnteChamber PYthon Parser interfacE v. 2020-10-24T12:16:34CEST (c) 2021 AWSdS |
=========================================================================================
ERROR: more than one residue detected '{'MOL', 'UNL'}'
ERROR: verify your input file 'li.mol2'. Aborting ...

By simply renaming the residue in the li.mol2 file, it gets a little farther. However, then I get the following:

=========================================================================================
| ACPYPE: AnteChamber PYthon Parser interfacE v. 2020-10-24T12:16:34CEST (c) 2021 AWSdS |
=========================================================================================
ERROR: Atoms with same coordinates in 'li.mol2'!
++++++++++start_quote+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
ATOM 1 Li 1.001 0.100 0.004
ATOM 2 Li1 1.001 0.100 0.004
ATOM 3 Li2 1.001 0.100 0.004
ATOM 4 Li3 1.001 0.100 0.004
++++++++++end_quote+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
ERROR: Use '-f' option if you want to proceed anyway. Aborting ...

This looks like it made it into a Li4 complex, but all the atoms are in the same location. With some other cases I was trying, it added 3 H atoms on top of an existing one, or added some F atoms on top of an existing one.

mquevill commented 3 years ago

I am using the versions from Conda when created with the following command: conda create -n acpype -c ambermd -c conda-forge acpype ambertools openbabel

Relevant packages:

# packages in environment at ~/.conda/envs/acpype:
#
# Name                    Version                   Build  Channel
acpype                    2020.10.24.12.16            py_0    conda-forge
ambertools                19.0                          0    ambermd
openbabel                 3.1.1            py39hc681964_1    conda-forge
python                    3.9.2           hffdb5ce_0_cpython    conda-forge

Let me know if there's any other information I can provide to assist in debugging this issue.

mquevill commented 3 years ago

It seems like this is a deficiency in the way that Antechamber reads in the .mol2 files. When writing a mol2 file, OpenBabel writes the ATOM section, then UNITY_ATOM_ATTR, then BOND. Antechamber does not expect the UNITY_ATOM_ATTR section, so it repeats the last atom for each line in the UNITY_ATOM_ATTR section, which messes up the input.

I would suggest passing the -xu flag for mol2 not to write the formal charge information. For instance, if a PDB is passed, OpenBabel will still write the UNITY_ATOM_ATTR section. This can also be used to parse existing mol2 files that might have the UNITY_ATOM_ATTR section. Removing charge information from the mol2 file will not be much of a loss, since it is already passed in via the -n flag, so ACPYPE already has that information.

Output of obabel -:"[Li+]" -o mol2

@<TRIPOS>MOLECULE
*****
 1 0 0 0 0
SMALL
GASTEIGER

@<TRIPOS>ATOM
      1 LI          0.0000    0.0000    0.0000 Li      1  UNL1        0.0000
@<TRIPOS>UNITY_ATOM_ATTR
1 1
charge 1
@<TRIPOS>BOND

Output of obabel -:"[Li+]" -o mol2 -xu

@<TRIPOS>MOLECULE
*****
 1 0 0 0 0
SMALL
GASTEIGER

@<TRIPOS>ATOM
      1 LI          0.0000    0.0000    0.0000 Li      1  UNL1        0.0000
@<TRIPOS>BOND
alanwilter commented 3 years ago

@mquevill Thanks for providing more details. I'm quite busy these days so sorry for not replying earlier. I still didn't have the time to look into this in full details but one thing I remember is that Antechamber may not work with Li, see FAQ.

Nevertheless, I tried your approach with O (oxygen):

>$ obabel -:"[O+]" -o mol2
@<TRIPOS>MOLECULE
*****
 1 0 0 0 0
SMALL
GASTEIGER

@<TRIPOS>ATOM
      1 O           0.0000    0.0000    0.0000 O.3     1  HOH1        1.0000
@<TRIPOS>UNITY_ATOM_ATTR
1 1
charge 1
@<TRIPOS>BOND
1 molecule converted

and, indeed, @<TRIPOS>UNITY_ATOM_ATTR trips a problem with Antechamber. I've read about Sybyl Mol2 format so, probably, Antechamber took a simple approach.

In the end, the onus is on the user when building the Mol2 file. I will post a question at Amber mailing list and add a note in my FAQ to warn about this situation.

And thanks for using ACPYPE!