djay0529 / mdanalysis

Automatically exported from code.google.com/p/mdanalysis
0 stars 0 forks source link

Add wrap method #190

Open GoogleCodeExporter opened 9 years ago

GoogleCodeExporter commented 9 years ago
Hello,

I'm currently a VMD user looking for better solutions to analyze my simulations 
using Python. I found that MDAnalysis is the best way to go! But it lacks some 
really important features that I used when scripting within VMD. The most 
important one (that delayed my migration to this package) is the ability to 
wrap systems into the original simulation cell, just like VMD's PBC Tools. 
Since I really like the way MDAnalysis works and want to use it, I wrote an 
wrap method that behave like `pbc wrap`. It is able to wrap orthorhombic and 
non-orthorhombic cells by moving atoms, residues, segments etc. One interesting 
group to be moved together is a fragment (a group of atoms that are bonded to 
the same molecule). So I also created an atom group named `fragment` and added 
it's builder to the `MDAnalysis.Universe`. Finally to avoid this `fragment` 
class being useful only when wrapping systems, I also added the `same` keyword 
into selections. This way one could, e.g., select whole lipid molecules using 
something like `segid DMPC and same fragment as around 3.5 protein`, or the 
complete protein residues that are in contact with water using `protein and 
same residue as around 3 resname TIP3`.

All these modifications can be found in 
https://code.google.com/p/mdanalysis-caio/source/list
The commits correspond to:
  * 1a82890f0c19: `wrap` method added to MDAnalysis.AtomGroup;
  * c03cbd2dfe72: `fragment` class an it's builder;
  * f160a65892a9: `same as` keyword added to selections.

I hope these codes help to improve MDAnalysis.
Please, let me know if you have any questions or if you want me to add/modify 
anything within this code.

Cheers

Original issue reported on code.google.com by caiobio...@gmail.com on 12 Jul 2014 at 3:47

GoogleCodeExporter commented 9 years ago
Hi caiobiounb,

Thank you for the contribution. I had a quick look at your proposed additions. 
They look all very interesting (wrap-functionality and fragments with addition 
to the selection parser). I have a number of questions/requests, though:

1) Can you briefly describe the new Fragments group, its use and implementation?
   - What is the use case? Can you give an example how it is used?
   - How is a fragment constructed (e.g. using bonded information by a topology reader) and when/where in the code?
   - Are fragments also constructed for e.g. water molecules and by how much does this slow down the initial universe construction?
   - What happens if no bonded information is available?

2) Can you briefly explain how AtomGroup.wrap() improves over 
AtomGroup.packintobox() (see Issue 153). Can you make use of the work that has 
already been done  for Issue 136 (and see also Issue 156)?

On the technical side: In order to integrate code we'll eventually also need 
the following:

1) Can you please rebase against the current develop branch instead of master? 
Otherwise merging becomes messy.

2) These are changes in the core so we'll need UnitTests. 

3) We'll also need documentation. You already have docstrings, which is great, 
but we'll also need something along the lines of your answers to the first 
question as part of the more general docs.

Cheers,
Oliver

Original comment by orbeckst on 14 Jul 2014 at 7:25

GoogleCodeExporter commented 9 years ago
wrap seems similar to packintobox, but then adds a lot of very nice options.  I 
think wrap() with no options is identical to packintobox(), so if we add wrap 
we can scrap packintobox.

I think it's important to add the inplace flag into this method, it handy for 
lots of methods to "non destructively" move atoms around.

Re: reusing code, the first half of the method to generate the cell seems to be 
the same as topology.core.triclinic_vectors
https://code.google.com/p/mdanalysis/source/browse/package/MDAnalysis/coordinate
s/core.py#280

In core.distances.applyPBC (src/numtools/distances.pyx and calc_distances.h in 
the source) there's some Cython that should be able to be used.
https://code.google.com/p/mdanalysis/source/browse/package/src/numtools/distance
s.pyx#615

applyPBC(coords, box) returns a copy of the coords shifted to all be within the 
defined box.  Box can be either orthogonal or triclinic.  I think it should be 
possible to do all the transformations as a translation then putting this into 
applyPBC.

Original comment by richardjgowers on 16 Jul 2014 at 10:40

GoogleCodeExporter commented 9 years ago
Hello
Thanks to your comments!

Fragments are well suited to select a group of atoms that are bonded to the
same molecule without knowing how the structure file is organized (it's not
important whether a molecule is grouped using segid, resid etc). This way,
if one selects just one atom, it's possible to get the entire molecule to
which this atom is bonded just using it's fragment attribute. By using a
selection like "resname TIP3 and same fragment as around 4 protein" you can
select not only the single atoms that are within 4 angstrom far from the
protein (given by the "around" keyword), but the complete water molecules.
Lastly, fragments are of great importance when wrapping, as stated below.
Fragments are constructed when Universe._init_topology is called by using
bond information, so any atoms bonded to the same molecule belong to unique
fragments. Also, isolated atoms (like ions) have their own fragments
containing just one atom. If no bond information is given then fragments
are not created.
The overhead of this process is small, since loading a PSF with 305,000
atoms took 12s, from which 1.5s corresponded to the fragments building.

When using AtomGroup.packintobox(), atoms that belong to the same molecule
(or segid, or residue etc) are moved indepedently, resulting in bonds
crossing the unitcell. With AtomGroup.wrap() it's possible to move group of
atoms together (residues, segids etc). When using it with fragments, entire
molecules are translated, thus there are no bonds crossing the simulation
cell. This is very useful when visualizing/constructing MD systems. It's
also possible to centralize the cell according to the center of any
selection of atoms, shift this center by user-defined values. Users may
also provide their own cell dimensions.

I'm very busy these days, but as soon as possible I'll provide UnitTests
and see how feasible is to use the code within applyPBC.

Here are the new list of commits, now rebased to the develop branch (
https://code.google.com/p/mdanalysis-caio/source/list)
  * *743cb0f5d2ae*: `wrap` method added to MDAnalysis.AtomGroup;
  *  *5c1f6497b3bb*: `fragment` class an it's builder;
  *  *3f1f6c938186*: `same as` keyword added to selections.

Original comment by caiobio...@gmail.com on 24 Jul 2014 at 1:05

GoogleCodeExporter commented 9 years ago
I'll take this on and finish it.

Original comment by richardjgowers on 21 Sep 2014 at 10:23

GoogleCodeExporter commented 9 years ago
Fragments added, but lazily built.  Wrap method still todo.

Original comment by richardjgowers on 16 Dec 2014 at 12:35

GoogleCodeExporter commented 9 years ago

Original comment by richardjgowers on 11 Jan 2015 at 1:19

GoogleCodeExporter commented 9 years ago

Original comment by orbeckst on 16 Feb 2015 at 7:28

GoogleCodeExporter commented 9 years ago
I don't really know how the selection methods work in MDA so I can't add the 
extra keyword.

The wrap method looks like it will work, but should have a lot of tests behind 
it to test all the combinations of keywords and geometry, so this will have to 
be post 0.9.

Original comment by richardjgowers on 16 Feb 2015 at 11:49

GoogleCodeExporter commented 9 years ago

Original comment by richardjgowers on 16 Feb 2015 at 1:27