Re: Report on a trajectory incompatibility between NAMD and CHARMM: Unfolding and RDFs

From: NAVEEN P MICHAUD AGRAWAL (nmichaud_at_jhu.edu)
Date: Sat Jul 30 2005 - 10:42:29 CDT

Looking over what I changed in dcdplugin.c, I think the reason it works is that catdcd uses VMD's file input plugins and can read both NAMD's dcd cell format and CHARMM's dcd cell format, but the dcd writing code uses CHARMM's format. So basically it just rewrites the cell dimensions in the correct format.

It seems like your approach only works when you have fixed dimensions (or you saved a record of the dimensions at every timestep you save the coordinates in another file besides the dcd). Otherwise you lose information on the cell dimension fluctuations.

Naveen

----- Original Message -----

From: Mauricio Carrillo Tripp <trippm@gmail.com>

Date: Saturday, July 30, 2005 11:07 am

Subject: Re: namd-l: Report on a trajectory incompatibility between NAMD and CHARMM: Unfolding and RDFs

> Hi Naveen, thanks for the reply and the tip. I guess I found a
> different way
> to
> correct the problem, without knowing exactly what it was. Although
> it would
> be
> nice to be able to use the dcd file directly from namd, without
> any steps in
> between.
> Or at least the manual should have a warning saying it won't work
> and a
> re-formatting
> is needed, instead of claiming namd uses charmm's dcd format, you
> know?
>
> On 7/29/05, NAVEEN P MICHAUD AGRAWAL <nmichaud@jhu.edu> wrote:
> >
> >
> > I had this same exact problem, and it was apparent when I viewed
> the
> > merged trajectories in VMD because I
> > used CHARMM to merge them. This issue was raised in the NAMD
> wiki:
> > http://www.ks.uiuc.edu/Research/
> > namd/wiki/index.cgi?NamdAndCHARMM. I hacked the dcdplugin.c file
> from the
> > newest version of catdcd to
> > take a trajectory and write out a new one with the correct cell
> > information. Tthis was back in January so I don't
> > remember exactly what i did, but I've included a tar of the
> source
> > directory (which includes the files i copied
> > from catdcd...hopefully I
> > am allowed to redistribute it). You can see my changes by doing
> a diff
> > between dcdcellfix.c and dcdplugin.c
> > (from the catdcd
> > program). To compile it just type
> >
> > gcc -I./ -DTEST_DCDPLUGIN dcdcellfix.c -o dcdcellfix -lm
> >
> > and then
> >
> > ./dcdcellfix [trajectory filename] [output trajectory]
> >
> > If you don't specify an output file it will default to test.trj.
> Hopefully
> > this will be helpful to anyone else having this
> > problem.
> >
> > Naveen
> >
> > ----- Original Message -----
> > From: Mauricio Carrillo Tripp <trippm@gmail.com>
> > Date: Friday, July 29, 2005 5:36 pm
> > Subject: namd-l: Report on a trajectory incompatibility between
> NAMD and
> > CHARMM: Unfolding and RDFs
> >
> > > Report on a trajectory incompatibility between NAMD and CHARMM.
> > >
> > > I describe a couple of problems we've had when applying Charmm's
> > > analysis
> > > tool set to a
> > > trajectory obtained from a NAMD simulation (using Charmm's FF) and
> > > the way
> > > we came
> > > around them. I apologize for the length of this posting but I hope
> > > that it
> > > will
> > > help other people who, either now or in the future, try to follow
> > > the same
> > > analysis procedure.
> > > What I describe might even be a bug in NAMD, and hopefully the
> > > developers
> > > will fix it in
> > > future versions. We currently use NAMD 2.5 for Linux-i686-TCP.
> > > We usually build the system using Charmm, creating the psf file
> > > both in
> > > Charmm and
> > > XPlore formats. We use the latter to run the simulation on NAMD,
> > > obtaining a
> > > DCD trajectory
> > > file. Finally, we go back to Charmm and analyze the DCD.
> > > The systems we work with are lipid membranes on different
> > > conditions and
> > > different sizes;
> > > this might be relevant for the problems described below since they
> > > are
> > > anisotropic.
> > >
> > > Unfolding:
> > >
> > > Some time ago I reported to this list a problem I had while trying
> > > to unfold
> > > a trajectory obtained from a
> > > NAMD simulation using Charmm's FF. The problem was that the
> > > coordinates got
> > > messed up
> > > during the process of unfolding the system using the Charmm's tool
> > > set. This
> > > was clearly seen when
> > > looking at the trajectory with VMD, the output would look
> something> > like
> > > this:
> > > http://chem.acad.wabash.edu/~trippm/tmp/unfold-
> > > namd.png<http://chem.acad.wabash.edu/%7Etrippm/tmp/unfold-
> namd.png>> >
> > > We tried a little trick (read below), and that fixed the problem:
> > > http://chem.acad.wabash.edu/~trippm/tmp/unfold-
> > > charmm.png<http://chem.acad.wabash.edu/%7Etrippm/tmp/unfold-
> charmm.png>> >
> > > RDFs:
> > >
> > > Recently we ran into some confusing results while trying to
> > > calculate a
> > > radial distribution
> > > function, also using the Charmm's tool set, on a trajectory file
> > > produced by
> > > NAMD. We used
> > > a couple of different commands, different options, different
> system> > size and
> > > different image
> > > transformations. All of them should have given the same result.
> > > However,
> > > they didn't.
> > > Most confusing was the fact that any of those results were the
> same> > to the
> > > result obtained from
> > > the same analysis made to a trajectory file produced by Charmm
> > > (which seemed
> > > to be the right answer).
> > > We didn't understand why would this be happening, until we
> > > remembered the
> > > 'small' issue we had before
> > > when trying to unfold the system. We applied the same 'trick' and
> > > everything
> > > worked just perfect!
> > > In the next figure, panel (a), are shown the results of the
> > > analysis done on
> > > the trajectory file produced by
> > > NAMD. It can be clearly seen that there is a disagreement between
> > > all of the
> > > methods used.
> > > Panel (b) shows all the results of the same methods applied to the
> > > 'fixed'
> > > trajectory file:
> > > http://chem.acad.wabash.edu/~trippm/tmp/rdf-
> > > comparison.png<
> > http://chem.acad.wabash.edu/%7Etrippm/tmp/rdf-comparison.png>
> > > As shown, only when 'fixing' the trajectory file all the methods
> > > used to
> > > calculate the RDF give the
> > > same result. Moreover, they are the same as the result obtained
> > > when
> > > analysing the trajectory
> > > produced by Charmm.
> > >
> > > The 'fix':
> > >
> > > There seems to be a difference between the way NAMD writes the
> > > cell's info
> > > to the DCD file and
> > > the way Charmm reads it (and writes it). This is why we tried to
> > > reset this
> > > info in the trajectory file.
> > > The first step is to read in the trajectory produced by NAMD into
> > > Charmm.
> > > Then, we dump only the
> > > coordinates into a temporary file. After defining the cell's
> > > lengths and
> > > geometry, we read in the
> > > coordinates from the temporary file and write the final
> 'fixed' DCD
> > > trajectory file. After doing this,
> > > all the analisys made with Charmm that involve periodic boundaries
> > > operations and the creation of
> > > images should work. This might not be the best or optimized
> > > solution, but it
> > > has worked for us so far.
> > >
> > > As I mentioned before, this problem might be only apparent for
> > > systems that
> > > are not isotropic, but for consistency's
> > > sake the trajectory file should be corrected anyway. At least
> until> > this bug
> > > is fixed in the code, shall it be one.
> > >
> > >
> > > For those who would like to try it for themselves, here's the
> > > Charmm script:
> > >
> >
> >
> ************************************************************************************************************************> **************************************
> > >
> > > ! Charmm's configuration file to RESET the Cell's info in the
> > > trajectory
> > > file
> > > ! to overcome an incompatibility between NAMD and Charmm
> > > ! by trippm@wabash.edu and fellers@wabash.edu
> > > ! Wabash College 290705
> > >
> > > set OPT 0 ! set this to [0] if just resetting the cell's info, [1]
> > > if
> > > unfolding the system
> > >
> > > set NS 500000 ! set this to the total number of steps in the
> > > simulationset NF 1000 ! set this to the actual number of
> frames in
> > > the trajectory file
> > >
> > > set boxx 49.00 ! set this to the cell's X length
> > > set boxy 49.00 ! set this to the cell's Y length
> > > set boxz 54.00 ! set this to the cell's Z length
> > >
> > > ! Read in the PSF file
> > > open unit 50 read card name my_system.psf
> > > read psf card unit 50
> > > close unit 50
> > >
> > > ! Read in the NAMD DCD trajectory file, and write it out to a temp
> > > fileopen unit 11 read file name my_system-namd.dcd
> > > open unit 1 write file name tmp.dcd
> > >
> > > ! Calculate the skip value
> > > set SK @NS
> > > divide SK by @NF
> > >
> > > merge coor firstu 11 nunit 1 skip @SK output 1 nfile @NF
> > >
> > > close unit 1
> > > close unit 11
> > >
> > > ! Define the cell's info, Charmm's format...
> > > crystal define tetragonal @boxx @boxy @boxz 90.0 90.0 90.0
> > > open unit 1 read card name image.cry
> > > crystal read unit 1 card
> > > close unit 1
> > >
> > > ! Read in the trajectory temp file, and write it out to a Charmm's
> > > formated
> > > file
> > > open unit 11 read file name tmp.dcd
> > > open unit 1 write file name my_system-charmm.dcd
> > >
> > > if OPT .eq. 0 then merge coor firstu 11 nunit 1 skip @SK
> output 1
> > > nfile @NF
> > > if OPT .eq. 1 then merge coor firstu 11 nunit 1 skip @SK
> output 1
> > > nfile @NF
> > > unfold
> > >
> > > stop
> > >
> > >
> >
> >
> ************************************************************************************************************************> **************************************
> > >
> > > --
> > >
> > > Mauricio Carrillo Tripp, PhD
> > > Department of Chemistry
> > > Wabash College
> > > trippm@wabash.edu
> > > http://trippm.bajacast.com/
> > >
> >
> >
> >
> >
>
>
> --
>
> Mauricio Carrillo Tripp, PhD
> Department of Chemistry
> Wabash College
> trippm@wabash.edu
> http://trippm.bajacast.com/
>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:39:45 CST