| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | module raffle__place_methods | ||
| 2 | !! Module containing the placement methods available within RAFFLE. | ||
| 3 | !! | ||
| 4 | !! This module contains procedures to query points for atom placement. | ||
| 5 | !! The available placement methods are: | ||
| 6 | !! - void: place the atom in the gridpoint with the largest void space | ||
| 7 | !! - rand: place the atom at a random gridpoint | ||
| 8 | !! - walk: place the atom using a random walk method | ||
| 9 | !! - growth: place the atom using a random walk, with last placement point | ||
| 10 | !! as the starting point | ||
| 11 | !! - min: place the atom at the gridpoint with the highest viability | ||
| 12 | use raffle__constants, only: real32, pi | ||
| 13 | use raffle__misc_linalg, only: inverse_3x3 | ||
| 14 | use raffle__geom_extd, only: extended_basis_type | ||
| 15 | use raffle__dist_calcs, only: & | ||
| 16 | get_min_dist, & | ||
| 17 | get_min_dist_between_point_and_species | ||
| 18 | use raffle__evaluator, only: evaluate_point | ||
| 19 | use raffle__distribs_container, only: distribs_container_type | ||
| 20 | implicit none | ||
| 21 | |||
| 22 | |||
| 23 | private | ||
| 24 | |||
| 25 | public :: & | ||
| 26 | place_method_void, place_method_rand, & | ||
| 27 | place_method_walk, place_method_growth, & | ||
| 28 | place_method_min | ||
| 29 | |||
| 30 | |||
| 31 | contains | ||
| 32 | |||
| 33 | !############################################################################### | ||
| 34 | 1 | function place_method_void( & | |
| 35 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | points, basis, viable & |
| 36 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | ) result(point) |
| 37 | !! VOID placement method. | ||
| 38 | !! | ||
| 39 | !! This method returns the gridpoint with the lowest neighbour density. | ||
| 40 | !! i.e. the point with the lowest density in the cell. | ||
| 41 | implicit none | ||
| 42 | |||
| 43 | ! Arguments | ||
| 44 | real(real32), dimension(:,:), intent(in) :: points | ||
| 45 | !! List of gridpoints to consider. | ||
| 46 | type(extended_basis_type), intent(inout) :: basis | ||
| 47 | !! Structure to add atom to. | ||
| 48 | logical, intent(out) :: viable | ||
| 49 | !! Boolean to indicate if point is viable. | ||
| 50 | real(real32), dimension(3) :: point | ||
| 51 | !! Point to add atom to. | ||
| 52 | |||
| 53 | ! Local variables | ||
| 54 | integer :: best_gridpoint | ||
| 55 | !! Index of best gridpoint. | ||
| 56 | |||
| 57 | |||
| 58 | 1 | viable = .false. | |
| 59 | |||
| 60 | !--------------------------------------------------------------------------- | ||
| 61 | ! find the gridpoint with the largest void space | ||
| 62 | !--------------------------------------------------------------------------- | ||
| 63 |
6/10✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 943 times.
✓ Branch 7 taken 1 times.
✓ Branch 8 taken 12 times.
✓ Branch 9 taken 931 times.
|
944 | best_gridpoint = maxloc(points(4,:), dim=1) |
| 64 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | if(best_gridpoint.eq.0) return |
| 65 | |||
| 66 | |||
| 67 | !--------------------------------------------------------------------------- | ||
| 68 | ! return the gridpoint with the largest void space | ||
| 69 | !--------------------------------------------------------------------------- | ||
| 70 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 times.
|
4 | point = points(1:3,best_gridpoint) |
| 71 | 1 | viable = .true. | |
| 72 | |||
| 73 | end function place_method_void | ||
| 74 | !############################################################################### | ||
| 75 | |||
| 76 | |||
| 77 | !############################################################################### | ||
| 78 | 1 | function place_method_rand( distribs_container, & | |
| 79 | bounds, & | ||
| 80 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | basis, species, radius_list, max_attempts, viable & |
| 81 | 1 | ) result(point) | |
| 82 | !! Random placement method. | ||
| 83 | !! | ||
| 84 | !! This method places the atom at a random gridpoint. | ||
| 85 | implicit none | ||
| 86 | |||
| 87 | ! Arguments | ||
| 88 | type(distribs_container_type), intent(in) :: distribs_container | ||
| 89 | !! Distribution function (gvector) container. | ||
| 90 | real(real32), dimension(2,3), intent(in) :: bounds | ||
| 91 | !! Bounds of the unit cell. | ||
| 92 | type(extended_basis_type), intent(inout) :: basis | ||
| 93 | !! Structure to add atom to. | ||
| 94 | integer, intent(in) :: species | ||
| 95 | !! Species index to add atom to. | ||
| 96 | real(real32), dimension(:), intent(in) :: radius_list | ||
| 97 | !! List of radii for each pair of elements. | ||
| 98 | integer, intent(in) :: max_attempts | ||
| 99 | !! Limit on number of attempts. | ||
| 100 | logical, intent(out) :: viable | ||
| 101 | !! Boolean to indicate if point is viable. | ||
| 102 | real(real32), dimension(3) :: point | ||
| 103 | !! Point to add atom to. | ||
| 104 | |||
| 105 | ! Local variables | ||
| 106 | integer :: i, is, js | ||
| 107 | !! Loop indices. | ||
| 108 | real(real32) :: rtmp1 | ||
| 109 | !! random number. | ||
| 110 | real(real32), dimension(3) :: rvec1 | ||
| 111 | !! random vector. | ||
| 112 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | integer, dimension(basis%nspec,basis%nspec) :: pair_index |
| 113 | |||
| 114 | |||
| 115 | 1 | viable = .false. | |
| 116 | |||
| 117 | !--------------------------------------------------------------------------- | ||
| 118 | ! get list of element pair indices | ||
| 119 | ! (i.e. the index for bond_info for each element pair) | ||
| 120 | !--------------------------------------------------------------------------- | ||
| 121 | 1 | i = 0 | |
| 122 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
|
2 | do is = 1, basis%nspec |
| 123 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
|
3 | do js = 1, basis%nspec |
| 124 | pair_index(is, js) = distribs_container%get_pair_index( & | ||
| 125 | basis%spec(is)%name, basis%spec(js)%name & | ||
| 126 | 2 | ) | |
| 127 | end do | ||
| 128 | end do | ||
| 129 | |||
| 130 | |||
| 131 | !--------------------------------------------------------------------------- | ||
| 132 | ! find a random gridpoint that is not too close to any other atom | ||
| 133 | !--------------------------------------------------------------------------- | ||
| 134 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | atom_loop: do i = 1, max_attempts |
| 135 | 1 | call random_number(rvec1) | |
| 136 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 times.
|
4 | point = bounds(1,:) + ( bounds(2,:) - bounds(1,:) ) * rvec1 |
| 137 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
|
2 | do js = 1, basis%nspec |
| 138 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | if( & |
| 139 | get_min_dist_between_point_and_species( & | ||
| 140 | basis, point, & | ||
| 141 | species = js & | ||
| 142 | ) .lt. max( & | ||
| 143 | distribs_container%cutoff_min(1), & | ||
| 144 | radius_list(pair_index(species,js)) * & | ||
| 145 | distribs_container%radius_distance_tol(1) & | ||
| 146 | ) & | ||
| 147 | 1 | )then | |
| 148 | ✗ | cycle atom_loop | |
| 149 | end if | ||
| 150 | end do | ||
| 151 | 1 | viable = .true. | |
| 152 | 1 | exit atom_loop | |
| 153 | end do atom_loop | ||
| 154 | |||
| 155 | 1 | end function place_method_rand | |
| 156 | !############################################################################### | ||
| 157 | |||
| 158 | |||
| 159 | !############################################################################### | ||
| 160 | 1 | function place_method_walk( distribs_container, & | |
| 161 | bounds, & | ||
| 162 | basis, species, & | ||
| 163 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | radius_list, max_attempts, & |
| 164 | step_size_coarse, step_size_fine, & | ||
| 165 | viable & | ||
| 166 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | ) result(point) |
| 167 | !! Random walk placement method. | ||
| 168 | !! | ||
| 169 | !! This method places the atom using a random walk method. | ||
| 170 | !! An initial point is chosen at random, and then points nearby are tested | ||
| 171 | !! to see if they are more suitable than the current point. If they are, | ||
| 172 | !! the query point is moved to the new point and the process is repeated. | ||
| 173 | !! The process is repeated, with each point being tested against a random | ||
| 174 | !! number. If the random number is less than the suitability of the point, | ||
| 175 | !! the atom is placed at that point. | ||
| 176 | implicit none | ||
| 177 | |||
| 178 | ! Arguments | ||
| 179 | type(distribs_container_type), intent(in) :: distribs_container | ||
| 180 | !! Distribution function (gvector) container. | ||
| 181 | real(real32), dimension(2,3), intent(in) :: bounds | ||
| 182 | !! Bounds of the unit cell. | ||
| 183 | type(extended_basis_type), intent(inout) :: basis | ||
| 184 | !! Structure to add atom to. | ||
| 185 | integer, intent(in) :: species | ||
| 186 | !! Species index to add atom to. | ||
| 187 | integer, intent(in) :: max_attempts | ||
| 188 | !! Limit on number of attempts. | ||
| 189 | real(real32), intent(in) :: step_size_coarse, step_size_fine | ||
| 190 | !! Step sizes for random walk. | ||
| 191 | logical, intent(out) :: viable | ||
| 192 | !! Boolean to indicate if point is viable. | ||
| 193 | real(real32), dimension(:), intent(in) :: radius_list | ||
| 194 | !! List of radii for each pair of elements. | ||
| 195 | real(real32), dimension(3) :: point | ||
| 196 | !! Point to add atom to. | ||
| 197 | |||
| 198 | ! Local variables | ||
| 199 | integer :: i, j | ||
| 200 | !! Loop indices. | ||
| 201 | integer :: nattempt, nstuck | ||
| 202 | !! Number of attempts and number of times stuck at same site | ||
| 203 | real(real32) :: rtmp1 | ||
| 204 | !! Random number. | ||
| 205 | real(real32), dimension(3) :: rvec1, abc | ||
| 206 | !! Random vector and lattice constants. | ||
| 207 | real(real32) :: crude_norm | ||
| 208 | !! Crude normalisation. | ||
| 209 | real(real32) :: site_value, test_value | ||
| 210 | !! Viability values. | ||
| 211 | real(real32), dimension(3) :: site_vector, test_vector | ||
| 212 | !! Vectors for gridpoints. | ||
| 213 | |||
| 214 | |||
| 215 | 1 | viable = .false. | |
| 216 | |||
| 217 | !--------------------------------------------------------------------------- | ||
| 218 | ! test a random point in the unit cell | ||
| 219 | !--------------------------------------------------------------------------- | ||
| 220 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 times.
|
4 | do i = 1, 3 |
| 221 |
5/6✓ Branch 0 taken 9 times.
✓ Branch 1 taken 3 times.
✓ Branch 2 taken 3 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
|
13 | abc(i) = norm2(basis%lat(i,:)) |
| 222 | end do | ||
| 223 | 1 | i = 0 | |
| 224 | 452 | random_loop : do | |
| 225 | 453 | i = i + 1 | |
| 226 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 453 times.
|
453 | if(i.gt.max_attempts) return |
| 227 | 453 | call random_number(site_vector) | |
| 228 |
2/2✓ Branch 0 taken 1359 times.
✓ Branch 1 taken 453 times.
|
1812 | site_vector = bounds(1,:) + ( bounds(2,:) - bounds(1,:) ) * site_vector |
| 229 | |||
| 230 | site_value = evaluate_point( distribs_container, & | ||
| 231 | site_vector, species, basis, radius_list & | ||
| 232 | 453 | ) | |
| 233 | 453 | call random_number(rtmp1) | |
| 234 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 452 times.
|
453 | if(rtmp1.lt.site_value) exit random_loop |
| 235 | |||
| 236 | end do random_loop | ||
| 237 | |||
| 238 | |||
| 239 | !--------------------------------------------------------------------------- | ||
| 240 | ! now do a random walk to find a suitable point to place the atom | ||
| 241 | !--------------------------------------------------------------------------- | ||
| 242 | 1 | nattempt = 0 | |
| 243 | 1 | nstuck = 0 | |
| 244 | 1 | crude_norm = 0.5_real32 | |
| 245 | 1 | i = 0 | |
| 246 | 24 | walk_loop : do | |
| 247 | 25 | i = i + 1 | |
| 248 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 25 times.
|
25 | if(i.gt.max_attempts) return |
| 249 | !------------------------------------------------------------------------ | ||
| 250 | ! if we have tried 10 times, then we need to reduce the step size | ||
| 251 | ! get the new test point and map it back into the unit cell | ||
| 252 | !------------------------------------------------------------------------ | ||
| 253 | 25 | call random_number(rvec1) | |
| 254 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 20 times.
|
25 | if(nattempt.ge.10) then |
| 255 | test_vector = site_vector + & | ||
| 256 |
2/2✓ Branch 0 taken 15 times.
✓ Branch 1 taken 5 times.
|
20 | ( rvec1 * 2._real32 - 1._real32 ) * step_size_fine / abc |
| 257 | else | ||
| 258 | test_vector = site_vector + & | ||
| 259 |
2/2✓ Branch 0 taken 60 times.
✓ Branch 1 taken 20 times.
|
80 | ( rvec1 * 2._real32 - 1._real32 ) * step_size_coarse / abc |
| 260 | end if | ||
| 261 |
3/4✓ Branch 0 taken 75 times.
✓ Branch 1 taken 25 times.
✓ Branch 2 taken 75 times.
✗ Branch 3 not taken.
|
100 | test_vector = test_vector - floor(test_vector) |
| 262 |
2/2✓ Branch 0 taken 75 times.
✓ Branch 1 taken 25 times.
|
100 | do j = 1, 3 |
| 263 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 75 times.
|
75 | if(test_vector(j).lt.bounds(1,j) .or. test_vector(j).ge.bounds(2,j)) & |
| 264 | 25 | cycle walk_loop | |
| 265 | end do | ||
| 266 | |||
| 267 | !------------------------------------------------------------------------ | ||
| 268 | ! evaluate the test point | ||
| 269 | !------------------------------------------------------------------------ | ||
| 270 | test_value = evaluate_point( distribs_container, & | ||
| 271 | test_vector, species, basis, radius_list & | ||
| 272 | 25 | ) | |
| 273 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 25 times.
|
25 | if(test_value.lt.1.E-6) cycle walk_loop |
| 274 | !------------------------------------------------------------------------ | ||
| 275 | ! if viability of test point is less than current point, then we | ||
| 276 | ! are stuck at the current point and need to try again | ||
| 277 | !------------------------------------------------------------------------ | ||
| 278 |
2/2✓ Branch 0 taken 24 times.
✓ Branch 1 taken 1 times.
|
25 | if(test_value.lt.site_value) then |
| 279 | 24 | nstuck = nstuck + 1 | |
| 280 |
2/2✓ Branch 0 taken 15 times.
✓ Branch 1 taken 9 times.
|
24 | if(nstuck.ge.10) then |
| 281 | 15 | nattempt = nattempt + 1 | |
| 282 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 15 times.
|
15 | if(crude_norm.lt.site_value) & |
| 283 | crude_norm = ( & | ||
| 284 | crude_norm + site_value/real(nattempt) & | ||
| 285 | ✗ | ) / 2._real32 | |
| 286 | |||
| 287 | ! if we have tried 10 times, and still no luck, then we need to | ||
| 288 | ! reduce the tolerance | ||
| 289 |
2/2✓ Branch 0 taken 6 times.
✓ Branch 1 taken 9 times.
|
15 | if(nattempt.ge.10) site_value = site_value / crude_norm |
| 290 | 15 | call random_number(rtmp1) | |
| 291 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 14 times.
|
15 | if(rtmp1.lt.site_value) exit walk_loop |
| 292 | end if | ||
| 293 | else | ||
| 294 | 1 | nstuck = 0 | |
| 295 | 1 | site_vector = test_vector | |
| 296 | 1 | site_value = test_value | |
| 297 | |||
| 298 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | if(nattempt.ge.10) test_value = test_value / crude_norm |
| 299 | 1 | call random_number(rtmp1) | |
| 300 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | if(rtmp1.lt.test_value) exit walk_loop |
| 301 | end if | ||
| 302 | |||
| 303 | end do walk_loop | ||
| 304 | |||
| 305 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 times.
|
4 | point = site_vector |
| 306 | 1 | viable=.true. | |
| 307 | |||
| 308 | end function place_method_walk | ||
| 309 | !############################################################################### | ||
| 310 | |||
| 311 | |||
| 312 | !############################################################################### | ||
| 313 | 1 | function place_method_growth( distribs_container, & | |
| 314 | prior_point, prior_species, & | ||
| 315 | bounds, & | ||
| 316 | basis, species, & | ||
| 317 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | radius_list, max_attempts, & |
| 318 | step_size_coarse, step_size_fine, & | ||
| 319 | viable & | ||
| 320 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
1 | ) result(point) |
| 321 | !! Random walk placement method. | ||
| 322 | !! | ||
| 323 | !! This method places the atom using a random walk method. | ||
| 324 | !! An initial point is chosen at random, and then points nearby are tested | ||
| 325 | !! to see if they are more suitable than the current point. If they are, | ||
| 326 | !! the query point is moved to the new point and the process is repeated. | ||
| 327 | !! The process is repeated, with each point being tested against a random | ||
| 328 | !! number. If the random number is less than the suitability of the point, | ||
| 329 | !! the atom is placed at that point. | ||
| 330 | implicit none | ||
| 331 | |||
| 332 | ! Arguments | ||
| 333 | type(distribs_container_type), intent(in) :: distribs_container | ||
| 334 | !! Distribution function (gvector) container. | ||
| 335 | real(real32), dimension(3), intent(in) :: prior_point | ||
| 336 | !! Point to start walk from. | ||
| 337 | integer, intent(in) :: prior_species | ||
| 338 | !! Species of last atom placed. | ||
| 339 | real(real32), dimension(2,3), intent(in) :: bounds | ||
| 340 | !! Bounds of the unit cell. | ||
| 341 | type(extended_basis_type), intent(inout) :: basis | ||
| 342 | !! Structure to add atom to. | ||
| 343 | integer, intent(in) :: species | ||
| 344 | !! Species index to add atom to. | ||
| 345 | integer, intent(in) :: max_attempts | ||
| 346 | !! Limit on number of attempts. | ||
| 347 | real(real32), intent(in) :: step_size_coarse, step_size_fine | ||
| 348 | !! Step sizes for random walk. | ||
| 349 | logical, intent(out) :: viable | ||
| 350 | !! Boolean to indicate if point is viable. | ||
| 351 | real(real32), dimension(:), intent(in) :: radius_list | ||
| 352 | !! List of radii for each pair of elements. | ||
| 353 | real(real32), dimension(3) :: point | ||
| 354 | !! Point to add atom to. | ||
| 355 | |||
| 356 | ! Local variables | ||
| 357 | integer :: i, j, idx | ||
| 358 | !! Loop indices. | ||
| 359 | integer :: nattempt, nstuck | ||
| 360 | !! Number of attempts and number of times stuck at same site | ||
| 361 | real(real32) :: rtmp1, min_radius | ||
| 362 | !! Random number and minimum radius. | ||
| 363 | real(real32), dimension(3) :: rvec1, abc | ||
| 364 | !! Random vector and lattice constants. | ||
| 365 | real(real32) :: crude_norm | ||
| 366 | !! Crude normalisation. | ||
| 367 | real(real32) :: site_value, test_value | ||
| 368 | !! Viability values. | ||
| 369 | real(real32), dimension(3) :: site_vector, test_vector | ||
| 370 | !! Vectors for gridpoints. | ||
| 371 | real(real32), dimension(3,3) :: inverse_lattice | ||
| 372 | |||
| 373 | |||
| 374 | 1 | viable = .false. | |
| 375 | |||
| 376 | !--------------------------------------------------------------------------- | ||
| 377 | ! get the lattice constants and the inverse lattice | ||
| 378 | !--------------------------------------------------------------------------- | ||
| 379 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 times.
|
4 | do i = 1, 3 |
| 380 |
5/6✓ Branch 0 taken 9 times.
✓ Branch 1 taken 3 times.
✓ Branch 2 taken 3 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
|
13 | abc(i) = norm2(basis%lat(i,:)) |
| 381 | end do | ||
| 382 | 1 | inverse_lattice = inverse_3x3(basis%lat) | |
| 383 | |||
| 384 | |||
| 385 | !--------------------------------------------------------------------------- | ||
| 386 | ! get the index of the pair of species | ||
| 387 | !--------------------------------------------------------------------------- | ||
| 388 | idx = distribs_container%get_pair_index( & | ||
| 389 | basis%spec(prior_species)%name, & | ||
| 390 | basis%spec(species)%name & | ||
| 391 | 1 | ) | |
| 392 | min_radius = & | ||
| 393 | max( & | ||
| 394 | distribs_container%cutoff_min(1), & | ||
| 395 | radius_list(idx) * distribs_container%radius_distance_tol(1) & | ||
| 396 | 1 | ) | |
| 397 | |||
| 398 | |||
| 399 | !--------------------------------------------------------------------------- | ||
| 400 | ! test a random point within a spherical shell around the prior point | ||
| 401 | !--------------------------------------------------------------------------- | ||
| 402 | 1 | i = 0 | |
| 403 | 436 | shell_loop: do | |
| 404 | 437 | i = i + 1 | |
| 405 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 437 times.
|
437 | if(i.gt.max_attempts) return |
| 406 | 437 | call random_number(rvec1) | |
| 407 | ! map rvec1(1) to ring between min_radius and min_radius + 1.0 | ||
| 408 | 437 | rvec1(1) = rvec1(1) + min_radius ! r | |
| 409 | 437 | rvec1(2) = rvec1(2) * 2._real32 * pi ! theta | |
| 410 | 437 | rvec1(3) = rvec1(3) * pi ! phi | |
| 411 | ! convert from spherical to cartesian | ||
| 412 | rvec1 = [ & | ||
| 413 | rvec1(1) * cos(rvec1(2)) * sin(rvec1(3)), & | ||
| 414 | rvec1(1) * sin(rvec1(2)) * sin(rvec1(3)), & | ||
| 415 | rvec1(1) * cos(rvec1(3)) & | ||
| 416 |
2/2✓ Branch 0 taken 1311 times.
✓ Branch 1 taken 437 times.
|
1748 | ] |
| 417 | ! convert from cartesian to direct | ||
| 418 |
2/2✓ Branch 1 taken 1311 times.
✓ Branch 2 taken 437 times.
|
1748 | rvec1 = matmul(rvec1, inverse_lattice) |
| 419 |
2/2✓ Branch 0 taken 1311 times.
✓ Branch 1 taken 437 times.
|
1748 | site_vector = prior_point + rvec1 |
| 420 |
2/2✓ Branch 0 taken 1213 times.
✓ Branch 1 taken 296 times.
|
1509 | do j = 1, 3 |
| 421 |
2/2✓ Branch 0 taken 141 times.
✓ Branch 1 taken 1072 times.
|
1213 | if(site_vector(j).lt.bounds(1,j) .or. site_vector(j).ge.bounds(2,j)) & |
| 422 | 437 | cycle shell_loop | |
| 423 | end do | ||
| 424 | ! now evaluate the point and check if it passes the initial criteria | ||
| 425 | site_value = evaluate_point( distribs_container, & | ||
| 426 | site_vector, species, basis, radius_list & | ||
| 427 | 296 | ) | |
| 428 | 296 | call random_number(rtmp1) | |
| 429 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 295 times.
|
296 | if(rtmp1.lt.site_value) exit shell_loop |
| 430 | |||
| 431 | end do shell_loop | ||
| 432 | |||
| 433 | |||
| 434 | !--------------------------------------------------------------------------- | ||
| 435 | ! now do a random walk to find a suitable point to place the atom | ||
| 436 | !--------------------------------------------------------------------------- | ||
| 437 | 1 | nattempt = 0 | |
| 438 | 1 | nstuck = 0 | |
| 439 | 1 | crude_norm = 0.5_real32 | |
| 440 | 1 | i = 0 | |
| 441 | 32 | walk_loop : do | |
| 442 | 33 | i = i + 1 | |
| 443 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 33 times.
|
33 | if(i.gt.max_attempts) return |
| 444 | !------------------------------------------------------------------------ | ||
| 445 | ! if we have tried 10 times, then we need to reduce the step size | ||
| 446 | ! get the new test point and map it back into the unit cell | ||
| 447 | !------------------------------------------------------------------------ | ||
| 448 | 33 | call random_number(rvec1) | |
| 449 |
2/2✓ Branch 0 taken 4 times.
✓ Branch 1 taken 29 times.
|
33 | if(nattempt.ge.10) then |
| 450 | test_vector = site_vector + & | ||
| 451 |
2/2✓ Branch 0 taken 12 times.
✓ Branch 1 taken 4 times.
|
16 | ( rvec1 * 2._real32 - 1._real32 ) * step_size_fine / abc |
| 452 | else | ||
| 453 | test_vector = site_vector + & | ||
| 454 |
2/2✓ Branch 0 taken 87 times.
✓ Branch 1 taken 29 times.
|
116 | ( rvec1 * 2._real32 - 1._real32 ) * step_size_coarse / abc |
| 455 | end if | ||
| 456 |
3/4✓ Branch 0 taken 99 times.
✓ Branch 1 taken 33 times.
✓ Branch 2 taken 99 times.
✗ Branch 3 not taken.
|
132 | test_vector = test_vector - floor(test_vector) |
| 457 |
2/2✓ Branch 0 taken 99 times.
✓ Branch 1 taken 30 times.
|
129 | do j = 1, 3 |
| 458 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 96 times.
|
99 | if(test_vector(j).lt.bounds(1,j) .or. test_vector(j).ge.bounds(2,j)) & |
| 459 | 33 | cycle walk_loop | |
| 460 | end do | ||
| 461 | |||
| 462 | !------------------------------------------------------------------------ | ||
| 463 | ! evaluate the test point | ||
| 464 | !------------------------------------------------------------------------ | ||
| 465 | test_value = evaluate_point( distribs_container, & | ||
| 466 | test_vector, species, basis, radius_list & | ||
| 467 | 30 | ) | |
| 468 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 30 times.
|
30 | if(test_value.lt.1.E-6) cycle walk_loop |
| 469 | !------------------------------------------------------------------------ | ||
| 470 | ! if viability of test point is less than current point, then we | ||
| 471 | ! are stuck at the current point and need to try again | ||
| 472 | !------------------------------------------------------------------------ | ||
| 473 |
2/2✓ Branch 0 taken 27 times.
✓ Branch 1 taken 3 times.
|
30 | if(test_value.lt.site_value) then |
| 474 | 27 | nstuck = nstuck + 1 | |
| 475 |
2/2✓ Branch 0 taken 14 times.
✓ Branch 1 taken 13 times.
|
27 | if(nstuck.ge.10) then |
| 476 | 14 | nattempt = nattempt + 1 | |
| 477 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
|
14 | if(crude_norm.lt.site_value) & |
| 478 | crude_norm = ( & | ||
| 479 | crude_norm + site_value/real(nattempt) & | ||
| 480 | ✗ | ) / 2._real32 | |
| 481 | |||
| 482 | ! if we have tried 10 times, and still no luck, then we need to | ||
| 483 | ! reduce the tolerance | ||
| 484 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 9 times.
|
14 | if(nattempt.ge.10) site_value = site_value / crude_norm |
| 485 | 14 | call random_number(rtmp1) | |
| 486 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 13 times.
|
14 | if(rtmp1.lt.site_value) exit walk_loop |
| 487 | end if | ||
| 488 | else | ||
| 489 | 3 | nstuck = 0 | |
| 490 | 3 | site_vector = test_vector | |
| 491 | 3 | site_value = test_value | |
| 492 | |||
| 493 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 3 times.
|
3 | if(nattempt.ge.10) test_value = test_value / crude_norm |
| 494 | 3 | call random_number(rtmp1) | |
| 495 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 3 times.
|
3 | if(rtmp1.lt.test_value) exit walk_loop |
| 496 | end if | ||
| 497 | |||
| 498 | end do walk_loop | ||
| 499 | |||
| 500 |
2/2✓ Branch 0 taken 3 times.
✓ Branch 1 taken 1 times.
|
4 | point = site_vector |
| 501 | 1 | viable=.true. | |
| 502 | |||
| 503 | end function place_method_growth | ||
| 504 | !############################################################################### | ||
| 505 | |||
| 506 | |||
| 507 | !############################################################################### | ||
| 508 | 32 | function place_method_min( & | |
| 509 | 32 | points, species, & | |
| 510 |
2/4✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 32 times.
✗ Branch 3 not taken.
|
32 | species_index_list, viable & |
| 511 |
1/2✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
|
32 | ) result(point) |
| 512 | !! Global minimum placement method. | ||
| 513 | !! | ||
| 514 | !! This method places the atom at the gridpoint with the highest | ||
| 515 | !! suitability. | ||
| 516 | implicit none | ||
| 517 | |||
| 518 | ! Arguments | ||
| 519 | logical, intent(out) :: viable | ||
| 520 | !! Boolean to indicate if point is viable. | ||
| 521 | integer, intent(in) :: species | ||
| 522 | !! Species index to add atom to. | ||
| 523 | integer, dimension(:), intent(in) :: species_index_list | ||
| 524 | !! List of species indices to add atoms to. | ||
| 525 | real(real32), dimension(:,:), intent(in) :: points | ||
| 526 | !! List of gridpoints to consider. | ||
| 527 | real(real32), dimension(3) :: point | ||
| 528 | !! Point to add atom to. | ||
| 529 | |||
| 530 | ! Local variables | ||
| 531 | integer :: species_index | ||
| 532 | !! Index of species in list. | ||
| 533 | integer :: best_gridpoint | ||
| 534 | !! Index of best gridpoint. | ||
| 535 | |||
| 536 | |||
| 537 | 32 | viable = .false. | |
| 538 | |||
| 539 | !--------------------------------------------------------------------------- | ||
| 540 | ! find the gridpoint with the highest viability | ||
| 541 | !--------------------------------------------------------------------------- | ||
| 542 |
2/4✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 32 times.
✗ Branch 3 not taken.
|
32 | species_index = findloc(species_index_list, species, 1) |
| 543 |
6/10✓ Branch 0 taken 32 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 32 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✓ Branch 6 taken 294322 times.
✓ Branch 7 taken 32 times.
✓ Branch 8 taken 512 times.
✓ Branch 9 taken 293810 times.
|
294354 | best_gridpoint = maxloc(points(4+species_index,:), dim=1) |
| 544 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 32 times.
|
32 | if(best_gridpoint.eq.0)then |
| 545 | ✗ | return | |
| 546 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 32 times.
|
32 | elseif(points(4+species,best_gridpoint).lt.1.E-6)then |
| 547 | ✗ | return | |
| 548 | end if | ||
| 549 | |||
| 550 | !--------------------------------------------------------------------------- | ||
| 551 | ! return the gridpoint with the highest viability | ||
| 552 | !--------------------------------------------------------------------------- | ||
| 553 |
2/2✓ Branch 0 taken 96 times.
✓ Branch 1 taken 32 times.
|
128 | point = points(1:3,best_gridpoint) |
| 554 | 32 | viable = .true. | |
| 555 | |||
| 556 | end function place_method_min | ||
| 557 | !############################################################################### | ||
| 558 | |||
| 559 | end module raffle__place_methods | ||
| 560 |