BraWl
Loading...
Searching...
No Matches
random_site Module Reference

Functions/Subroutines

pure integer function, dimension(4), public simple_cubic_random_site (setup)
 Function to get a random site on the simple cubic lattice.
 
pure integer function, dimension(4), public simple_cubic_random_nbr (setup, site)
 Function to get a random neighbour of a site on the simple cubic lattice.
 
pure integer function, dimension(4), public bcc_random_site (setup)
 Function to get a random site on the bcc lattice.
 
pure integer function, dimension(4), public bcc_random_nbr (setup, site)
 Function to get a random neighbour of a site on the bcc lattice.
 
pure integer function, dimension(4), public fcc_random_site (setup)
 Function to get a random site on the fcc lattice.
 
pure integer function, dimension(4), public fcc_random_nbr (setup, site)
 Function to get a random neighbour of a site on the fcc lattice.
 
subroutine, public pair_swap (config, idx1, idx2)
 Function to swap a pair of lattice site occupancies.
 

Variables

integer, dimension(3, 6), parameter sc_nbrs = reshape((/ 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, -1, 0, -1, 0, -1, 0, 0 /), (/3, 6/))
 
integer, dimension(3, 8), parameter bcc_nbrs = reshape((/ 1, 1, 1, 1, 1, -1, 1, -1, 1, 1, -1, -1, -1, 1, 1, -1, 1, -1, -1, -1, 1, -1, -1, -1 /), (/3, 8/))
 
integer, dimension(3, 12), parameter fcc_nbrs = reshape((/ 0, 1, 1, 0, 1, -1, 0, -1, 1, 0, -1, -1, 1, 1, 0, 1, -1, 0, 1, 0, 1, 1, 0, -1, -1, 1, 0, -1, -1, 0, -1, 0, 1, -1, 0, -1 /), (/3, 12/))
 

Function/Subroutine Documentation

◆ bcc_random_nbr()

pure integer function, dimension(4), public random_site::bcc_random_nbr ( class(run_params), intent(in) setup,
integer, dimension(4), intent(in) site )

Function to get a random neighbour of a site on the bcc lattice.

Author
C. D. Woodgate
Date
2019-2025
Parameters
setupDerived type containing simulation parameters
siteSite for which to get the neighbour
Returns
Indices of a random neighbour of the input site

Definition at line 152 of file random_site.f90.

153
154 class(run_params), intent(in) :: setup
155 integer, dimension(4), intent(in) :: site
156 integer, dimension(4) :: nbr
157 integer :: n
158
159 ! Generate random integer corresponding to neighbour
160 n = floor(8.0_real64*genrand())+1
161
162 ! Make this the neighbour
163 nbr(2:) = site(2:) + bcc_nbrs(:,n)
164
165 ! Wrap coordinates into box
166 nbr(1) = 1
167 nbr(2) = modulo(nbr(2)-1, 2*setup%n_1) + 1
168 nbr(3) = modulo(nbr(3)-1, 2*setup%n_2) + 1
169 nbr(4) = modulo(nbr(4)-1, 2*setup%n_3) + 1
170

◆ bcc_random_site()

pure integer function, dimension(4), public random_site::bcc_random_site ( class(run_params), intent(in) setup)

Function to get a random site on the bcc lattice.

Author
C. D. Woodgate
Date
2019-2025
Parameters
setupDerived type containing simulation parameters
Returns
Indices of a random site on the bcc lattice

Definition at line 127 of file random_site.f90.

128
129 class(run_params), intent(in) :: setup
130 integer, dimension(4) :: site
131
132 site(1) = 1
133 site(4) = floor(2.0_real64*genrand()*real(setup%n_3, real64)) +1
134 site(2) = 2 * floor(genrand()*real(setup%n_1, real64)) &
135 + 2 - modulo(site(4), 2)
136 site(3) = 2 * floor(genrand()*real(setup%n_2, real64)) &
137 + 2 - modulo(site(4), 2)
138

◆ fcc_random_nbr()

pure integer function, dimension(4), public random_site::fcc_random_nbr ( class(run_params), intent(in) setup,
integer, dimension(4), intent(in) site )

Function to get a random neighbour of a site on the fcc lattice.

Author
C. D. Woodgate
Date
2019-2025
Parameters
setupDerived type containing simulation parameters
siteSite for which to get the neighbour
Returns
Indices of a random neighbour of the input site

Definition at line 206 of file random_site.f90.

207
208 class(run_params), intent(in) :: setup
209 integer, dimension(4), intent(in) :: site
210 integer, dimension(4) :: nbr
211 integer :: n
212
213 ! Generate random integer corresponding to neighbour
214 n = floor(12.0_real64*genrand())+1
215
216 ! Make this the neighbour
217 nbr(2:) = site(2:) + fcc_nbrs(:,n)
218
219 ! Wrap coordinates into box
220 nbr(1) = 1
221 nbr(2) = modulo(nbr(2)-1, 2*setup%n_1) + 1
222 nbr(3) = modulo(nbr(3)-1, 2*setup%n_2) + 1
223 nbr(4) = modulo(nbr(4)-1, 2*setup%n_3) + 1
224

◆ fcc_random_site()

pure integer function, dimension(4), public random_site::fcc_random_site ( class(run_params), intent(in) setup)

Function to get a random site on the fcc lattice.

Author
C. D. Woodgate
Date
2019-2025
Parameters
setupDerived type containing simulation parameters
Returns
Indices of a random site on the fcc lattice

Definition at line 182 of file random_site.f90.

