VMD-L Mailing List

From: Paweł Kędzierski (pawel.kedzierski_at_pwr.edu.pl)
Date: Wed Mar 05 2025 - 02:06:02 CST

Dear All,

There is a related bug in Molefacture 2 affecting FFTK users. When you
do not make any changes but just want to save the topology file, it is
not saved without any message.

The Josh suggestion cited below may therefore not work in general.
Actually what I observe that in these workflows when saving PSF+PDB
files for a new ligand does actually work, Molefacture ALSO saves the
top file with the same prefix. So perhaps the primary cause of both bugs
is that sometimes Molefacture "forgets" to create the the TOP file which
is needed by psfgen to create the PDB+PSF pair.

Now the problem I observe:
When using FFTK, one ends up with PDB and PSF files. I use Molefacture
to covert the PSF to a TOP file format to be used with the regular
system preparation procedure.
The steps to reproduce:
- load a molecule from PSF and PDB files,
- open Molefacture
- Molefacture menu File-> New molecule -> From selection, type "all"
- Molefacture menu File-> Write top/str file
Neither str nor top format works, nothing is written. Saving xyz or mol2
formats does work.
When I tried to save psf+pdb pair of files, vmd actually crashed (tested
with the distributed binary version 1.9.4a57 for Linux) with these messages:

psfgen) building segment O1
psfgen) reading residues from pdb file O1.pdb
psfgen) unknown residue type MSA
psfgen) extracted 1 residues from pdb file
psfgen) setting patch for first residue to NONE
psfgen) setting patch for last residue to NONE
psfgen) enabling angle autogeneration
psfgen) enabling dihedral autogeneration
Info: generating structure...
psfgen) unknown residue type MSA
double free or corruption (out)

I have found a workaround, namely if you introduce any change in the
Molefacture editor, even something like oxidation state of any atom
(irrelevant for the TOP file), then you can save the TOP file and the
PSF+PDB files normally. But it is clearly a bug.

As for the topic cited, I think I remember that the issue with all-zero
PDB coords only happend to me when the molecule was read from a PDB
file. When I start with the ligand read from a MOL2 file, the problem
does not appear, but again I have to do some editing first to be able to
save the PSF+PDB files for FFTK.

Best regards,
Paweł

