Search code examples
gfortrando-loopschemistry

FORTRAN counter loop returns multiple iterations of the same value


First of all I am a complete novice to FORTRAN. With that said I am attempting to "build" a box, then randomly generate x, y, z coordinates for 100 atoms. From there, the goal is to calculate the distance between each atom, which becomes the value "r" of the Lennard-Jones potential energy equation. Then calculate the LJ potential, and finally sum the potential of the entire box. A previous question that I had asked about this project is here. The problem is that I get the same calculated value over and over and over again. My code is below.

    program energytot
    implicit none

    integer, parameter :: n = 100
    integer :: i, j, k, seed(12)
    double precision :: sigma, r, epsilon, lx, ly, lz
    double precision, dimension(n) :: x, y, z, cx, cy, cz
    double precision, dimension(n*(n+1)/2) :: dx, dy, dz, LJx, LJy, LJz
    sigma = 4.1
    epsilon = 1.7
    !Box length with respect to the axis
    lx = 15
    ly = 15
    lz = 15

    do i=1,12
        seed(i)=i+3
    end do
    !generate n random numbers for x, y, z
    call RANDOM_SEED(PUT = seed)
    call random_number(x)
    call random_number(y)
    call random_number(z)
    !convert random numbers into x, y, z coordinates
    cx = ((2*x)-1)*(lx*0.5)
    cy = ((2*y)-1)*(lx*0.5)
    cz = ((2*z)-1)*(lz*0.5)

    do j=1,n-1
        do k=j+1,n
            dx = ABS((cx(j) - cx(k)))
            LJx = 4 * epsilon * ((sigma/dx(j))**12 - (sigma/dx(j))**6)
            dy = ABS((cy(j) - cy(k)))
            LJy = 4 * epsilon * ((sigma/dy(j))**12 - (sigma/dy(j))**6)
            dz = ABS((cz(j) - cz(k)))
            LJz = 4 * epsilon * ((sigma/dz(j))**12 - (sigma/dz(j))**6)
        end do
    end do

    print*, dx
    end program energytot

Solution

  • What exactly is your question? What do you want your code to do, and what does it do instead?

    If you're having problems with the final print statement print*, dx, try this instead:

      print *, 'dx = '
      do i = 1, n * (n + 1) / 2
        print *, dx(i)
      end do
    

    It seems that dx is too big to be printed without a loop.

    Also, it looks like you're repeatedly assigning the array dx (and other arrays in the loop) to a single value. Try this instead:

      i = 0
      do j=1,n-1
          do k=j+1,n
              i = i + 1
              dx(i) = ABS((cx(j) - cx(k)))
          end do
      end do
    

    This way, each value cx(j) - cx(k) gets saved to a different element of dx, instead of overwriting previously saved values.