From: Joshua Adelman (jla65_at_pitt.edu)
Date: Tue Mar 12 2013 - 14:42:25 CDT

Hi John,

If I read more than one frame in, I get the same results as if I had read N frames in, where N > 2. So I would get the same results for the energy of frame 0 with any of the commands that you wrote. I see this problem consistently for several different dcd files that I've tested on.

I haven't dug through the tcl code, but the documentation makes it seem like the unit cell should be read from the dcd:

"However, for proper initialization of a periodic system NAMD must also be fed an initial guess for the periodic cell, which will then be replaced by the data from the trajectory as NAMDenergy runs. "

Josh

On Mar 12, 2013, at 3:29 PM, John Stone wrote:

>
> Joshua,
> Do you also get different results if you change it to:
> animate read dcd $dcd_file beg 0 end 1 waitfor all
> or
> animate read dcd $dcd_file beg 0 end 2 waitfor all
> or
> animate read dcd $dcd_file beg 0 end 3 waitfor all
>
> If so, I would guess that this has to do with the definition of the
> periodic cell, and that namdenergy may be taking the periodic cell
> from the last frame (PBC info is in the DCD trajectory, not the PSF file).
> This is just a guess, as I didn't write any of that code.
>
> Cheers,
> John Stone
> vmd_at_ks.uiuc.edu
>
> On Tue, Mar 12, 2013 at 12:53:01PM -0400, Joshua Adelman wrote:
>> I was wondering if anyone had any insight into a strange result I'm getting when using the namdenergy plugin. I'm using VMD 1.9.1 and version 1.4 of the namdenergy package. I'm calculating the electrostatic energy for a series of snapshots, and realized that the results differ depending on whether all of the frames were loaded in at once or if a separate call to namdenergy was done for each frame.
>>
>> So for example, if my script looks something like (read in 10 frames from a dcd file):
>>
>> proc get_energy {sys_psf dcd_file out_name tmp_file xsc_file ff1 ff2 ff3} {
>> package require namdenergy
>> mol new $sys_psf
>> animate read dcd $dcd_file beg 0 end 10 waitfor all
>> set all [atomselect top "all"]
>> namdenergy -switch 8.0 -cutoff 10.0 -elec -sel $all -extsys $xsc_file -pme -par $ff1 -par $ff2 -par $ff3 -ofile $out_name -tempname $tmp_file -debug
>> }
>>
>> My results look like:
>>
>> Frame Time Elec Total
>> 0 0 -94731.9 -94731.9
>> 1 1 -94733.5 -94733.5
>> 2 2 -94574.4 -94574.4
>> 3 3 -94674.7 -94674.7
>> 4 4 -94597.5 -94597.5
>> 5 5 -94239.8 -94239.8
>> 6 6 -94677.2 -94677.2
>> 7 7 -94478.9 -94478.9
>> 8 8 -94750.7 -94750.7
>> 9 9 -94442.5 -94442.5
>> 10 10 -94824.9 -94824.9
>>
>> However, if I keep everything else the same (force field files, xsc, etc, but I change the animate read statement to:
>>
>> animate read dcd $dcd_file beg 0 end 0 waitfor all
>>
>> the energy returned for that first frame is:
>> Frame Time Elec Total
>> 0 0 -87170.6 -87170.6
>>
>> The same energy is obtained if I make a pdb from the first frame and load it as:
>> mol addfile fr0.pdb
>>
>> I've checked that the same PME grid size is being used in all cases, and that the {a b c alpha beta gamma} values are the same for the pdb frame and dcd frames. I have also noticed that if I use a different xsc file for the single frame calculations that I different energy, but it doesn't effect the results if I'm processing a dcd file. I thought that the xsc file was just there for setup, but the actual box size was determined from the unit cell info stored in the dcd file or in the CRYST line in the pdb file.
>>
>> Any insight into what might be going on here would be most appreciated.
>>
>> Josh
>>
>>
>
> --
> NIH Center for Macromolecular Modeling and Bioinformatics
> Beckman Institute for Advanced Science and Technology
> University of Illinois, 405 N. Mathews Ave, Urbana, IL 61801
> http://www.ks.uiuc.edu/~johns/ Phone: 217-244-3349
> http://www.ks.uiuc.edu/Research/vmd/