Page 1 of 2 1 2 >
Topic Options
#403 - 10/23/03 06:51 AM self-diffusion-coeff.inp
lennart Offline

Forum Member

Registered: 09/25/03
Posts: 4688
Loc: ~ 59N, 15E
*FILENAME: self-diffusion-coeff.inp
*PURPOSE: compute Mean Square Displacement vs time to get diffusion coeffcient D for water
* via the Einstein relation: 6*D is the slope of MSD vs t at long t.
*AUTHOR: Lennart Nilsson, Karolinska Institutet, October 2003
*
!unix environment variable CHM_HOME points to CHARMM installation directory
read rtf card name $CHM_HOME/toppar/top_all22_prot.inp
read para card name $CHM_HOME/toppar/par_all22_prot.inp
read psf card name my_psf.psf

open unit 51 read unform name my_traj1.cor
open unit 52 read unform name my_traj2.cor

open unit 21 write form name msd.dat
write title unit 21
* TIME MSD
*

coor anal select type oh2 end - ! what atoms to look at, these are the TIP3 oxygens
firstu 51 nunit 2 skip 50 - ! trajectory specification
imsd 21 - ! flag to do the MSD analysis
rspin 0.0 rspout 999.9 - ! we are interested in ALL waters
ncors 200 - ! compute MSD out to NCORS*SKIP steps
xbox 40.0 ybox 45.0 zbox 45.0 ! and we did use PBC, simple rectangular box
! if trajectory is from constant pressure simulation
! actual box size is taken from trajectory
_________________________
Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden

Top
#404 - 11/03/04 07:56 PM Re: self-diffusion-coeff.inp [Re: lennart]
chuan Offline
Forum Member

Registered: 05/17/04
Posts: 222
Hi Prof Nilsson,

Quote:

! if trajectory is from constant pressure simulation
! actual box size is taken from trajectory




Is it referring to @FILE-box-size.str? Or is there some way of getting the sizes instantaneously? from traj?
_________________________
CHARMM 30b1 driven by 1/ Xeon (32 bits) 2/ Redhat 7.3 (32 bits) with a Quadrics-modified 2.4-18-5 kernel 3/ Chuan, with 95% of the mentorship coming from great scientists frequenting this forum. 4/ Gracious support from the forum.

Top
#405 - 11/04/04 02:33 AM Re: self-diffusion-coeff.inp [Re: chuan]
lennart Offline

Forum Member

Registered: 09/25/03
Posts: 4688
Loc: ~ 59N, 15E
the statement "if trajectory is from constant pressure simulation
actual box size is taken from trajectory" means exactly what it says, namely that if you have a trajectory that was performed using constant pressure (in which case the box size may vary), the box size at each step (which in this case is stored in the trajectory) is taken from the trajectory. you do not have to do anything about it.
_________________________
Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden

Top
#406 - 10/06/05 02:29 PM Re: self-diffusion-coeff.inp [Re: lennart]
dmr Offline
Forum Member

Registered: 04/02/04
Posts: 9
Loc: University of Wisconsin, Madis...
Can this be modified for computing the self-diffusion within RHDO pbc? If so, would you be so kind as to reveal the modifications?

Thanks a bunch!

Top
#407 - 10/06/05 03:53 PM Re: self-diffusion-coeff.inp [Re: dmr]
lennart Offline

Forum Member

Registered: 09/25/03
Posts: 4688
Loc: ~ 59N, 15E
The current implementation of PBC corrections in the COOR ANAL code is very simple, and only applies to orthorhombic lattices.

The problems with computing diffusion coefficients using the Einstein relation in PBC are: 1/ mean square displacements (MSD) at long times are underestimated due to the limited distance particles can travel since they are confined to the box by the recentering, and 2/ you can get spuriously high MSDs when molecules are moved across the box by the recentering. It is the second effect that is corrected by the PBC code in COOR ANAL (the first effect is usually not a problem because you can get a good estimate of the slope of MSD vs t before it sets in).

You can try approximately three different approaches:
1/ Implement PBC corrections that can also handle RHDO in subroutine corfuncb of source/solana.src (you can see how the same thing is done in source/manip/hbanal.src in the most recent version of CHARMM). When done, and satisfied that it works, please submit changes. I will try to find the time to put this in for the next release.

2/ Undo the recentering events using MERGE UNFOLD on your trajectory before the analysis.

3/ Restrict the MSD calculations to molecules that are not in the vicinity of the RHDO borders (to avoid problems of the second kind): Set RSPOut to a value well inside your RHDO. This is the simplest way, but you probably have enough data to get a reasonable estimate of the diffusion coefficient anyway; whether this is the case is highly system-dependent.
_________________________
Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden

Top
#408 - 10/06/05 04:16 PM Re: self-diffusion-coeff.inp [Re: lennart]
dmr Offline
Forum Member

Registered: 04/02/04
Posts: 9
Loc: University of Wisconsin, Madis...
Thank you for the helpful response! I will consider implementation, if other research things go smoothly as to release my time constraints.

-Demian

Top
#409 - 10/06/05 05:33 PM Re: self-diffusion-coeff.inp [Re: dmr]
rmv Online   content

Forum Member

Registered: 09/17/03
Posts: 8259
Loc: 39 03 48 N, 77 06 54 W
The MERGE COOR UNFOLD operation is also currently restricted to orthonormal shapes (CUBIC, TETRagonal, ORTHorhombic).
_________________________
Rick Venable
computational chemist


Top
#410 - 07/12/06 02:02 PM Re: self-diffusion-coeff.inp [Re: lennart]
leo Offline
Forum Member

Registered: 10/04/04
Posts: 22
Dear All,

In the script it shows that the trajectory was saved every 50 steps. My question is how often we should save the trajectory for diffusion coefficient calculation of TIP3 water. Is there any rule or something?

What is more, I understand that only the linear part of MSD should be used for diffusion coefficient calculation. However, is there any recommend value of ncors to calculate the diffusion of water base on your experience?

Thanks,
Leo

Top
#411 - 07/12/06 05:08 PM Re: self-diffusion-coeff.inp [Re: leo]
lennart Offline

Forum Member

Registered: 09/25/03
Posts: 4688
Loc: ~ 59N, 15E
please read the relevant original literature!
for water you can get a good estimate of Dtransl for the 1-20 ps part of MSD(t)
_________________________
Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden

Top
#412 - 07/12/06 06:57 PM Re: self-diffusion-coeff.inp [Re: lennart]
leo Offline
Forum Member

Registered: 10/04/04
Posts: 22
Thanks, Lennart. I do appreciate your help.

-Leo

Top
Page 1 of 2 1 2 >

Moderator:  chmgr, John Legato, petrella