Why subtract 0.0?
> 67 call pbc(x,y,z,boxl)
> 68 end do
> 69 !- - - - - - - - - - - - - - - - - - - - - - - -
> 70 ! KEEPING MINIMUM DIST. BETWEEN ANY TWO
> 71 ! ATOM GREATER THEN RCUT
> 72 !- - - - - - - - - - - - - - - - - - - - - - - -
rcut2 = rcut**2
> 73 do i=1,nap-1
> 74 do j=i+1,nap
> 75 dx=x(i)-x(j)
> 76 dy=y(i)-y(j)
> 77 dz=z(i)-z(j)
> 78
> 79 dist=dsqrt(dx*dx+dy*dy+dz*dz)
Depending on the value of nap, it may be prudent to do
dist2 = dx*dx+dy*dy+dz*d
> 80 tr=tr+1
> 81 ! write(*,*) tr
> 82 if (dist<rcut)go to 10
if (dist2 < rcut2) goto 10
> 83 end do
> 84 end do
> 85 write(*,'(1x,"Success after",1x,i0,1x,"try")')tr
> 86 call mindist(x,y,z)
> 87 write(*,*) &
> 88
> "==========================================================="
> 89 write(*,*) ""
> 90
> 91 !-----------------------------------------------
> 92 ! INITIAL CONFIGURATION DONE
> 93 !-----------------------------------------------
> while compiling with ifort, its converging after 160k steps, but even
> after waiting for 1000k steps with gfortran, its not
> converging......do you think its problem with ran? can you suggest a
> better way?
What is seed? What is boxl()?
program a
integer i, seed
seed = 0
do i = 1, 4
print *, ran(seed)
end do
print *
seed = 1
do i = 1, 4
print *, ran(seed)
end do
seed = 42
do i = 1, 4
print *, ran(seed)
end do
end program a
RTFM?
> while compiling with ifort, its converging after 160k steps, but even
> after waiting for 1000k steps with gfortran, its not
> converging......do you think its problem with ran? can you suggest a
> better way?
Are you certain that it necessarily will converge? Are you getting the
same random sequence from the two compilers?
What have you done to try to characterize the nonconvergence? Do you
even have the slightest idea where the two compilers are diverging? I
would strongly recommend that you study that, using write statements or
a debugger.
Nothing like deleting all context from your follow up.
When I asked 'what is seed?', I literally meant give
us the value. Did you actually RTFM, which tells you
how to use RAN() in gfortran? When I asked 'What is
boxl?' Either show us the code or at least tell us
what boxl() does.
>> while compiling with ifort, its converging after 160k steps, but even
>> after waiting for 1000k steps with gfortran, its not
>> converging......do you think its problem with ran? can you suggest a
>> better way?
> Are you certain that it necessarily will converge? Are you getting the
> same random sequence from the two compilers?
I agree. A program could be very sensitive to roundoff such that
it doesn't converge or doesn't seem to converge. In some cases,
the result can alternate between values differing in the low order
but, and so will 'fail to converge' if compared for equality.
If it diverges (goes to infinity) then that is very different.
-- glen
> differing in the low order but ...
Some form of fuzzy logic?
Regards,
Mike Metcalf
>>differing in the low order but ...
> Some form of fuzzy logic?
http://en.wikipedia.org/wiki/Periodic_point
http://en.wikipedia.org/wiki/Limit_set
When using something like Newton-Raphson, it is non
unusual that the limit, even for a scalar, has a
period greater than one. For multi-dimensional systems
the result can be much more complicated.
Since the OP didn't give many details, it is hard to know,
but it might be that on one system it does converge to
a fixed point (period one cycle), and on others it doesn't.
-- glen
It looks like it's just trying to initialize a particle simulation
keeping a minimum distance between particles using periodic boundary
conditions. That's something where the ease of converging it is going
to depend on the particle density. If it needs 100k steps to succeed on
ifort, it may be dense enough that it's difficult to find a good
configuration, and the problem isn't due to one compiler or another but
rather just that the random sequences aren't the same (which comes back
to my question about whether the OP is certain that the random sequences
are the same, and what's been done to characterize the divergence
between the two compilers).
(Incidentally, that also suggests that it may be appropriate to find
another way to create an initial configuration... e.g. start out with
some sort of regular configuration and evolve it.)
I didn't even really look at the code in any detail initially, because
the basic steps to dealing with the problem are independent of the
details of the code. I will note, however, that if performance is any
concern, the OP should *definitely* be comparing distsq to rcutsq
instead of dist to rcut, as the sqrt taken to get dist will be very
expensive and rcutsq can be precomputed.
Oh, and one other thought... you shouldn't need to be doing PBC on your
initial positions. Fix your random assignments so that they land in the
box correctly.
Actually, a third thought... any particular reason why you do the entire
configuration at once and throw the whole thing out if there's an
incursion? You could also build it up incrementally, redoing an
individual particle if it's inside rcut.
Regards,
Mike Metcalf
>Actually, a third thought... any particular reason why you do the entire
>configuration at once and throw the whole thing out if there's an
>incursion? You could also build it up incrementally, redoing an
>individual particle if it's inside rcut.
That would be a much bigger improvement than a few sqrts.
-- greg
!===========================================
! This is the main driver routine
! of the reverse monte carlo program
program rmonte
!===========================================
implicit none
integer::i,j
integer:: seed,tr,rmcs
integer::ntyp,nap,ntot,it1,it2,atnm
real(selected_real_kind(8))::rmin,rmax,boxl,dr
real(selected_real_kind(8))::rr,rdf
real(selected_real_kind(8)),dimension(100)::x,y,z
real(selected_real_kind(8)),dimension(100):: gr,rrl
real(selected_real_kind(8)):: dx,dy,dz,dist,rcut
character(2)::elmnt
!
real::t1,t2
t1=secnds(t1)
!
seed=12345
BOXL=20
RCUT=1.5
DR=5.0
RMCS=1000
ntot=100
tr=1
write(*,*) &
"===============Creating Initial Configuration=============="
write(*,&
'(1x,"Minimum distance between any two particle is",1x,f6.4)'),
&
rcut
!-----------------------------------------------
! INITIALISE POSITION
!-----------------------------------------------
10 do i=1,ntot
x(i)=boxl*(ran(seed)-0.0)
y(i)=boxl*(ran(seed)-0.0)
z(i)=boxl*(ran(seed)-0.0)
call pbc(x,y,z,boxl)
end do
!- - - - - - - - - - - - - - - - - - - - - - - -
! KEEPING MINIMUM DIST. BETWEEN ANY TWO
! ATOM GREATER THEN RCUT
!- - - - - - - - - - - - - - - - - - - - - - - -
do i=1,ntot-1
do j=i+1,ntot
dx=x(i)-x(j)
dy=y(i)-y(j)
dz=z(i)-z(j)
dist=dsqrt(dx*dx+dy*dy+dz*dz)
tr=tr+1
! write(*,*) tr
if (dist<rcut)go to 10
end do
end do
write(*,'(1x,"Success after",1x,i0,1x,"try")')tr
call mindist(x,y,z)
write(*,*) &
"==========================================================="
write(*,*) ""
!-----------------------------------------------
! INITIAL CONFIGURATION DONE
!-----------------------------------------------
t2=secnds(t1)
write(*,'("Elapsed Time =",1x,f10.4)'), t2
end
!===========================================
! Subroutine to apply periodic
! boundary value
!===========================================
subroutine pbc(xx,yy,zz,boxl)
implicit none
real(selected_real_kind(4))::xx,yy,zz,boxl
if(xx.gt. boxl)xx=xx-boxl
if(xx.lt.-boxl)xx=xx+boxl
if(yy.gt. boxl)yy=yy-boxl
if(yy.lt.-boxl)yy=yy+boxl
if(zz.gt. boxl)zz=zz-boxl
if(zz.lt.-boxl)zz=zz+boxl
return
end
!==============================================
! Subroutine to calculate minimum
! distance between two atom
subroutine mindist(x,y,z)
!==============================================
implicit none
integer::i,j,p1,p2
real(selected_real_kind(8)),dimension(100):: x,y,z
real(selected_real_kind(8)):: dist,dx,dy,dz,dmini,dmin
write(*, &
'(" Calculating minimum distance between atoms....",$)')
dx=x(1)-x(2)
dy=y(1)-y(2)
dz=z(1)-z(2)
dmini=dsqrt(dx*dx+dy*dy+dz*dz)
do i=1,99
do j=i+1,100
dx=x(i)-x(j)
dy=y(i)-y(j)
dz=z(i)-z(j)
dmin=dsqrt(dx*dx+dy*dy+dz*dz)
if (dmini>dmin)then
dmini=dmin
p1=i
p2=j
end if
end do
end do
write(*,*) "DONE"
write &
(*,'(1x,"Shortest distance between to atom is",2x,f8.6)')dmini
write(*,'(1x,"between particle pair",1x,i0," - ",i0)') p1,p2
return
end
while running with ifort, its prompting:
$ ./a.out
===============Creating Initial Configuration==============
Minimum distance between any two particle is 1.5000
Success after 1680724 try
Calculating minimum distance between atoms.... DONE
Shortest distance between to atom is 1.626997
between particle pair 40 - 92
===========================================================
Elapsed Time = 0.0623
but while running with gfortran, its never success...why?
To quote my previous suggestions:
Are you certain that it necessarily will converge? Are you getting the
same random sequence from the two compilers?
What have you done to try to characterize the nonconvergence? Do you
Dick Hendrickson