Topic Options
#16295 - 11/19/07 10:11 AM Cell Dimensions
adewar Offline
Forum Member

Registered: 07/27/06
Posts: 6
I have produced a trajectory file and would like to extract the cell dimensions from each frame. I think I have to set up the PBC's as per the simulation, but I do not know how to extract the X, Y and Z cell dimensions from each individual frame, how would I go about this and is it possible?

Here is what I have so far;

Code:
trajectory query unit 51
trajectory iread 51 begin ?start skip ?skip

set i 1

label loop
trajectory read
SET A 23.00000
SET B 23.00000
SET C 23.00000
SHAKE BOND
CRYSTAL DEFINE CUBIC @A @B @C 90. 90. 90.
CRYSTAL BUILD CUTOFF 7.0 NOPERATIONS 0

set list @i

write coor card name frame-@list.crd
* @LX @LY @LZ
*
increment i
if i le ?nfile goto loop



With Thanks

Alastair

Top
#16296 - 11/19/07 10:19 AM Re: Cell Dimensions [Re: adewar]
rmv Offline

Forum Member

Registered: 09/17/03
Posts: 8379
Loc: 39 03 48 N, 77 06 54 W
The easy way is to use the CELL time series type in CORREL; see correl.doc

It can also be done with TRAJ READ, but CRYSTAL should only be set up once, before the loop; after the TRAJ READ command, the CHARMM internal variables ?XTLA etc. will contain the unit cell info.

Top
#16297 - 11/19/07 02:07 PM Re: Cell Dimensions [Re: rmv]
adewar Offline
Forum Member

Registered: 07/27/06
Posts: 6
Dr Venable,

I have changed the input to the following:

Code:
 
trajectory iread 51 begin ?start skip ?skip

set i 1
SET A 23.00000
SET B 23.00000
SET C 23.00000
CRYSTAL DEFINE CUBIC @A @B @C 90. 90. 90.
!CRYSTAL BUILD CUTOFF 7.0 NOPERATIONS 0

label loop

trajectory read
set list @i
set LX ?xlta

write coor card name methane_water-@list.crd
*
increment i
if i le ?nfile goto loop




Is this the sequence of commands as you suggested? I have also tried to use the set command to get the cell dimensions, but get an error from RDCMND: cannot substitute energy “?XLTA”. I suspect I'm using the wrong command and/or omitting crucial commands here, but do not know what I should use due to my limited knowledge of charmm. I've also had a look for examples of what I'm trying to do, but cannot find any specific examples.

I've also had a read of correl.doc and think I can see how to get the cell dimensions via this route using the following:

Code:
 
CORREL MAXSERIES 3 MAXTIMESTEP 3
ENTER LXCD CELL A
ENTER LYCD CELL B
ENTER LZCD CELL C
trajectory iread firstu 51 nunit 1 begin ?start skip ?skip




The correl function as it is set-up above reads the entire trajectory and calculates the average X, Y and Z from all the time steps. Is there anyway of getting charmm to output the X, Y and Z from each individual frame, as opposed to an average of all time steps? And how can I convert the LXCD (and other) functions to a simple variable, which I can then write to the frame-@list.crd files?

With Thanks in advance

Alastair

Top
#16298 - 11/19/07 02:24 PM Re: Cell Dimensions [Re: adewar]
lennart Offline

Forum Member

Registered: 09/25/03
Posts: 4742
Loc: ~ 59N, 15E
1/ Correct spelling is important.
2/ Yes you can output timseries data from the CORREL module; see examples in the Script Archive, and in correl.doc - you apparently missed the examples showing the normal usage of this module, in which the resulting timeseries are written to disk.
3/ See subst.doc
_________________________
Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden

Top
#16299 - 11/19/07 02:33 PM Re: Cell Dimensions [Re: adewar]
rmv Offline

Forum Member

Registered: 09/17/03
Posts: 8379
Loc: 39 03 48 N, 77 06 54 W
After the TRAJ command (omit the IREAD), write out the time series to a file; for this case, a simple way is

open unit 4 write card name cellabc.dat
write all dumb time unit 4
* dummy
*


This will result in a text file with one line per time point and 4 data columns:

Time CellA CellB CellC

You'll probably want to increase MAXTimesteps a good deal, and I recommend the use of NOUPDATE when starting CORREL; the exception is when you want to compute an energy time series, which isn't very often (at least for me).

Top

Moderator:  John Legato, lennart