* Move particles (fish) in periodic box according to external force field (current) * Inputs : * d = grid dimension * fish = Indicates fish presence at the grid points * tsteps = number of time steps * Outputs: * plot of fish positions after each time step * * program fish3 include '/usr/cm5/cmf/include/cm/timer-fort.h' include '/usr/cm5/cmf/include/cm/CMF_defs.h' * Define the size of the display window integer, parameter :: m = 256 integer, parameter :: disp_interval = 1 * Max. Grid Size integer, parameter :: MAXG = 64 integer fish(MAXG,MAXG), nbrcnt(MAXG,MAXG) CMF$ layout fish(:news,:news),nbrcnt(:news,:news) integer d, i, tsteps, stride, limit, off real etime0, etime1 call CMF_DESCRIBE_ARRAY(fish) print*, 'Input...' print*, 'Give grid size (max 64) : ' read (*,*) d call cmf_randomize(100) call cmf_random(fish,2) tsteps = 200 * * Initialize the timer call CM_timer_clear(0) call CM_timer_clear(1) call CM_timer_start(0) call disp_init(m, d, stride, limit, off) * * Loop over time steps do i=1, tsteps print *,'i = ',i * * Get neighbor count call CM_timer_start(1) nbrcnt(2:d-1,2:d-1) = fish(1:d-2,2:d-1) + fish(1:d-2,1:d-2) + ! fish(1:d-2,3:d) + fish(2:d-1,1:d-2) + ! fish(2:d-1,3:d) + fish(3:d,2:d-1) + ! fish(3:d,1:d-2) + fish(3:d,3:d) * print *,nbrcnt(2:d-1,2:d-1) where ((nbrcnt(2:d-1,2:d-1) .eq. 3) .or. ( (fish(2:d-1,2:d-1) ! .eq. 1) .and. (nbrcnt(2:d-1,2:d-1) .eq. 2))) fish(2:d-1,2:d-1) = 1 elsewhere fish(2:d-1,2:d-1) = 0 end where call updateboundaries(fish, d) * print *,fish(1:d,1:d) call CM_timer_stop(1) * * Display output every disp_interval steps if (mod(i-1,disp_interval) .eq. 0) then call show_fish(fish, d, stride, limit, off, m) endif enddo * * Read the timer call CM_timer_stop(0) etime0 = CM_timer_read_elapsed(0) etime1 = CM_timer_read_elapsed(1) print *, ' Computation time =', etime1, ' sec' print *, ' Total Elapsed time =', etime0, ' sec' stop end subroutine disp_init(m, d, stride, limit, off) integer m, d, stride, limit, off * * Set parameters for display stride = m/d limit = d*stride off = 0 off = stride/2 * print *,'m = ',m * print *,'d = ',d * print *,'limit = ',limit * print *,'stride = ',stride * print *,'off = ',off end subroutine show_fish(fish, d, stride, limit, off, m) integer fish(:,:), d, stride, limit, off real show(m,m) show = 0 * * Map fish grid onto show show(off+1:off+limit:stride,off+1:off+limit:stride) = fish(1:d,1:d) * * The following enlarges each pixel for better visibility show(:,1:m-1) = show(:,1:m-1) + show(:,2:m) show(1:m-1,:) = show(1:m-1,:) + show(2:m,:) call display (show,m,m,MINVAL(show),MAXVAL(show)+0.0001) end subroutine updateboundaries(fish, d) integer fish(:,:), d fish(1,:) = fish(d-1,:) fish(d,:) = fish(2,:) fish(:,1) = fish(:,d-1) fish(:,d) = fish(:,2) end