183
184 class(run_params), intent(in) :: setup
185 integer, dimension(4) :: site
186
187 site(1) = 1
188 site(4) = floor(2.0_real64*genrand()*real(setup%n_3, real64)) + 1
189 site(2) = floor(2.0_real64*genrand()*real(setup%n_1, real64)) + 1
190 site(3) = 2 * floor(genrand()*real(setup%n_2, real64)) + 1 &
191 + modulo((site(2)-modulo(site(4),2)), 2)
192

◆ pair_swap()

subroutine, public random_site::pair_swap ( integer(array_int), dimension(:,:,:,:) config,
integer, dimension(4), intent(in) idx1,
integer, dimension(4), intent(in) idx2 )

Function to swap a pair of lattice site occupancies.

Author
C. D. Woodgate
Date
2019-2025
Parameters
configSystem configuration
idx1Indices of first lattice site
idx2Indices of second lattice site
Returns
None

Definition at line 238 of file random_site.f90.

239
240 integer(array_int), dimension(:,:,:,:) :: config
241 integer, dimension(4), intent(in) :: idx1, idx2
242 integer(array_int) :: species1, species2
243
244 species1 = config(idx1(1), idx1(2), idx1(3), idx1(4))
245 species2 = config(idx2(1), idx2(2), idx2(3), idx2(4))
246
247 config(idx1(1), idx1(2), idx1(3), idx1(4)) = species2
248 config(idx2(1), idx2(2), idx2(3), idx2(4)) = species1
249
Here is the caller graph for this function:

◆ simple_cubic_random_nbr()

pure integer function, dimension(4), public random_site::simple_cubic_random_nbr ( class(run_params), intent(in) setup,
integer, dimension(4), intent(in) site )

Function to get a random neighbour of a site on the simple cubic lattice.

Author
C. D. Woodgate
Date
2019-2025
Parameters
setupDerived type containing simulation parameters
siteSite for which to get the neighbour
Returns
Indices of a random neighbour of the input site

Definition at line 97 of file random_site.f90.

98
99 class(run_params), intent(in) :: setup
100 integer, dimension(4), intent(in) :: site
101 integer, dimension(4) :: nbr
102 integer :: n
103
104 ! Generate random integer corresponding to neighbour
105 n = floor(6.0_real64*genrand())+1
106
107 ! Make this the neighbour
108 nbr(2:) = site(2:) + sc_nbrs(:,n)
109
110 ! Wrap coordinates into box
111 nbr(1) = 1
112 nbr(2) = modulo(nbr(2)-1, 2*setup%n_1) + 1
113 nbr(3) = modulo(nbr(3)-1, 2*setup%n_2) + 1
114 nbr(4) = modulo(nbr(4)-1, 2*setup%n_3) + 1
115

◆ simple_cubic_random_site()

pure integer function, dimension(4), public random_site::simple_cubic_random_site ( class(run_params), intent(in) setup)

Function to get a random site on the simple cubic lattice.

Author
C. D. Woodgate
Date
2019-2025
Parameters
setupDerived type containing simulation parameters
Returns
Indices of a random site on the simple cubic lattice

Definition at line 74 of file random_site.f90.

75
76 class(run_params), intent(in) :: setup
77 integer, dimension(4) :: site
78
79 site(1) = 1
80 site(2) = floor(genrand()*2.0_real64*real(setup%n_1)) +1
81 site(3) = floor(genrand()*2.0_real64*real(setup%n_2)) +1
82 site(4) = floor(genrand()*2.0_real64*real(setup%n_3)) +1
83

Variable Documentation

◆ bcc_nbrs

integer, dimension(3,8), parameter random_site::bcc_nbrs = reshape((/ 1, 1, 1, 1, 1, -1, 1, -1, 1, 1, -1, -1, -1, 1, 1, -1, 1, -1, -1, -1, 1, -1, -1, -1 /), (/3, 8/))
private

Definition at line 38 of file random_site.f90.

38 integer, parameter, dimension(3,8) :: &
39 bcc_nbrs = reshape((/ 1, 1, 1, &
40 1, 1, -1, &
41 1, -1, 1, &
42 1, -1, -1, &
43 -1, 1, 1, &
44 -1, 1, -1, &
45 -1, -1, 1, &
46 -1, -1, -1 /), (/3, 8/))

◆ fcc_nbrs

integer, dimension(3,12), parameter random_site::fcc_nbrs = reshape((/ 0, 1, 1, 0, 1, -1, 0, -1, 1, 0, -1, -1, 1, 1, 0, 1, -1, 0, 1, 0, 1, 1, 0, -1, -1, 1, 0, -1, -1, 0, -1, 0, 1, -1, 0, -1 /), (/3, 12/))
private

Definition at line 49 of file random_site.f90.

49 integer, parameter, dimension(3,12) :: &
50 fcc_nbrs = reshape((/ 0, 1, 1, &
51 0, 1, -1, &
52 0, -1, 1, &
53 0, -1, -1, &
54 1, 1, 0, &
55 1, -1, 0, &
56 1, 0, 1, &
57 1, 0, -1, &
58 -1, 1, 0, &
59 -1, -1, 0, &
60 -1, 0, 1, &
61 -1, 0, -1 /), (/3, 12/))

◆ sc_nbrs

integer, dimension(3,6), parameter random_site::sc_nbrs = reshape((/ 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, -1, 0, -1, 0, -1, 0, 0 /), (/3, 6/))
private

Definition at line 29 of file random_site.f90.

29 integer, parameter, dimension(3,6) :: &
30 sc_nbrs = reshape((/ 0, 0, 1, &
31 0, 1, 0, &
32 1, 0, 0, &
33 0, 0, -1, &
34 0, -1, 0, &
35 -1, 0, 0 /), (/3, 6/))