Metropolis Code in a Data Parallel Language

C----------------------------------------------------------------------C
C                                                                      C
C    Parallel Metropolis algorithm for 2-dimensional Ising model.      C
C                      CMFortran version.                              C
C----------------------------------------------------------------------C


      SUBROUTINE metupdate( black, spin )
      LOGICAL black(,)
      INTEGER spin(,) new_spin(,)
      INTEGER newE(,), oldE(,), deltaE(,)
      INTEGER checkerboard
      REAL rand(,)
C
C  Generate random numbers
C
      CALL RANDOM(rand)
C
C  New trial spin (= - old spin)
C
      new_spin = -spin

C
C  Red/black checkerboard update
C
      DO checkerboard = 1,2
C
C  Calculate old and new energy (new energy = - old energy)
C  by summing over neighboring spins.
C
        WHERE (black)
          oldE = - spin * ( CSHIFT(spin,1,1) + CSHIFT(spin,1,-1) +
     &                      CSHIFT(spin,2,1) + CSHIFT(spin,2,-1) )
          newE = -oldE
          deltaE = newE - oldE
      END WHERE
C
C  Metropolis accept/reject 
C 
        WHERE (black.AND.((deltaE.LE.0).OR.(EXP(-beta*deltaE).GT.rand)))
          spin = new_spin
        END WHERE
C
C  Change from black to red sites in checkerboard update
C
        black = .NOT.black

      END DO

      RETURN
      END