Question on Molecular dynamics programming

[caseyjay]

Dear all,

I am not sure if this is the right place to ask the following. (Moderators, please move this thread to the right sub-forum if necessary.)

The question I am having is related to molecular dynamics. In the book "Molecular Dynamics Simulation: Elementary Methods" by J. M. Haile, a FORTRAN program is provided page 459. A subroutine called SSX2SC is listed and is described as

"the routine to be called only once during a run. At the start of the initialization it is used to scale the initial forces to the program units for acceleration. The scale factor is STPSQH = $$(\Delta t)^2 / 2$$."

In the FORTRAN program, it goes on to equate:

X2(I)=FX(I)*STPSQH
Y2(I)=FY(I)*STPSQH
Z2(I)=FZ(I)*STPSQH

where X2, Y2 and Z2 are accelerations in x, y and z directions respectively and FX, FY and FZ are forces in x, y and z directions respectively.

The program is treating the problem using reduced units such that mass m = 1, time =$$\sigma\sqrt{m/\epsilon}$$, energy $$\epsilon$$ and length $$\sigma$$.


What I do not understand is:

1. What is the meaning of "to scale the initial forces to the program units for acceleration"?
2. Why multiply force by STPSQH?

Thank you very much and I hope someone can shed some light to this.

Casey

 

[nbo10]

Can you fix the LaTeX?

 

[caseyjay]

Hi nbo10,

I am not sure if there is something wrong with LATEX or not. Anyway I have duplicate the earlier question over here but this time without using LATEX in it.

The question I am having is related to molecular dynamics. In the book "Molecular Dynamics Simulation: Elementary Methods" by J. M. Haile, a FORTRAN program is provided page 459. A subroutine called SSX2SC is listed and is described as

"the routine to be called only once during a run. At the start of the initialization it is used to scale the initial forces to the program units for acceleration. The scale factor is STPSQH = [(delta t)^2]/2."

In the FORTRAN program, it goes on to equate:

X2(I)=FX(I)*STPSQH
Y2(I)=FY(I)*STPSQH
Z2(I)=FZ(I)*STPSQH

where X2, Y2 and Z2 are accelerations in x, y and z directions respectively and FX, FY and FZ are forces in x, y and z directions respectively.

The program is treating the problem using reduced units such that mass m = 1, time =sigma*sqrt(m/episilon), energy is epsilon and length is sigma.

What I do not understand is:

1. What is the meaning of "to scale the initial forces to the program units for acceleration"?
2. Why multiply force by STPSQH?

Thank you very much and I hope someone can shed some light to this.

Casey

 

[alxm]

The program is treating the problem using reduced units such that mass m = 1, time =sigma*sqrt(m/episilon), energy is epsilon and length is sigma.

What I do not understand is:

1. What is the meaning of "to scale the initial forces to the program units for acceleration"?
2. Why multiply force by STPSQH?
Casey

Seems to me you answered your own question;
F = ma, so if m=1 then F = a in the program's internal units, where STPSQH is the scaling factor to convert to those units.

 

[caseyjay]

Seems to me you answered your own question;
F = ma, so if m=1 then F = a in the program's internal units, where STPSQH is the scaling factor to convert to those units.

Hi alxm,

Thanks for your reply, but why multiply it with STPSQH [(delta t)^2]/2?

What do you mean by "scaling factor" and to convert to "what" units?

 

[t!m]

The factor, STPSQH (probably referring to something like "[time]Step SQuared Halved") comes about from the integrator you are using. Many are derived using some kind of short time Taylor expansion to update the position of the particles, i.e.
x(dt)=x(0)+x'(0)dt+x''(0)(dt)^2/2+...
Clearly x'(0) is the velocity and x''(0) is the acceleration (here equal to the force due to m=1) and the (dt)^2/2 term is the one you are asking about. The variables X2, Y2, Z2, would be more appropriately called "the coordinate contribution due to acceleration", because they clearly have units of position, after the "scaling factor" is applied. In this way, the coordinates can be updated by something like:
X(I) = X1(I) + X2(I)