michellab / Sire

Sire Molecular Simulations Framework
http://siremol.org
GNU General Public License v3.0
95 stars 26 forks source link

Format of solvent pointers #371

Closed kexul closed 2 years ago

kexul commented 2 years ago

Dear Sire developers: Recently I found that pymol cannot load the topology file generated by Sire when using the command: load somd.prm7, somd, 0, top, which raised an error:

ObjectMolecule: Attempting to read Amber7 topology file.
 ObjectMolecule-Error: Unrecognized file format (can't find "%FORMAT(3I8)").
ObjectMoleculeTOPStr2CoordSet-Error: failed

It was caused by these lines in the topology file

%FLAG SOLVENT_POINTERS
%FORMAT(10I8)
    1806    1817    1807

which should be written by this line: https://github.com/michellab/Sire/blob/3a664356096964546df3108c889c6b7117bd7709/corelib/src/libs/SireIO/amberprm.cpp#L3338

I found that Amber suggests a format of 3I8 in their manual for SOLVENT_POINTERS, should we use that format instead?

lohedges commented 2 years ago

Thanks for this. As you say, the online AMBER format guide (here) implies that 3I8 should be used. I'm surprised that pymol doesn't just try to read the data irrespective of the format, since it is just an integer regardless.

I'll update this now.

lohedges commented 2 years ago

Actually, I believe that the AMBER format page might be incorrect and pymol should actually read 10I8 as we use. Indeed, Parmed write this format too. See here (search for SOLVENT_POINTERS.)

lohedges commented 2 years ago

Maybe not, since the actual Parmed code says otherwise.

kexul commented 2 years ago

I'm surprised that pymol doesn't just try to read the data irrespective of the format, since it is just an integer regardless.

When I modify the extension of the file from .prm7 to .parm7, the topology file could be loaded by pymol directly by load somd.parm7 , interesting :)

lohedges commented 2 years ago

It looks like AMBER might have change the format between the parm/top and prm7 format. Looking more closely at the Parmed code I see that they also have another parser where they use 10I8.

kexul commented 2 years ago

It looks like AMBER might have change the format between the parm/top and prm7 format.

Thanks! Just tried load somd.prm7, somd, 0, parm7 and it works! 🤗

lohedges commented 2 years ago

Using ParmEd and saving to both the .parmtop and .parm7 extension gives a SOLVENT_POINTERS record using the 3I8 format. I'm not sure the reason for the inconsistency in their code, or what would trigger the writing of 10I8 records instead. (The parser using 10I8 has been updated much more recently.

The full AMBER specification manual also states 3I8, so I'm going to go with that as the canonical reference.