32 #include "KIM_API_status.h" 33 #define THIS_FILE_NAME __FILE__ 34 #define TRUEFALSE(TRUTH) merge(1,0,(TRUTH)) 46 use,
intrinsic :: iso_c_binding
52 integer(c_int),
pointer :: neighborlist(:,:)
53 real(c_double),
pointer :: rijlist(:,:,:)
64 integer(c_int) function get_neigh(pkim,mode,request,part,numnei,pnei1part, &
69 type(c_ptr),
intent(in) :: pkim
70 integer(c_int),
intent(in) :: mode
71 integer(c_int),
intent(in) :: request
72 integer(c_int),
intent(out) :: part
73 integer(c_int),
intent(out) :: numnei
74 type(c_ptr),
intent(out) :: pnei1part
75 type(c_ptr),
intent(out) :: prij
78 integer(c_int),
parameter :: dim = 3
79 integer(c_int),
save :: iterval = 0
81 integer(c_int) parttoreturn
82 integer(c_int),
pointer :: numberofparticles;
type(c_ptr) :: pnparts
84 integer(c_int) ier, idum
87 pnparts = kim_api_get_data(pkim,
"numberOfParticles", ier)
88 if (ier.lt.kim_status_ok)
then 89 idum = kim_api_report_error(__line__, this_file_name, &
90 "kim_api_get_data", ier)
93 call c_f_pointer(pnparts, numberofparticles)
96 pneighobject = kim_api_get_data(pkim,
"neighObject", ier)
97 if (ier.lt.kim_status_ok)
then 98 idum = kim_api_report_error(__line__, this_file_name, &
99 "kim_api_get_data", ier)
102 call c_f_pointer(pneighobject, neighobject)
104 n =
size(neighobject%neighborList, 2)
108 if (request.eq.0)
then 110 get_neigh = kim_status_neigh_iter_init_ok
112 elseif (request.eq.1)
then 113 iterval = iterval + 1
114 if (iterval.gt.n)
then 115 get_neigh = kim_status_neigh_iter_past_end
118 parttoreturn = iterval
121 idum = kim_api_report_error(__line__, this_file_name, &
122 "Invalid request in get_neigh", &
123 kim_status_neigh_invalid_request)
124 get_neigh = kim_status_neigh_invalid_request
127 elseif (mode.eq.1)
then 128 if ( (request.gt.n) .or. (request.lt.1))
then 129 idum = kim_api_report_error(__line__, this_file_name, &
130 "Invalid part ID in get_neigh", &
131 kim_status_particle_invalid_id)
132 get_neigh = kim_status_particle_invalid_id
135 parttoreturn = request
138 idum = kim_api_report_error(__line__, this_file_name, &
139 "Invalid mode in get_neigh", &
140 kim_status_neigh_invalid_mode)
141 get_neigh = kim_status_neigh_invalid_mode
149 numnei = neighobject%neighborList(1,part)
152 pnei1part = c_loc(neighobject%neighborList(2,part))
155 if (
associated(neighobject%RijList))
then 156 prij = c_loc(neighobject%RijList(1,1,part))
181 boxSideLengths, neighObject)
182 use,
intrinsic :: iso_c_binding
186 integer(c_int),
parameter :: cd = c_double
189 logical,
intent(in) :: half
190 integer(c_int),
intent(in) :: numberOfParticles
191 real(c_double),
dimension(3,numberOfParticles), &
193 real(c_double),
intent(in) :: rcut
194 real(c_double),
dimension(3),
intent(in) :: boxSideLengths
195 type(neighObject_type),
intent(inout) :: neighObject
198 integer(c_int) i, j, a
205 do i=1,numberofparticles
207 do j=1,numberofparticles
208 dx(:) = coords(:, j) - coords(:, i)
209 where (abs(dx) > 0.5_cd*boxsidelengths)
210 dx = dx - sign(boxsidelengths,dx)
212 r2 = dot_product(dx, dx)
213 if (r2.le.rcut2)
then 215 if ( (j .gt. i) .or. ((.not. half) .AND. (i.ne.j)) )
then 218 neighobject%neighborList(a,i) = j
224 neighobject%neighborList(1,i) = a-1
238 use,
intrinsic :: iso_c_binding
244 logical,
intent(in) :: half
245 integer(c_int),
intent(in) :: numberOfParticles
246 real(c_double),
dimension(3,numberOfParticles), &
248 real(c_double),
intent(in) :: cutoff
249 type(neighObject_type),
intent(inout) :: neighObject
252 integer(c_int) i, j, a
255 real(c_double) cutoff2
259 do i=1,numberofparticles
261 do j=1,numberofparticles
262 dx(:) = coords(:, j) - coords(:, i)
263 r2 = dot_product(dx, dx)
264 if (r2.le.cutoff2)
then 266 if ( (j .gt. i) .OR. ((.not. half) .AND. (i.ne.j)) )
then 268 neighobject%neighborList(a,i) = j
273 neighobject%neighborList(1,i) = a-1
287 use,
intrinsic :: iso_c_binding
293 logical,
intent(in) :: half
294 integer(c_int),
intent(in) :: numberOfParticles
295 real(c_double),
dimension(3,numberOfParticles), &
297 real(c_double),
intent(in) :: cutoff
298 type(neighObject_type),
intent(inout) :: neighObject
301 integer(c_int) i, j, a
304 real(c_double) cutoff2
308 do i=1,numberofparticles
310 do j=1,numberofparticles
311 dx(:) = coords(:, j) - coords(:, i)
312 r2 = dot_product(dx, dx)
313 if (r2.le.cutoff2)
then 314 if ((half .and. i.lt.j) .or. (.not.half .and. i.ne.j))
then 317 neighobject%neighborList(a,i) = j
318 neighobject%RijList(:,a-1,i) = dx
323 neighobject%neighborList(1,i) = a-1
348 coords, MiddlePartId)
349 use,
intrinsic :: iso_c_binding
351 integer(c_int),
parameter :: cd = c_double
354 real(c_double),
intent(in) :: FCCspacing
355 integer(c_int),
intent(in) :: nCellsPerSide
356 logical,
intent(in) :: periodic
357 real(c_double),
intent(out) :: coords(3,*)
358 integer(c_int),
intent(out) :: MiddlePartId
362 real(c_double) FCCshifts(3,4)
363 real(c_double) latVec(3)
364 integer(c_int) a, i, j, k, m
368 fccshifts(1,1) = 0.0_cd
369 fccshifts(2,1) = 0.0_cd
370 fccshifts(3,1) = 0.0_cd
371 fccshifts(1,2) = 0.5_cd*fccspacing
372 fccshifts(2,2) = 0.5_cd*fccspacing
373 fccshifts(3,2) = 0.0_cd
374 fccshifts(1,3) = 0.5_cd*fccspacing
375 fccshifts(2,3) = 0.0_cd
376 fccshifts(3,3) = 0.5_cd*fccspacing
377 fccshifts(1,4) = 0.0_cd
378 fccshifts(2,4) = 0.5_cd*fccspacing
379 fccshifts(3,4) = 0.5_cd*fccspacing
384 latvec(1) = (i-1)*fccspacing
386 latvec(2) = (j-1)*fccspacing
388 latvec(3) = (k-1)*fccspacing
391 coords(:,a) = latvec + fccshifts(:,m)
392 if ((i.eq.ncellsperside/2+1).and.(j.eq.ncellsperside/2+1) .and. &
393 (k.eq.ncellsperside/2+1) .and. (m.eq.1))
then 394 coords(:,1) = latvec + fccshifts(:,m)
399 if (.not. periodic)
then 402 latvec(1) = ncellsperside*fccspacing
403 latvec(2) = (i-1)*fccspacing
404 latvec(3) = (j-1)*fccspacing
405 a = a+1; coords(:,a) = latvec
406 a = a+1; coords(:,a) = latvec + fccshifts(:,4)
408 latvec(1) = (i-1)*fccspacing
409 latvec(2) = ncellsperside*fccspacing
410 latvec(3) = (j-1)*fccspacing
411 a = a+1; coords(:,a) = latvec
412 a = a+1; coords(:,a) = latvec + fccshifts(:,3)
414 latvec(1) = (i-1)*fccspacing
415 latvec(2) = (j-1)*fccspacing
416 latvec(3) = ncellsperside*fccspacing
417 a = a+1; coords(:,a) = latvec
418 a = a+1; coords(:,a) = latvec + fccshifts(:,2)
421 if (.not. periodic)
then 423 latvec(1) = (i-1)*fccspacing
424 latvec(2) = ncellsperside*fccspacing
425 latvec(3) = ncellsperside*fccspacing
426 a = a+1; coords(:,a) = latvec
427 latvec(1) = ncellsperside*fccspacing
428 latvec(2) = (i-1)*fccspacing
429 latvec(3) = ncellsperside*fccspacing
430 a = a+1; coords(:,a) = latvec
431 latvec(1) = ncellsperside*fccspacing
432 latvec(2) = ncellsperside*fccspacing
433 latvec(3) = (i-1)*fccspacing
434 a = a+1; coords(:,a) = latvec
437 if (.not. periodic)
then 439 a = a+1; coords(:,a) = ncellsperside*fccspacing
462 use,
intrinsic :: iso_c_binding
467 integer(c_int),
parameter :: cd = c_double
469 integer(c_int),
parameter :: ncellsperside = 2
470 integer(c_int),
parameter :: dim = 3
471 integer(c_int),
parameter :: aspecies = 1
473 real(c_double),
parameter :: cutpad = 0.75_cd
474 real(c_double),
parameter :: fccspacing = 5.260_cd
475 real(c_double),
parameter :: min_spacing = 0.8*fccspacing
476 real(c_double),
parameter :: max_spacing = 1.2*fccspacing
477 real(c_double),
parameter :: spacing_incr = 0.025*fccspacing
478 real(c_double) :: current_spacing
479 real(c_double) :: force_norm
481 integer(c_int),
parameter :: &
482 n = 4*(ncellsperside)**3 + 6*(ncellsperside)**2 + 3*(ncellsperside) + 1
483 integer(c_int),
parameter :: sizeone = 1
490 character(len=KIM_KEY_STRING_LENGTH) :: testkimfile =
"descriptor.kim" 491 character(len=KIM_KEY_STRING_LENGTH) :: modelname
492 character(len=KIM_KEY_STRING_LENGTH) :: nbc_method
498 integer(c_int),
pointer :: numberofparticles;
type(c_ptr) :: pnparts
499 integer(c_int),
pointer :: numcontrib;
type(c_ptr) :: pnumcontrib
500 integer(c_int),
pointer :: numberofspecies;
type(c_ptr) :: pnofspecies
501 integer(c_int),
pointer :: particlespecies(:);
type(c_ptr) :: pparticlespecies
502 real(c_double),
pointer :: cutoff;
type(c_ptr) :: pcutoff
503 real(c_double),
pointer :: energy;
type(c_ptr) :: penergy
504 real(c_double),
pointer :: coords(:,:);
type(c_ptr) :: pcoor
505 real(c_double),
pointer :: forces(:,:);
type(c_ptr) :: pforces
506 real(c_double),
pointer :: boxsidelengths(:);
type(c_ptr) :: pboxsidelengths
509 integer(c_int) ier, idum
517 print
'("Please enter a valid KIM model name: ")' 523 print *,
'This is Test : ex_test_Ar_fcc_cluster_fortran' 526 print
'("Results for KIM Model : ",A)', trim(modelname)
529 ier = kim_api_file_init(pkim, testkimfile, modelname)
530 if (ier.lt.kim_status_ok)
then 531 idum = kim_api_report_error(__line__, this_file_name, &
532 "kim_api_file_init", ier)
537 ier = kim_api_get_nbc_method(pkim, nbc_method)
538 if (ier.lt.kim_status_ok)
then 539 idum = kim_api_report_error(__line__, this_file_name, &
540 "kim_api_get_nbc_method", ier)
543 if (index(nbc_method,
"NEIGH_RVEC_H").eq.1)
then 545 elseif (index(nbc_method,
"NEIGH_PURE_H").eq.1)
then 547 elseif (index(nbc_method,
"NEIGH_RVEC_F").eq.1)
then 549 elseif (index(nbc_method,
"NEIGH_PURE_F").eq.1)
then 551 elseif (index(nbc_method,
"MI_OPBC_H").eq.1)
then 553 elseif (index(nbc_method,
"MI_OPBC_F").eq.1)
then 555 elseif (index(nbc_method,
"CLUSTER").eq.1)
then 558 ier = kim_status_fail
559 idum = kim_api_report_error(__line__, this_file_name, &
560 "Unknown NBC method", ier)
565 call kim_api_allocate(pkim, n, aspecies, ier)
566 if (ier.lt.kim_status_ok)
then 567 idum = kim_api_report_error(__line__, this_file_name, &
568 "kim_api_allocate", ier)
573 if (nbc.lt.6)
allocate(neighobject%neighborList(n+1, n))
574 if (nbc.eq.0 .or. nbc.eq.2)
then 575 allocate(neighobject%RijList(dim,n+1, n))
579 ier = kim_api_set_data(pkim,
"neighObject", sizeone, c_loc(neighobject))
580 if (ier.lt.kim_status_ok)
then 581 idum = kim_api_report_error(__line__, this_file_name, &
582 "kim_api_set_data", ier)
588 ier = kim_api_set_method(pkim,
"get_neigh", sizeone, c_funloc(
get_neigh))
589 if (ier.lt.kim_status_ok)
then 590 idum = kim_api_report_error(__line__, this_file_name, &
591 "kim_api_set_method", ier)
597 ier = kim_api_model_init(pkim)
598 if (ier.lt.kim_status_ok)
then 599 idum = kim_api_report_error(__line__, this_file_name, &
600 "kim_api_model_init", ier)
607 call kim_api_getm_data(pkim, ier, &
608 "numberOfParticles", pnparts, 1, &
609 "numberContributingParticles", pnumcontrib, truefalse((nbc.eq.0).or.(nbc.eq.1).or.(nbc.eq.4)), &
610 "numberOfSpecies", pnofspecies, 1, &
611 "particleSpecies", pparticlespecies, 1, &
612 "coordinates", pcoor, 1, &
613 "cutoff", pcutoff, 1, &
614 "boxSideLengths", pboxsidelengths, truefalse((nbc.eq.4).or.(nbc.eq.5)), &
615 "energy", penergy, 1, &
616 "forces", pforces, 1)
617 if (ier.lt.kim_status_ok)
then 618 idum = kim_api_report_error(__line__, this_file_name, &
619 "kim_api_getm_data", ier)
622 call c_f_pointer(pnparts, numberofparticles)
623 if ((nbc.eq.0).or.(nbc.eq.1).or.(nbc.eq.4))
call c_f_pointer(pnumcontrib, &
625 call c_f_pointer(pnofspecies, numberofspecies)
626 call c_f_pointer(pparticlespecies, particlespecies, [n])
627 call c_f_pointer(pcoor, coords, [dim,n])
628 call c_f_pointer(pcutoff, cutoff)
629 if ((nbc.eq.4).or.(nbc.eq.5))
call c_f_pointer(pboxsidelengths, &
630 boxsidelengths, [dim])
631 call c_f_pointer(penergy, energy)
632 call c_f_pointer(pforces, forces, [dim,n])
636 numberofparticles = n
637 if ((nbc.eq.0).or.(nbc.eq.1).or.(nbc.eq.4)) numcontrib = n
638 numberofspecies = aspecies
639 particlespecies(:) = kim_api_get_species_code(pkim,
"Ar", ier)
640 if (ier.lt.kim_status_ok)
then 641 idum = kim_api_report_error(__line__, this_file_name, &
642 "kim_api_get_species_code", ier)
647 if (nbc.eq.4 .or. nbc.eq.5) boxsidelengths(:) = 600.0_cd
650 print
'(3A20)',
"Energy",
"Force Norm",
"Lattice Spacing" 652 current_spacing = min_spacing
653 do while (current_spacing < max_spacing)
661 elseif (nbc.eq.1)
then 664 elseif (nbc.eq.2)
then 667 elseif (nbc.eq.3)
then 670 elseif (nbc.eq.4)
then 672 boxsidelengths, neighobject)
673 elseif (nbc.eq.5)
then 675 boxsidelengths, neighobject)
681 ier = kim_api_model_compute(pkim)
682 if (ier.lt.kim_status_ok)
then 683 idum = kim_api_report_error(__line__, this_file_name, &
684 "kim_api_model_compute", ier)
692 force_norm = force_norm + forces(j,i)*forces(j,i)
695 force_norm = sqrt(force_norm)
699 print
'(3ES20.10)', energy, force_norm, current_spacing
701 current_spacing = current_spacing + spacing_incr
705 if (nbc.lt.6)
deallocate( neighobject%neighborList )
706 if (nbc.eq.0 .or. nbc.eq.2)
then 707 deallocate( neighobject%RijList )
710 ier = kim_api_model_destroy(pkim)
711 if (ier.lt.kim_status_ok)
then 712 idum = kim_api_report_error(__line__, this_file_name, &
713 "kim_api_model_destroy", ier)
716 call kim_api_free(pkim, ier)
717 if (ier.lt.kim_status_ok)
then 718 idum = kim_api_report_error(__line__, this_file_name, &
subroutine neigh_pure_cluster_neighborlist(half, numberOfParticles, coords, cutoff, neighObject)
subroutine create_fcc_configuration(FCCspacing, nCellsPerSide, periodic, coords, MiddlePartId)
subroutine mi_opbc_cluster_neighborlist(half, numberOfParticles, coords, rcut, boxSideLengths, neighObject)
program ex_test_ar_fcc_cluster_fortran
subroutine neigh_rvec_cluster_neighborlist(half, numberOfParticles, coords, cutoff, neighObject)
integer(c_int) function, public get_neigh(pkim, mode, request, part, numnei, pnei1part, pRij)