W dniu 4.03.2025 o 19:11, Josh Vermaas pisze:
> Hi Mike,
>
> Yeah, there have been changes to molefacture since 2015, and I think
> that none of the currently listed maintainers are still actively
> working on VMD, so it wouldn't surprise me that there is a bug
> *somewhere*. I can reproduce that a psf/pdb written from molefacture
> based on the tutorial gives you zeros for all the coordinates, but the
> xyz files work!
>
> Now, I can sort of see where this goes wrong. Its in the
> molefacture/molefacture_internals.tcl file, in the write_psfpdbfiles
> procedure, lines 2756-2800 is what is actually being used to build the
> psf/pdb. Note that by default, VMD's new and improved molefacture
> tries to be clever, and will load topology files that already are
> distributed with VMD, including top_all36_cgenff.rtf (which has a ETOH
> residue defined in it). To do what you want it to, you'll need to
> create another topology file for the molecule you made in molefacture
> first (File->write top/str file), which is appended to the list of
> topology files in your psfcontext. Only after you write this top/str
> file, will psfgen recognize the atoms and put them in the right place.
>
> Hope this helps!
>
> -Josh
>
> On 3/4/25 12:23 PM, CHARLES MCCALLUM wrote:
>> I really appreciate Josh’s willingness to help (Go State!), but I
>> think probably I need help from an author of the ffTK tutorial from
>> 2015.  I believe the point of this tutorial is to
>> 1) show how to parameterize a new, smallish molecule
>> 2) use ethanol, from scratch, because the results may be compared to
>> the established ETOH, standard parameterization.
>>
>> After Josh’s reply, I re-did the tutorial, and tried a couple things.
>>
>> First, I made the resname “ETO”, but that returned no coordinates.
>>  This made me realize Molefacture IS taking coordinates from the
>> established ETOH.  I then re-ran and made sure the oxygen was
>> attached to the C1 carbon (though the OH hydrogen ended up being
>> index 8/atom 9 in the Molefacture table), and used “ETOH”.  This
>> returned a full set of coordinates, so I guess we can move on, but
>> I’m a bit skeptical about the rest of the workflow.
>>
>> REMARK original generated coordinate pdb file
>> ATOM      1  C1  ETOHA   1 0.000   0.000   0.000  1.00  0.00      A    C
>> ATOM      2  O1  ETOHA   1  -0.495   1.437  -0.000  1.00  0.00      A
>>    O
>> ATOM      3  HO1 ETOHA   1  -1.902   1.714  -0.000  1.00  0.00      A
>>    H
>> ATOM      4  H11 ETOHA   1  -0.326  -0.473   0.819  1.00  0.00      A
>>    H
>> ATOM      5  H12 ETOHA   1  -0.326  -0.473  -0.819  1.00  0.00      A
>>    H
>> ATOM      6  C2  ETOHA   1 1.500   0.000   0.000  1.00  0.00      A    C
>> ATOM      7  H21 ETOHA   1 1.872  -0.103   1.042  1.00  0.00      A    H
>> ATOM      8  H22 ETOHA   1 1.899  -0.852  -0.592  1.00  0.00      A    H
>> ATOM      9  H23 ETOHA   1 1.909   0.945  -0.419  1.00  0.00      A    H
>> END
>>
>>
>> Why does the tutorial explicitly use a different order?  Somehow this
>> worked, as the example files have a full set with this different order:
>>
>> REMARK original generated coordinate pdb file
>> ATOM      1  C1  ETOHL   1 0.000   0.000   0.000  1.00  0.00      L    C
>> ATOM      2  H11 ETOHL   1  -0.326   0.946  -0.000  1.00  0.00      L
>>    H
>> ATOM      3  H12 ETOHL   1  -0.326  -0.473   0.819  1.00  0.00      L
>>    H
>> ATOM      4  H13 ETOHL   1  -0.326  -0.473  -0.819  1.00  0.00      L
>>    H
>> ATOM      5  C2  ETOHL   1 1.500   0.000   0.000  1.00  0.00      L    C
>> ATOM      6  H21 ETOHL   1 1.872  -0.103   1.042  1.00  0.00      L    H
>> ATOM      7  H22 ETOHL   1 1.899  -0.852  -0.592  1.00  0.00      L    H
>> ATOM      8  O2  ETOHL   1 2.059   1.292  -0.573  1.00  0.00      L    O
>> ATOM      9  HO2 ETOHL   1 1.186   2.397  -0.844  1.00  0.00      L    H
>> END
>>
>> What am I missing here?
>>
>> Best,
>> Mike
>>
>>> On Mar 4, 2025, at 07:42, CHARLES MCCALLUM <cmccallum_at_me.com> wrote:
>>>
>>> Josh, thanks.
>>>
>>> However, I understand how psfgen does things, as we’ve been doing
>>> that for years.  However, I’m trying to do things through
>>> Molefacture.  Why isn’t following the tutorial working?  The goal is
>>> toe parameterize some new molecules we want to use.  If I’m trying
>>> to teach people an easier way, front to back, to parameterize new
>>> molecules, I’d like to be able to stick to one tool, and have the
>>> tutorial explain things.
>>>
>>> Best
>>> Mike
>>>
>>>
>>>> On Mar 4, 2025, at 07:18, Josh Vermaas <vermaasj_at_msu.edu> wrote:
>>>>
>>>> Hi Charles,
>>>>
>>>> If you are building a psf/pdb from psfgen, the order of the atoms
>>>> will be based on the order for the residue in the topology file
>>>> (usually .top or .rtf). Here is the relevant section from
>>>> top_all36_cgenff.rtf:
>>>>
>>>> RESI ETOH          0.00 ! C2H6O, Ethanol, adm jr.
>>>> GROUP
>>>> ATOM C1   CG321    0.05 !  H21  H11 H12
>>>> ATOM O1   OG311   -0.65 !     \   \  /
>>>> ATOM HO1  HGP1     0.42 ! H22--C2--C1
>>>> ATOM H11  HGA2     0.09 !     /      \
>>>> ATOM H12  HGA2     0.09 !  H23 O1--HO1
>>>> GROUP
>>>> ATOM C2   CG331   -0.27
>>>> ATOM H21  HGA3     0.09
>>>> ATOM H22  HGA3     0.09
>>>> ATOM H23  HGA3     0.09
>>>>
>>>> You'll notice that your order matches here! The other thing that I
>>>> notice is that you have a few atoms where the occupancy column is
>>>> -1. That means that the input pdb for psfgen did not have an atom
>>>> that matched the name that psfgen expected, and so nothing was
>>>> loaded in. If you are doing this in a standard psfgen script, I
>>>> think these two lines before you write out your psf/pdb will do the
>>>> trick:
>>>>
>>>> regenerate angles dihedrals
>>>> guesscoord
>>>>
>>>>
>>>> -Josh
>>>>
>>>> On 3/3/25 9:34 PM, CHARLES MCCALLUM wrote:
>>>>> Hello, we recently starting learning to us Molefacture to
>>>>> parameterize, and we are running into some difficulties.
>>>>> Using the tutorial
>>>>> (https://www.ks.uiuc.edu/Training/Tutorials/science/ffTK/fftk-tutorial.pdf
>>>>> <https://urldefense.com/v3/__https://www.ks.uiuc.edu/Training/Tutorials/science/ffTK/fftk-tutorial.pdf__;!!HXCxUKc!05ug8JAR-3nPHBSss6MjnihamaigCYNB2fkqMBS_p9TkaoRVanfXqPDOuSym8SGF7WCqcAUQt445pzE$>)
>>>>> for ethanol, we are finding some coordinates are not being written
>>>>> out in the PDB file.  We have edited the table as the tutorial notes.
>>>>>
>>>>> The console reports PSFGEN output, and all seems ok until it
>>>>> reports “failed to set coordinates” for H13, O2 and HO2 (for some
>>>>> reason I am not able to copy and paste the tkConsole output).  The
>>>>> funny thing is, the output PDB contains more missing coordinates,
>>>>> and it has some labels wrong or missing, and it seems to have
>>>>> changed the order (shouldn’t these be in the order of the index in
>>>>> the table?):
>>>>>
>>>>> REMARK original generated coordinate pdb file
>>>>> ATOM    1  C1  ETOHL   1       0.000   0.000   0.000  1.00  0.00
>>>>>  L    C
>>>>> ATOM    2  O1  ETOHL   1       0.000   0.000   0.000 -1.00  0.00  L
>>>>> ATOM    3  HO1 ETOHL   1       0.000   0.000   0.000 -1.00  0.00  L
>>>>> ATOM    4  H11 ETOHL   1      -0.326   0.946  -0.000  1.00  0.00
>>>>>  L    H
>>>>> ATOM    5  H12 ETOHL   1      -0.326  -0.473   0.819  1.00  0.00
>>>>>  L    H
>>>>> ATOM    6  C2  ETOHL   1       1.500   0.000   0.000  1.00  0.00
>>>>>  L    C
>>>>> ATOM    7  H21 ETOHL   1       1.872  -0.103   1.042  1.00  0.00
>>>>>  L    H
>>>>> ATOM    8  H22 ETOHL   1       1.899  -0.852  -0.592  1.00  0.00
>>>>>  L    H
>>>>> ATOM    9  H23 ETOHL   1       0.000   0.000   0.000 -1.00  0.00  L
>>>>> END
>>>>>
>>>>> Writing out an xyz file does not cause problems:
>>>>> 9
>>>>>  generated by VMD
>>>>>  C       0.000000        0.000000        0.000000
>>>>>  H      -0.325568        0.945519       -0.000000
>>>>>  H      -0.325568       -0.472759        0.818843
>>>>>  H      -0.325568       -0.472759       -0.818843
>>>>>  C       1.500000        0.000000        0.000000
>>>>>  H       1.872000       -0.103000        1.042000
>>>>>  H       1.899000       -0.852000       -0.592000
>>>>>  O       2.059218        1.292080       -0.572891
>>>>>  H       1.186122        2.396645       -0.844260
>>>>>
>>>>> We are also seeing occasional problems with the table deleting our
>>>>> edits, leaving blanks.  When this happens, we have found we have
>>>>> to restart VMD in order to use Molefacture again.
>>>>>
>>>>> VMD 1.9.a57-x86_64-Rev12 on OS X 15.3.1
>>>>>
>>>>> Any help would be appreciated!
>>>>>
>>>>> Mike
>>>>>
>>>>> --
>>>>> Mike McCallum
>>>>> Morada/Stockton CA, USA
>>>>> cmccallum at me ^dot^ com
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>
>>>> --
>>>> Josh Vermaas
>>>>
>>>> vermaasj_at_msu.edu
>>>> Assistant Professor, Plant Research Laboratory and Biochemistry and Molecular Biology
>>>> Michigan State University
>>>> vermaaslab.github.io
>>>
>>> --
>>> Mike McCallum
>>> Morada/Stockton CA, USA
>>> cmccallum at me ^dot^ com
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>
>> --
>> Mike McCallum
>> Morada/Stockton CA, USA
>> cmccallum at me ^dot^ com
>>
>>
>>
>>
>>
>>
>>
>
> --
> Josh Vermaas
>
> vermaasj_at_msu.edu
> Assistant Professor, Plant Research Laboratory and Biochemistry and Molecular Biology
> Michigan State University
> vermaaslab.github.io 

-- 
Room 3.03 A-3, tel. +48 71 320 3200
email:Pawel.Kedzierski_at_pwr.edu.pl
Institute of Advanced Materials I01W03D10
Faculty of Chemistry
Wrocław University of Science and Technology
Wrocław, Poland