29 use,
intrinsic :: iso_c_binding
35 recursive subroutine my_error(message)
37 character(len=*, kind=c_char),
intent(in) :: message
39 print *,
"* Error : ", trim(message)
45 character(len=*, kind=c_char),
intent(in) :: message
47 print *,
"* Warning : ", trim(message)
61 use,
intrinsic :: iso_c_binding
66 real(c_double) :: cutoff
67 integer(c_int) :: number_of_particles
68 type(c_ptr) :: neighbor_list_pointer
79 recursive subroutine get_neigh(data_object, number_of_neighbor_lists, &
80 cutoffs, neighbor_list_index, request, &
81 numnei, pnei1part, ierr) bind(c)
86 type(c_ptr),
value,
intent(in) :: data_object
87 integer(c_int),
value,
intent(in) :: number_of_neighbor_lists
88 real(c_double),
intent(in) :: cutoffs(*)
89 integer(c_int),
value,
intent(in) :: neighbor_list_index
90 integer(c_int),
value,
intent(in) :: request
91 integer(c_int),
intent(out) :: numnei
92 type(c_ptr),
intent(out) :: pnei1part
93 integer(c_int),
intent(out) :: ierr
96 integer(c_int) numberOfParticles
98 integer(c_int),
pointer :: neighborList(:, :)
100 call c_f_pointer(data_object, neighobject)
101 numberofparticles = neighobject%number_of_particles
102 call c_f_pointer(neighobject%neighbor_list_pointer, neighborlist, &
103 [numberofparticles + 1, numberofparticles])
105 if (number_of_neighbor_lists /= 1)
then 111 if (cutoffs(1) > neighobject%cutoff)
then 112 call my_warning(
"neighbor list cutoff too small for model cutoff")
117 if (neighbor_list_index /= 1)
then 123 if ((request > numberofparticles) .or. (request < 1))
then 125 call my_warning(
"Invalid part ID in get_neigh")
131 numnei = neighborlist(1, request)
134 pnei1part = c_loc(neighborlist(2, request))
154 half, numberOfParticles, coords, cutoff, neighObject)
155 use,
intrinsic :: iso_c_binding
160 logical,
intent(in) :: half
161 integer(c_int),
intent(in) :: numberOfParticles
162 real(c_double),
dimension(3, numberOfParticles), &
164 real(c_double),
intent(in) :: cutoff
168 integer(c_int) i, j, a
171 real(c_double) cutoff2
172 integer(c_int),
pointer :: neighborList(:, :)
174 call c_f_pointer(neighobject%neighbor_list_pointer, neighborlist, &
175 [numberofparticles + 1, numberofparticles])
177 neighobject%cutoff = cutoff
181 do i = 1, numberofparticles
183 do j = 1, numberofparticles
184 dx(:) = coords(:, j) - coords(:, i)
185 r2 = dot_product(dx, dx)
186 if (r2 <= cutoff2)
then 188 if ((j > i) .OR. ((.not. half) .AND. (i /= j)))
then 190 neighborlist(a, i) = j
195 neighborlist(1, i) = a - 1
220 periodic, coords, MiddlePartId)
221 use,
intrinsic :: iso_c_binding
223 integer(c_int),
parameter :: cd = c_double
226 real(c_double),
intent(in) :: FCCspacing
227 integer(c_int),
intent(in) :: nCellsPerSide
228 logical,
intent(in) :: periodic
229 real(c_double),
intent(out) :: coords(3, *)
230 integer(c_int),
intent(out) :: MiddlePartId
234 real(c_double) FCCshifts(3, 4)
235 real(c_double) latVec(3)
236 integer(c_int) a, i, j, k, m
240 fccshifts(1, 1) = 0.0_cd
241 fccshifts(2, 1) = 0.0_cd
242 fccshifts(3, 1) = 0.0_cd
243 fccshifts(1, 2) = 0.5_cd * fccspacing
244 fccshifts(2, 2) = 0.5_cd * fccspacing
245 fccshifts(3, 2) = 0.0_cd
246 fccshifts(1, 3) = 0.5_cd * fccspacing
247 fccshifts(2, 3) = 0.0_cd
248 fccshifts(3, 3) = 0.5_cd * fccspacing
249 fccshifts(1, 4) = 0.0_cd
250 fccshifts(2, 4) = 0.5_cd * fccspacing
251 fccshifts(3, 4) = 0.5_cd * fccspacing
255 do i = 1, ncellsperside
256 latvec(1) = (i - 1) * fccspacing
257 do j = 1, ncellsperside
258 latvec(2) = (j - 1) * fccspacing
259 do k = 1, ncellsperside
260 latvec(3) = (k - 1) * fccspacing
263 coords(:, a) = latvec + fccshifts(:, m)
264 if ((i == ncellsperside / 2 + 1) &
265 .and. (j == ncellsperside / 2 + 1) &
266 .and. (k == ncellsperside / 2 + 1) &
270 coords(:, 1) = latvec + fccshifts(:, m)
275 if (.not. periodic)
then 278 latvec(1) = ncellsperside * fccspacing
279 latvec(2) = (i - 1) * fccspacing
280 latvec(3) = (j - 1) * fccspacing
281 a = a + 1; coords(:, a) = latvec
282 a = a + 1; coords(:, a) = latvec + fccshifts(:, 4)
284 latvec(1) = (i - 1) * fccspacing
285 latvec(2) = ncellsperside * fccspacing
286 latvec(3) = (j - 1) * fccspacing
287 a = a + 1; coords(:, a) = latvec
288 a = a + 1; coords(:, a) = latvec + fccshifts(:, 3)
290 latvec(1) = (i - 1) * fccspacing
291 latvec(2) = (j - 1) * fccspacing
292 latvec(3) = ncellsperside * fccspacing
293 a = a + 1; coords(:, a) = latvec
294 a = a + 1; coords(:, a) = latvec + fccshifts(:, 2)
297 if (.not. periodic)
then 299 latvec(1) = (i - 1) * fccspacing
300 latvec(2) = ncellsperside * fccspacing
301 latvec(3) = ncellsperside * fccspacing
302 a = a + 1; coords(:, a) = latvec
303 latvec(1) = ncellsperside * fccspacing
304 latvec(2) = (i - 1) * fccspacing
305 latvec(3) = ncellsperside * fccspacing
306 a = a + 1; coords(:, a) = latvec
307 latvec(1) = ncellsperside * fccspacing
308 latvec(2) = ncellsperside * fccspacing
309 latvec(3) = (i - 1) * fccspacing
310 a = a + 1; coords(:, a) = latvec
313 if (.not. periodic)
then 315 a = a + 1; coords(:, a) = ncellsperside * fccspacing
338 use,
intrinsic :: iso_c_binding
345 integer(c_int),
parameter :: cd = c_double
347 integer(c_int),
parameter :: nCellsPerSide = 2
348 integer(c_int),
parameter :: DIM = 3
350 real(c_double),
parameter :: cutpad = 0.75_cd
351 real(c_double),
parameter :: FCCspacing = 5.260_cd
352 real(c_double),
parameter :: min_spacing = 0.8 * fccspacing
353 real(c_double),
parameter :: max_spacing = 1.2 * fccspacing
354 real(c_double),
parameter :: spacing_incr = 0.025 * fccspacing
355 real(c_double) :: current_spacing
356 real(c_double) :: force_norm
358 character(len=256, kind=c_char) :: modelname
360 integer(c_int),
parameter :: N = 4 * (ncellsperside)**3 &
361 + 6 * (ncellsperside)**2 &
362 + 3 * (ncellsperside) + 1
365 integer(c_int),
allocatable,
target :: neighborList(:, :)
367 type(kim_model_handle_type) :: model_handle
368 type(kim_compute_arguments_handle_type) :: compute_arguments_handle
369 real(c_double) :: influence_distance
370 integer(c_int) :: number_of_neighbor_lists
371 real(c_double) :: cutoff
372 real(c_double) :: cutoffs(1)
374 model_will_not_request_neighbors_of_noncontributing_particles(1)
375 integer(c_int),
target :: particle_species_codes(n)
376 integer(c_int),
target :: particle_contributing(n)
377 real(c_double),
target :: energy
378 real(c_double),
target :: coords(dim, n)
379 real(c_double),
target :: forces(dim, n)
380 integer(c_int) i, j, ierr, ierr2
382 integer(c_int) species_is_supported
383 integer(c_int) species_code
384 integer(c_int) requested_units_accepted
385 integer(c_int) number_of_model_routine_names
386 type(kim_model_routine_name_type) model_routine_name
387 integer(c_int) present
388 integer(c_int) required
389 type(kim_supported_extensions_type),
target :: supported_extensions
390 character(len=KIM_MAX_EXTENSION_ID_LENGTH, kind=c_char) :: id_string
396 supported_extensions%number_of_supported_extensions = 0
399 print
'("Please enter a valid KIM model name: ")' 400 read (*, *) modelname
404 call kim_model_create(kim_numbering_one_based, &
406 kim_energy_unit_ev, &
408 kim_temperature_unit_k, &
411 requested_units_accepted, &
418 if (requested_units_accepted == 0)
then 419 call my_error(
"Must adapt to model units")
423 call kim_get_number_of_model_routine_names(number_of_model_routine_names)
424 do i = 1, number_of_model_routine_names
425 call kim_get_model_routine_name(i, model_routine_name, ierr)
427 call my_error(
"kim_get_model_routine_name")
429 call kim_is_routine_present(model_handle, model_routine_name,
present, &
432 call my_error(
"kim_is_routine_present")
435 if ((
present == 1) .and. (required == 1))
then 436 if (.not. ((model_routine_name == kim_model_routine_name_create) &
437 .or. (model_routine_name == &
438 kim_model_routine_name_compute_arguments_create) &
439 .or. (model_routine_name == kim_model_routine_name_compute) &
440 .or. (model_routine_name == kim_model_routine_name_refresh) &
441 .or. (model_routine_name == &
442 kim_model_routine_name_compute_arguments_destroy) &
443 .or. (model_routine_name == kim_model_routine_name_destroy))) &
445 call my_error(
"Unknown required ModelRoutineName found.")
452 call kim_get_species_support_and_code( &
453 model_handle, kim_species_name_ar, species_is_supported, species_code, ierr)
454 if ((ierr /= 0) .or. (species_is_supported /= 1))
then 455 call my_error(
"Model does not support Ar")
459 call kim_is_routine_present( &
460 model_handle, kim_model_routine_name_extension,
present, required, ierr)
462 call my_error(
"Unable to get Extension present/required.")
464 if (
present /= 0)
then 466 c_loc(supported_extensions), ierr)
468 call my_error(
"Error returned from kim_model_extension().")
470 write (*,
'(A,I2,A)')
"Model Supports ", &
471 supported_extensions%number_of_supported_extensions, &
473 do i = 1, supported_extensions%number_of_supported_extensions
474 call kim_c_char_array_to_string( &
475 supported_extensions%supported_extension_id(:, i), id_string)
476 write (*,
'(A,I2,A,A,A,A,I2)')
" supportedExtensionID[", i,
'] = "', &
477 trim(id_string),
'" ', &
478 "which has required = ", &
479 supported_extensions%supported_extension_required(i)
487 call kim_compute_arguments_create( &
488 model_handle, compute_arguments_handle, ierr)
490 call my_error(
"kim_model_compute_arguments_create")
495 call kim_set_argument_pointer( &
496 compute_arguments_handle, kim_compute_argument_name_number_of_particles, &
499 call kim_set_argument_pointer( &
500 compute_arguments_handle, &
501 kim_compute_argument_name_particle_species_codes, &
502 particle_species_codes, ierr2)
504 call kim_set_argument_pointer( &
505 compute_arguments_handle, &
506 kim_compute_argument_name_particle_contributing, particle_contributing, &
509 call kim_set_argument_pointer( &
510 compute_arguments_handle, kim_compute_argument_name_coordinates, coords, &
513 call kim_set_argument_pointer( &
514 compute_arguments_handle, kim_compute_argument_name_partial_energy, &
517 call kim_set_argument_pointer( &
518 compute_arguments_handle, kim_compute_argument_name_partial_forces, &
522 call my_error(
"set_argument_pointer")
527 allocate (neighborlist(n + 1, n))
528 neighobject%neighbor_list_pointer = c_loc(neighborlist)
529 neighobject%number_of_particles = n
533 call kim_set_callback_pointer( &
534 compute_arguments_handle, kim_compute_callback_name_get_neighbor_list, &
535 kim_language_name_fortran, c_funloc(
get_neigh), c_loc(neighobject), ierr)
537 call my_error(
"set_callback_pointer")
540 call kim_get_influence_distance(model_handle, influence_distance)
541 call kim_get_number_of_neighbor_lists(model_handle, &
542 number_of_neighbor_lists)
543 if (number_of_neighbor_lists /= 1)
then 544 call my_error(
"too many neighbor lists")
546 call kim_get_neighbor_list_values( &
547 model_handle, cutoffs, &
548 model_will_not_request_neighbors_of_noncontributing_particles, ierr)
550 call my_error(
"get_neighbor_list_values")
557 particle_species_codes(i) = species_code
562 particle_contributing(i) = 1
568 print *,
'This is Test : ex_test_Ar_fcc_cluster_fortran' 571 print
'("Results for KIM Model : ",A)', trim(modelname)
574 print
'(3A20)',
"Energy",
"Force Norm",
"Lattice Spacing" 576 current_spacing = min_spacing
577 do while (current_spacing < max_spacing)
583 (cutoff + cutpad), neighobject)
586 call kim_compute(model_handle, compute_arguments_handle, ierr)
588 call my_error(
"kim_api_model_compute")
595 force_norm = force_norm + forces(j, i) * forces(j, i)
598 force_norm = sqrt(force_norm)
602 print
'(3ES20.10)', energy, force_norm, current_spacing
604 current_spacing = current_spacing + spacing_incr
608 deallocate (neighborlist)
610 call kim_compute_arguments_destroy( &
611 model_handle, compute_arguments_handle, ierr)
613 call my_error(
"compute_arguments_destroy")
615 call kim_model_destroy(model_handle)
recursive subroutine, public get_neigh(data_object, number_of_neighbor_lists, cutoffs, neighbor_list_index, request, numnei, pnei1part, ierr)
recursive subroutine my_warning(message)
recursive subroutine neigh_pure_cluster_neighborlist(half, numberOfParticles, coords, cutoff, neighObject)
recursive subroutine my_error(message)
recursive subroutine create_fcc_configuration(FCCspacing, nCellsPerSide, periodic, coords, MiddlePartId)
character(len= *, kind=c_char), parameter, public kim_supported_extensions_id
program ex_test_ar_fcc_cluster_fortran
The only standard extension defined by the KIM API.