6.5 Program 6.1, Fortran Code
99
vc = sqrt( mvcx**2 + mvcy**2 + mvcz**2)
t = dt*step
if(step.GT.7.E5) then
step1 = step - 7.E5
t2 = dt*step1
write(10,*) t2, mvcx
endif
c mvcx = sumvcx/step1
c mvcx2 =sumvcx2 /step1
c———————-FORNES—————DRXSQ = (rx(1)-rx(2))**2
DRYSQ = (ry(1)-ry(2))**2
DRZSQ = (rz(1)-rz(2))**2
R = SQRT(DRXSQ + DRYSQ + DRZSQ)
RSQ = R*R
SUMR = SUMR + R
SUMR2 = SUMR2 + RSQ
do i = 1, n
rx0(i) = rx(i)
ry0(i) = ry(i)
rz0(i) = rz(i)
end do
rcx0 = (rx0(1)+rx0(2))*0.5
rcy0 = (ry0(1)+ry0(2))*0.5
rcz0 = (rz0(1)+rz0(2))*0.5
rx210 = rx0(2)-rx0(1)
ry210 = ry0(2)-ry0(1)
rz210 = rz0(2)-rz0(1)
c SUZUKI EFFICIENCY
if(step.GT.7.E5.AND.step.LT.700503.) then
step1 = step - 7.E5
t1 = dt*step1
Frat = Fo*(cos(DPI*rcx) -0.5*cos(2*DPI*rcx))
sumFrat = sumFrat + Frat*vcx
endif
100 continue
RMEAN = SUMR/real(NSTEP)
RMEANSQ = SUMR2/real(NSTEP)
DELTAR = SQRT(abs(RMEANSQ - RMEAN**2))
DEFF = acdiff/real(nstep)/6/DT
DEFFX = acdiffx/real(nstep)/2/DT
PECLET = ABS(VCX)/DEFFX
VCX2EQ = DEFFX*gamma/2.
c write(*,*) VCX2EQ, mvcx2
c ========================================================