Creating new molecule

24 views
Skip to first unread message

Matthias Ernst

unread,
Apr 11, 2013, 1:22:46 PM4/11/13
to mdnalysis-...@googlegroups.com
Hi,

I started using MDAnalysis today and right now, I like it a lot.
Just one problem: I want to create a new series of atomic coordinates
and write them to a trajectory. Unfortunately, everything I found in the
documentation based on manipulating existing structures but I want to
create something new.
How can I do this?

I tried the following:
#format for it: [atomname atomindex resname resindex x y z]
a=MDAnalysis.core.AtomGroup.Atom(it[1], it[0], it[0], it[2], it[3], "A",
0.0, 0.0)
a.position(numpy.array(it[4:7]))

but this results in the following error:
Traceback (most recent call last):
File "test.py", line 244, in <module>
main()
File "test.py", line 197, in main
a.position(numpy.array(it[4:7]))
File
"/usr/local/lib/python2.6/dist-packages/MDAnalysis-0.7.7-py2.6-linux-x86_64.egg/MDAnalysis/core/AtomGroup.py",
line 219, in position
return self.universe.coord[self.number] # internal numbering starts at 0
File
"/usr/local/lib/python2.6/dist-packages/MDAnalysis-0.7.7-py2.6-linux-x86_64.egg/MDAnalysis/core/AtomGroup.py",
line 252, in universe
if not self.__universe == None: return self.__universe
AttributeError: _Atom__universe

How can I create a new set of atoms (not starting with a new universe
from an existing file: should be creating a whole new one!) and get them
into a (Gromacs xtc) trajectory file?

Thanks a lot,
Matthias

Oliver Beckstein

unread,
Apr 11, 2013, 10:14:50 PM4/11/13
to mdnalysis-...@googlegroups.com
Hi Matthias,

On 11 Apr, 2013, at 10:22, Matthias Ernst wrote:

> Hi,
>
> I started using MDAnalysis today and right now, I like it a lot.

Glad to hear it.

> Just one problem: I want to create a new series of atomic coordinates
> and write them to a trajectory. Unfortunately, everything I found in the
> documentation based on manipulating existing structures but I want to
> create something new.

Well, I hate to admit it but you stumbled into a major gap of MDAnalysis, at least in the currently released version.

Release 0.8 is supposed to have better model building capabilities.

Andy Somogyi has added a lot of code related to topology manipulation that is currently on the feature-psfgen branch (see http://code.google.com/p/mdanalysis/wiki/DevelopmentBranch on how to use the development branch; you can switch to a feature branch with 'git checkout -n feature-psfgen origin/feature-psfgen'). However, Andy knows more about how to use the new stuff than I do.

> How can I do this?

There may be a way to manually build the Universe._psf dictionary by hand; have a look at the topology parsers for how they do it. If you can get an "empty" Universe you could add Atom instances to _psf['_atoms'] and then run Universe._build_segments() again.

Maybe other people have better ideas...

>
> I tried the following:
> #format for it: [atomname atomindex resname resindex x y z]
> a=MDAnalysis.core.AtomGroup.Atom(it[1], it[0], it[0], it[2], it[3], "A",
> 0.0, 0.0)
> a.position(numpy.array(it[4:7]))

The above fails because coordinates are not stored with an atom (i.e. with the topology) but in Timestep objects that are part of a trajectory reader object. The Universe binds topology and trajectory together and without it you can't do much.

Oliver

>
> but this results in the following error:
> Traceback (most recent call last):
> File "test.py", line 244, in <module>
> main()
> File "test.py", line 197, in main
> a.position(numpy.array(it[4:7]))
> File
> "/usr/local/lib/python2.6/dist-packages/MDAnalysis-0.7.7-py2.6-linux-x86_64.egg/MDAnalysis/core/AtomGroup.py",
> line 219, in position
> return self.universe.coord[self.number] # internal numbering starts at 0
> File
> "/usr/local/lib/python2.6/dist-packages/MDAnalysis-0.7.7-py2.6-linux-x86_64.egg/MDAnalysis/core/AtomGroup.py",
> line 252, in universe
> if not self.__universe == None: return self.__universe
> AttributeError: _Atom__universe
>
> How can I create a new set of atoms (not starting with a new universe
> from an existing file: should be creating a whole new one!) and get them
> into a (Gromacs xtc) trajectory file?
>
> Thanks a lot,
> Matthias
>

--
Oliver Beckstein * orbe...@gmx.net
skype: orbeckst * orbe...@gmail.com

Matthias Ernst

unread,
Apr 12, 2013, 8:46:52 AM4/12/13
to mdnalysis-...@googlegroups.com
Hi Oliver,

thanks for your mail.
I had a gro writer on my own which I used to write a dummy structure
file. With this I could construct a universe and then manipulate the
coordinates. For me, this works reasonably well but is not the cleanest
solution).

I have another question concerning this issue:
when I write an xtc trajectory file, I can set for example the
universe.trajectory.ts.frame value and time. But when I write the xtc
file and have a look at it with gmxdump, every step is numbered "1"
although the time is set correctly. No matter if I set time, frame,
delta or whatever and independent, this always is "1". Shouldn't the
step increase with each written frame? This occurs no matter if I use
write or write_next_timestep function for trajectory.

Best regards,
Matthias

Oliver Beckstein

unread,
Apr 12, 2013, 12:01:04 PM4/12/13
to mdnalysis-...@googlegroups.com

On 12 Apr, 2013, at 05:46, Matthias Ernst wrote:

> I had a gro writer on my own which I used to write a dummy structure file. With this I could construct a universe and then manipulate the coordinates. For me, this works reasonably well but is not the cleanest solution).

I completely agree with you.

> I have another question concerning this issue:
> when I write an xtc trajectory file, I can set for example the universe.trajectory.ts.frame value and time. But when I write the xtc file and have a look at it with gmxdump, every step is numbered "1" although the time is set correctly. No matter if I set time, frame, delta or whatever and independent, this always is "1". Shouldn't the step increase with each written frame?

I suppose yes, so it might be a bug.

> This occurs no matter if I use write or write_next_timestep function for trajectory.

Can you file an issue together with code that reproduces the problem? (See http://code.google.com/p/mdanalysis/wiki/ReportingProblems for details.)

Thanks,
Oliver
Reply all
Reply to author
Forward
0 new messages