From: Peter Freddolino (
Date: Thu May 24 2018 - 21:57:32 CDT

Dear David,
Molefacture is necessarily sort of a beast because of the amount of
overhead needed to allow vmd molecules to undergo arbitrary edits. Also, it
is heavily dependent on hooks associated with the gui, such that when I've
had to hack something together to use its internals from the command line,
I have often regretted it.

I think molefacture is probably overkill for your purpose. One key
question: what is the downstream application that you're going to use this
for? If it is for most any md simulation, you should just make use of
whatever topology generation utilities match your md package, and do the
hydrogen addition that way. If for some other purpose, I suspect that you'd
be better off just using text manipulation (e.g., cat) to add the hydrogen
atom to each file, and then a tcl script in vmd to position it (using the
vector commands to find a suitable position, and $sel move to move it). You
could get some good ideas on the needed geometry calculation
from ::Molefacture::calc_tetrahedral_geo


On Thu, May 24, 2018 at 12:58 PM, David Stevens <> wrote:

> Dear all,
> In order to perform a specific type of calculation, I need to take a
> large set of PDB files and protonate a single oxygen atom. I have
> decided to use molefacture, so I don't have to worry about calculating
> torsions, bends, and bond lengths. I can do all of the following manually,
> but if I could automate this through the command line, it would save me
> numerous afternoons worth of work.
> The main question I have is, how can I ensure that the oxygen of
> question is selected/picked?
> I have figured out how to call Molefacture through the command line:
> The beginning of the script is a for-loop that opens the PDB file and then
> performs the following:
> package require molefacture
> ::Molefacture::molefacture_gui_aux "index 3871" - this is the
> particular atom that I want to protonate - but this opens up the
> Molefacture gui (is there a way to not open the Molefacture gui?)
> ::Molefacture::add_hydrogen 0 (MolID number - but this command does not
> work, because the oxygen that I want to protonate isn't being automatically
> selected/highlight from the above command. Is there another command I
> should use to pick the atom through the command line?
> ::Molefacture::apply_changes_to_parent_mol
> set sel [atomselect top all]
> $sel writepdb O3H_prot.pdb.$i
> End of For-loop
> Thank you for your help,
> David