Skip to content

Commit

Permalink
Minor updates
Browse files Browse the repository at this point in the history
  • Loading branch information
Michael-P-Allen committed Aug 20, 2022
1 parent d420850 commit b83c924
Show file tree
Hide file tree
Showing 2 changed files with 7 additions and 7 deletions.
2 changes: 1 addition & 1 deletion cluster.f90
Original file line number Diff line number Diff line change
Expand Up @@ -114,7 +114,7 @@ PROGRAM cluster

IF ( ALL ( done > 0 ) ) EXIT ! Loop until all done

i = MINLOC ( done, dim = 1 ) ! Find first zero (FINDLOC is not implemented in gfortran at the time of writing)
i = FINDLOC ( done, 0, dim = 1 ) ! Find first zero (FINDLOC was implemented in gfortran v9)
cluster_id = cluster_id + 1
WRITE ( unit=output_unit, fmt='(a,i5,a)', advance='no' ) 'Cluster ', cluster_id, ' = '
j = i
Expand Down
12 changes: 6 additions & 6 deletions ewald_module.f90
Original file line number Diff line number Diff line change
Expand Up @@ -210,8 +210,8 @@ FUNCTION pot_k_pm_ewald ( nk, n, r, q, kappa ) RESULT ( pot )

INTEGER :: ix, iy, iz
REAL :: kx, ky, kz, fx, fxy, fxyz, g, k_sq, dr
real, DIMENSION(:,:,:), ALLOCATABLE :: rho_sq ! Square modulus of charge density (0:sc-1,0:sc-1,0:sc-1)
real, dimension(:), allocatable :: kmesh ! k-values in wraparound convention (0:sc-1)
REAL, DIMENSION(:,:,:), ALLOCATABLE :: rho_sq ! Square modulus of charge density (0:sc-1,0:sc-1,0:sc-1)
REAL, DIMENSION(:), ALLOCATABLE :: kmesh ! k-values in wraparound convention (0:sc-1)

REAL, PARAMETER :: pi = 4.0*ATAN(1.0), twopi = 2.0*pi, fourpi = 4.0*pi, dk = twopi

Expand All @@ -220,7 +220,7 @@ FUNCTION pot_k_pm_ewald ( nk, n, r, q, kappa ) RESULT ( pot )
! Use nk to determine mesh dimension
sc = 2*nk
ALLOCATE ( fft_inp(0:sc-1,0:sc-1,0:sc-1), fft_out(0:sc-1,0:sc-1,0:sc-1) )
allocate ( rho_sq(0:sc-1,0:sc-1,0:sc-1), kmesh(0:sc-1) )
ALLOCATE ( rho_sq(0:sc-1,0:sc-1,0:sc-1), kmesh(0:sc-1) )

! Assign charge density to complex array ready for FFT
! Assume r in unit box with range (-0.5,0.5)
Expand All @@ -231,10 +231,10 @@ FUNCTION pot_k_pm_ewald ( nk, n, r, q, kappa ) RESULT ( pot )
CALL fftw_destroy_plan ( fft_plan ) ! Release plan

fft_out = fft_out / REAL(sc**3) ! Incorporate scaling by number of grid points
rho_sq = real ( fft_out*conjg(fft_out) ) ! Convert to square modulus of charge density
rho_sq = REAL ( fft_out*CONJG(fft_out) ) ! Convert to square modulus of charge density

kmesh = [ ( real(ix)*dk, ix = 0, sc-1 ) ] ! Set up k-components
kmesh(nk:) = kmesh(nk:) - real(sc)*dk ! in wraparound convention
kmesh = [ ( REAL(ix)*dk, ix = 0, sc-1 ) ] ! Set up k-components
kmesh(nk:) = kmesh(nk:) - REAL(sc)*dk ! in wraparound convention
pot = 0.0

! Triple loop over xyz grid points (uses wraparound indexing)
Expand Down

0 comments on commit b83c924

Please sign in to comment.