32 use,
intrinsic :: iso_c_binding
37 subroutine my_error(message, line, file)
39 character(len=*, kind=c_char),
intent(in) :: message
40 integer,
intent(in) :: line
41 character(len=*, kind=c_char),
intent(in) :: file
43 print *,
"* Error : '", trim(message),
"' ", line,
":", &
50 character(len=*, kind=c_char),
intent(in) :: message
51 integer,
intent(in) :: line
52 character(len=*, kind=c_char),
intent(in) :: file
54 print *,
"* Error : '", trim(message),
"' ", line,
":", &
70 use,
intrinsic :: iso_c_binding
75 real(c_double) :: cutoff
76 integer(c_int) :: number_of_particles
77 integer(c_int),
pointer :: neighborList(:,:)
88 subroutine get_neigh(data_object, number_of_cutoffs, cutoffs, &
89 neighbor_list_index, request, numnei, pnei1part, ierr) bind(c)
94 type(c_ptr),
value,
intent(in) :: data_object
95 integer(c_int),
value,
intent(in) :: number_of_cutoffs
96 real(c_double),
intent(in) :: cutoffs(number_of_cutoffs)
97 integer(c_int),
value,
intent(in) :: neighbor_list_index
98 integer(c_int),
value,
intent(in) :: request
99 integer(c_int),
intent(out) :: numnei
100 type(c_ptr),
intent(out) :: pnei1part
101 integer(c_int),
intent(out) :: ierr
104 integer(c_int),
parameter :: DIM = 3
105 integer(c_int) numberOfParticles
106 type(neighObject_type),
pointer :: neighObject
108 call c_f_pointer(data_object, neighobject)
110 if (number_of_cutoffs /= 1)
then 111 call my_warning(
"invalid number of cutoffs", __line__, __file__)
116 if (cutoffs(1) > neighobject%cutoff)
then 117 call my_warning(
"neighbor list cutoff too small for model cutoff", &
123 if (neighbor_list_index /= 1)
then 124 call my_warning(
"wrong list index", __line__, __file__)
129 numberofparticles = neighobject%number_of_particles
131 if ( (request.gt.numberofparticles) .or. (request.lt.1))
then 133 call my_warning(
"Invalid part ID in get_neigh", &
140 numnei = neighobject%neighborList(1,request)
143 pnei1part = c_loc(neighobject%neighborList(2,request))
165 use,
intrinsic :: iso_c_binding
170 logical,
intent(in) :: half
171 integer(c_int),
intent(in) :: numberOfParticles
172 real(c_double),
dimension(3,numberOfParticles), &
174 real(c_double),
intent(in) :: cutoff
175 type(neighObject_type),
intent(inout) :: neighObject
178 integer(c_int) i, j, a
181 real(c_double) cutoff2
183 neighobject%cutoff = cutoff
187 do i=1,numberofparticles
189 do j=1,numberofparticles
190 dx(:) = coords(:, j) - coords(:, i)
191 r2 = dot_product(dx, dx)
192 if (r2.le.cutoff2)
then 194 if ( (j .gt. i) .OR. ((.not. half) .AND. (i.ne.j)) )
then 196 neighobject%neighborList(a,i) = j
201 neighobject%neighborList(1,i) = a-1
226 coords, MiddlePartId)
227 use,
intrinsic :: iso_c_binding
229 integer(c_int),
parameter :: cd = c_double
232 real(c_double),
intent(in) :: FCCspacing
233 integer(c_int),
intent(in) :: nCellsPerSide
234 logical,
intent(in) :: periodic
235 real(c_double),
intent(out) :: coords(3,*)
236 integer(c_int),
intent(out) :: MiddlePartId
240 real(c_double) FCCshifts(3,4)
241 real(c_double) latVec(3)
242 integer(c_int) a, i, j, k, m
246 fccshifts(1,1) = 0.0_cd
247 fccshifts(2,1) = 0.0_cd
248 fccshifts(3,1) = 0.0_cd
249 fccshifts(1,2) = 0.5_cd*fccspacing
250 fccshifts(2,2) = 0.5_cd*fccspacing
251 fccshifts(3,2) = 0.0_cd
252 fccshifts(1,3) = 0.5_cd*fccspacing
253 fccshifts(2,3) = 0.0_cd
254 fccshifts(3,3) = 0.5_cd*fccspacing
255 fccshifts(1,4) = 0.0_cd
256 fccshifts(2,4) = 0.5_cd*fccspacing
257 fccshifts(3,4) = 0.5_cd*fccspacing
262 latvec(1) = (i-1)*fccspacing
264 latvec(2) = (j-1)*fccspacing
266 latvec(3) = (k-1)*fccspacing
269 coords(:,a) = latvec + fccshifts(:,m)
270 if ((i.eq.ncellsperside/2+1).and.(j.eq.ncellsperside/2+1) .and. &
271 (k.eq.ncellsperside/2+1) .and. (m.eq.1))
then 272 coords(:,1) = latvec + fccshifts(:,m)
277 if (.not. periodic)
then 280 latvec(1) = ncellsperside*fccspacing
281 latvec(2) = (i-1)*fccspacing
282 latvec(3) = (j-1)*fccspacing
283 a = a+1; coords(:,a) = latvec
284 a = a+1; coords(:,a) = latvec + fccshifts(:,4)
286 latvec(1) = (i-1)*fccspacing
287 latvec(2) = ncellsperside*fccspacing
288 latvec(3) = (j-1)*fccspacing
289 a = a+1; coords(:,a) = latvec
290 a = a+1; coords(:,a) = latvec + fccshifts(:,3)
292 latvec(1) = (i-1)*fccspacing
293 latvec(2) = (j-1)*fccspacing
294 latvec(3) = ncellsperside*fccspacing
295 a = a+1; coords(:,a) = latvec
296 a = a+1; coords(:,a) = latvec + fccshifts(:,2)
299 if (.not. periodic)
then 301 latvec(1) = (i-1)*fccspacing
302 latvec(2) = ncellsperside*fccspacing
303 latvec(3) = ncellsperside*fccspacing
304 a = a+1; coords(:,a) = latvec
305 latvec(1) = ncellsperside*fccspacing
306 latvec(2) = (i-1)*fccspacing
307 latvec(3) = ncellsperside*fccspacing
308 a = a+1; coords(:,a) = latvec
309 latvec(1) = ncellsperside*fccspacing
310 latvec(2) = ncellsperside*fccspacing
311 latvec(3) = (i-1)*fccspacing
312 a = a+1; coords(:,a) = latvec
315 if (.not. periodic)
then 317 a = a+1; coords(:,a) = ncellsperside*fccspacing
340 use,
intrinsic :: iso_c_binding
346 integer(c_int),
parameter :: cd = c_double
348 integer(c_int),
parameter :: ncellsperside = 2
349 integer(c_int),
parameter :: dim = 3
351 real(c_double),
parameter :: cutpad = 0.75_cd
352 real(c_double),
parameter :: fccspacing = 5.260_cd
353 real(c_double),
parameter :: min_spacing = 0.8*fccspacing
354 real(c_double),
parameter :: max_spacing = 1.2*fccspacing
355 real(c_double),
parameter :: spacing_incr = 0.025*fccspacing
356 real(c_double) :: current_spacing
357 real(c_double) :: force_norm
359 character(len=256, kind=c_char) :: modelname
361 integer(c_int),
parameter :: &
362 n = 4*(ncellsperside)**3 + 6*(ncellsperside)**2 + 3*(ncellsperside) + 1
366 type(kim_model_handle_type) :: model_handle
367 type(kim_compute_arguments_handle_type) :: compute_arguments_handle
368 real(c_double) :: influence_distance
369 integer(c_int) :: number_of_cutoffs
370 real(c_double) :: cutoff
371 real(c_double) :: cutoffs(1)
372 integer(c_int),
target :: particle_species_codes(n)
373 integer(c_int),
target :: particle_contributing(n)
374 real(c_double),
target :: energy
375 real(c_double),
target :: coords(dim, n)
376 real(c_double),
target :: forces(dim, n)
377 integer(c_int) i, j, ierr, ierr2
379 integer(c_int) species_is_supported
380 integer(c_int) species_code
381 integer(c_int) requested_units_accepted
389 print
'("Please enter a valid KIM model name: ")' 395 print *,
'This is Test : ex_test_Ar_fcc_cluster_fortran' 398 print
'("Results for KIM Model : ",A)', trim(modelname)
402 call kim_model_create(kim_numbering_one_based, &
404 kim_energy_unit_ev, &
406 kim_temperature_unit_k, &
409 requested_units_accepted, &
412 call my_error(
"kim_api_create", __line__, __file__)
416 if (requested_units_accepted == 0)
then 417 call my_error(
"Must adapt to model units", __line__, __file__)
422 call kim_model_get_species_support_and_code(model_handle, &
423 kim_species_name_ar, species_is_supported, species_code, ierr)
424 if ((ierr /= 0) .or. (species_is_supported /= 1))
then 425 call my_error(
"Model does not support Ar", __line__, __file__)
432 call kim_model_compute_arguments_create( &
433 model_handle, compute_arguments_handle, ierr)
435 call my_error(
"kim_model_compute_arguments_create", __line__, __file__)
440 call kim_compute_arguments_set_argument_pointer(compute_arguments_handle, &
441 kim_compute_argument_name_number_of_particles, n, ierr2)
443 call kim_compute_arguments_set_argument_pointer(compute_arguments_handle, &
444 kim_compute_argument_name_particle_species_codes, particle_species_codes, &
447 call kim_compute_arguments_set_argument_pointer(compute_arguments_handle, &
448 kim_compute_argument_name_particle_contributing, particle_contributing, &
451 call kim_compute_arguments_set_argument_pointer(compute_arguments_handle, &
452 kim_compute_argument_name_coordinates, coords, ierr2)
454 call kim_compute_arguments_set_argument_pointer(compute_arguments_handle, &
455 kim_compute_argument_name_partial_energy, energy, ierr2)
457 call kim_compute_arguments_set_argument_pointer(compute_arguments_handle, &
458 kim_compute_argument_name_partial_forces, forces, ierr2)
461 call my_error(
"set_argument_pointer", __line__, __file__)
466 allocate(neighobject%neighborList(n+1,n))
467 neighobject%number_of_particles = n
471 call kim_compute_arguments_set_callback_pointer(compute_arguments_handle, &
472 kim_compute_callback_name_get_neighbor_list, kim_language_name_fortran, &
473 c_funloc(
get_neigh), c_loc(neighobject), ierr)
475 call my_error(
"set_callback_pointer", __line__, __file__)
478 call kim_model_get_influence_distance(model_handle, influence_distance)
479 call kim_model_get_number_of_neighbor_list_cutoffs(model_handle, &
481 if (number_of_cutoffs /= 1)
then 482 call my_error(
"too many cutoffs", __line__, __file__)
484 call kim_model_get_neighbor_list_cutoffs(model_handle, cutoffs, ierr)
486 call my_error(
"get_cutoffs", __line__, __file__)
493 particle_species_codes(i) = species_code
498 particle_contributing(i) = 1
502 print
'(3A20)',
"Energy",
"Force Norm",
"Lattice Spacing" 504 current_spacing = min_spacing
505 do while (current_spacing < max_spacing)
514 call kim_model_compute(model_handle, compute_arguments_handle, ierr)
516 call my_error(
"kim_api_model_compute", __line__, __file__)
523 force_norm = force_norm + forces(j,i)*forces(j,i)
526 force_norm = sqrt(force_norm)
530 print
'(3ES20.10)', energy, force_norm, current_spacing
532 current_spacing = current_spacing + spacing_incr
536 deallocate( neighobject%neighborList )
538 call kim_model_compute_arguments_destroy(&
539 model_handle, compute_arguments_handle, ierr)
541 call my_error(
"compute_arguments_destroy", __line__, __file__)
543 call kim_model_destroy(model_handle)
subroutine my_error(message, line, file)
subroutine neigh_pure_cluster_neighborlist(half, numberOfParticles, coords, cutoff, neighObject)
subroutine create_fcc_configuration(FCCspacing, nCellsPerSide, periodic, coords, MiddlePartId)
program ex_test_ar_fcc_cluster_fortran
subroutine my_warning(message, line, file)
integer(c_int) function, public get_neigh(pkim, mode, request, part, numnei, pnei1part, pRij)