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(:,:)
63 integer(c_int) function get_neigh(pkim,mode,request,part,numnei,pnei1part, &
68 type(c_ptr),
intent(in) :: pkim
69 integer(c_int),
intent(in) :: mode
70 integer(c_int),
intent(in) :: request
71 integer(c_int),
intent(out) :: part
72 integer(c_int),
intent(out) :: numnei
73 type(c_ptr),
intent(out) :: pnei1part
74 type(c_ptr),
intent(out) :: pRij
77 integer(c_int),
parameter :: DIM = 3
79 integer(c_int) partToReturn
80 integer(c_int),
pointer :: numberOfParticles;
type(c_ptr) :: pnParts
81 type(neighObject_type),
pointer :: neighObject;
type(c_ptr) :: pneighObject
82 integer(c_int) ier, idum
85 pnparts = kim_api_get_data(pkim,
"numberOfParticles", ier)
86 if (ier.lt.kim_status_ok)
then 87 idum = kim_api_report_error(__line__, this_file_name, &
88 "kim_api_get_data", ier)
91 call c_f_pointer(pnparts, numberofparticles)
94 pneighobject = kim_api_get_data(pkim,
"neighObject", ier)
95 if (ier.lt.kim_status_ok)
then 96 idum = kim_api_report_error(__line__, this_file_name, &
97 "kim_api_get_data", ier)
100 call c_f_pointer(pneighobject, neighobject)
102 n =
size(neighobject%neighborList, 2)
106 if ( (request.gt.n) .or. (request.lt.1))
then 107 idum = kim_api_report_error(__line__, this_file_name, &
108 "Invalid part ID in get_neigh", &
109 kim_status_particle_invalid_id)
110 get_neigh = kim_status_particle_invalid_id
113 parttoreturn = request
116 idum = kim_api_report_error(__line__, this_file_name, &
117 "Invalid mode in get_neigh", &
118 kim_status_neigh_invalid_mode)
119 get_neigh = kim_status_neigh_invalid_mode
127 numnei = neighobject%neighborList(1,part)
130 pnei1part = c_loc(neighobject%neighborList(2,part))
154 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
165 type(neighObject_type),
intent(inout) :: neighObject
168 integer(c_int) i, j, a
171 real(c_double) cutoff2
175 do i=1,numberofparticles
177 do j=1,numberofparticles
178 dx(:) = coords(:, j) - coords(:, i)
179 r2 = dot_product(dx, dx)
180 if (r2.le.cutoff2)
then 182 if ( (j .gt. i) .OR. ((.not. half) .AND. (i.ne.j)) )
then 184 neighobject%neighborList(a,i) = j
189 neighobject%neighborList(1,i) = a-1
214 coords, MiddlePartId)
215 use,
intrinsic :: iso_c_binding
217 integer(c_int),
parameter :: cd = c_double
220 real(c_double),
intent(in) :: FCCspacing
221 integer(c_int),
intent(in) :: nCellsPerSide
222 logical,
intent(in) :: periodic
223 real(c_double),
intent(out) :: coords(3,*)
224 integer(c_int),
intent(out) :: MiddlePartId
228 real(c_double) FCCshifts(3,4)
229 real(c_double) latVec(3)
230 integer(c_int) a, i, j, k, m
234 fccshifts(1,1) = 0.0_cd
235 fccshifts(2,1) = 0.0_cd
236 fccshifts(3,1) = 0.0_cd
237 fccshifts(1,2) = 0.5_cd*fccspacing
238 fccshifts(2,2) = 0.5_cd*fccspacing
239 fccshifts(3,2) = 0.0_cd
240 fccshifts(1,3) = 0.5_cd*fccspacing
241 fccshifts(2,3) = 0.0_cd
242 fccshifts(3,3) = 0.5_cd*fccspacing
243 fccshifts(1,4) = 0.0_cd
244 fccshifts(2,4) = 0.5_cd*fccspacing
245 fccshifts(3,4) = 0.5_cd*fccspacing
250 latvec(1) = (i-1)*fccspacing
252 latvec(2) = (j-1)*fccspacing
254 latvec(3) = (k-1)*fccspacing
257 coords(:,a) = latvec + fccshifts(:,m)
258 if ((i.eq.ncellsperside/2+1).and.(j.eq.ncellsperside/2+1) .and. &
259 (k.eq.ncellsperside/2+1) .and. (m.eq.1))
then 260 coords(:,1) = latvec + fccshifts(:,m)
265 if (.not. periodic)
then 268 latvec(1) = ncellsperside*fccspacing
269 latvec(2) = (i-1)*fccspacing
270 latvec(3) = (j-1)*fccspacing
271 a = a+1; coords(:,a) = latvec
272 a = a+1; coords(:,a) = latvec + fccshifts(:,4)
274 latvec(1) = (i-1)*fccspacing
275 latvec(2) = ncellsperside*fccspacing
276 latvec(3) = (j-1)*fccspacing
277 a = a+1; coords(:,a) = latvec
278 a = a+1; coords(:,a) = latvec + fccshifts(:,3)
280 latvec(1) = (i-1)*fccspacing
281 latvec(2) = (j-1)*fccspacing
282 latvec(3) = ncellsperside*fccspacing
283 a = a+1; coords(:,a) = latvec
284 a = a+1; coords(:,a) = latvec + fccshifts(:,2)
287 if (.not. periodic)
then 289 latvec(1) = (i-1)*fccspacing
290 latvec(2) = ncellsperside*fccspacing
291 latvec(3) = ncellsperside*fccspacing
292 a = a+1; coords(:,a) = latvec
293 latvec(1) = ncellsperside*fccspacing
294 latvec(2) = (i-1)*fccspacing
295 latvec(3) = ncellsperside*fccspacing
296 a = a+1; coords(:,a) = latvec
297 latvec(1) = ncellsperside*fccspacing
298 latvec(2) = ncellsperside*fccspacing
299 latvec(3) = (i-1)*fccspacing
300 a = a+1; coords(:,a) = latvec
303 if (.not. periodic)
then 305 a = a+1; coords(:,a) = ncellsperside*fccspacing
328 use,
intrinsic :: iso_c_binding
333 integer(c_int),
parameter :: cd = c_double
335 integer(c_int),
parameter :: ncellsperside = 2
336 integer(c_int),
parameter :: dim = 3
337 integer(c_int),
parameter :: aspecies = 1
339 real(c_double),
parameter :: cutpad = 0.75_cd
340 real(c_double),
parameter :: fccspacing = 5.260_cd
341 real(c_double),
parameter :: min_spacing = 0.8*fccspacing
342 real(c_double),
parameter :: max_spacing = 1.2*fccspacing
343 real(c_double),
parameter :: spacing_incr = 0.025*fccspacing
344 real(c_double) :: current_spacing
345 real(c_double) :: force_norm
347 integer(c_int),
parameter :: &
348 n = 4*(ncellsperside)**3 + 6*(ncellsperside)**2 + 3*(ncellsperside) + 1
349 integer(c_int),
parameter :: sizeone = 1
356 character(len=KIM_KEY_STRING_LENGTH) :: testkimfile =
"descriptor.kim" 357 character(len=KIM_KEY_STRING_LENGTH) :: modelname
358 character(len=KIM_KEY_STRING_LENGTH) :: nbc_method
361 integer(c_int),
pointer :: numberofparticles;
type(c_ptr) :: pnparts
362 integer(c_int),
pointer :: numberofspecies;
type(c_ptr) :: pnofspecies
363 integer(c_int),
pointer :: particlespecies(:);
type(c_ptr) :: pparticlespecies
364 real(c_double),
pointer :: cutoff;
type(c_ptr) :: pcutoff
365 real(c_double),
pointer :: energy;
type(c_ptr) :: penergy
366 real(c_double),
pointer :: coords(:,:);
type(c_ptr) :: pcoor
367 real(c_double),
pointer :: forces(:,:);
type(c_ptr) :: pforces
370 integer(c_int) ier, idum
378 print
'("Please enter a valid KIM model name: ")' 384 print *,
'This is Test : ex_test_Ar_fcc_cluster_fortran' 387 print
'("Results for KIM Model : ",A)', trim(modelname)
390 ier = kim_api_file_init(pkim, testkimfile, modelname)
391 if (ier.lt.kim_status_ok)
then 392 idum = kim_api_report_error(__line__, this_file_name, &
393 "kim_api_file_init", ier)
398 ier = kim_api_get_nbc_method(pkim, nbc_method)
399 if (ier.lt.kim_status_ok)
then 400 idum = kim_api_report_error(__line__, this_file_name, &
401 "kim_api_get_nbc_method", ier)
404 if (index(nbc_method,
"NEIGH_PURE_F").ne.1)
then 405 ier = kim_status_fail
406 idum = kim_api_report_error(__line__, this_file_name, &
407 "Unknown NBC method", ier)
412 call kim_api_allocate(pkim, n, aspecies, ier)
413 if (ier.lt.kim_status_ok)
then 414 idum = kim_api_report_error(__line__, this_file_name, &
415 "kim_api_allocate", ier)
421 allocate(neighobject%neighborList(n+1,n))
425 ier = kim_api_set_data(pkim,
"neighObject", sizeone, c_loc(neighobject))
426 if (ier.lt.kim_status_ok)
then 427 idum = kim_api_report_error(__line__, this_file_name, &
428 "kim_api_set_data", ier)
432 ier = kim_api_set_method(pkim,
"get_neigh", sizeone, c_funloc(
get_neigh))
433 if (ier.lt.kim_status_ok)
then 434 idum = kim_api_report_error(__line__, this_file_name, &
435 "kim_api_set_method", ier)
440 ier = kim_api_model_init(pkim)
441 if (ier.lt.kim_status_ok)
then 442 idum = kim_api_report_error(__line__, this_file_name, &
443 "kim_api_model_init", ier)
450 call kim_api_getm_data(pkim, ier, &
451 "numberOfParticles", pnparts, 1, &
452 "numberOfSpecies", pnofspecies, 1, &
453 "particleSpecies", pparticlespecies, 1, &
454 "coordinates", pcoor, 1, &
455 "cutoff", pcutoff, 1, &
456 "energy", penergy, 1, &
457 "forces", pforces, 1)
458 if (ier.lt.kim_status_ok)
then 459 idum = kim_api_report_error(__line__, this_file_name, &
460 "kim_api_getm_data", ier)
463 call c_f_pointer(pnparts, numberofparticles)
464 call c_f_pointer(pnofspecies, numberofspecies)
465 call c_f_pointer(pparticlespecies, particlespecies, [n])
466 call c_f_pointer(pcoor, coords, [dim,n])
467 call c_f_pointer(pcutoff, cutoff)
468 call c_f_pointer(penergy, energy)
469 call c_f_pointer(pforces, forces, [dim,n])
473 numberofparticles = n
474 numberofspecies = aspecies
475 particlespecies(:) = kim_api_get_species_code(pkim,
"Ar", ier)
476 if (ier.lt.kim_status_ok)
then 477 idum = kim_api_report_error(__line__, this_file_name, &
478 "kim_api_get_species_code", ier)
483 print
'(3A20)',
"Energy",
"Force Norm",
"Lattice Spacing" 485 current_spacing = min_spacing
486 do while (current_spacing < max_spacing)
495 ier = kim_api_model_compute(pkim)
496 if (ier.lt.kim_status_ok)
then 497 idum = kim_api_report_error(__line__, this_file_name, &
498 "kim_api_model_compute", ier)
506 force_norm = force_norm + forces(j,i)*forces(j,i)
509 force_norm = sqrt(force_norm)
513 print
'(3ES20.10)', energy, force_norm, current_spacing
515 current_spacing = current_spacing + spacing_incr
519 deallocate( neighobject%neighborList )
521 ier = kim_api_model_destroy(pkim)
522 if (ier.lt.kim_status_ok)
then 523 idum = kim_api_report_error(__line__, this_file_name, &
524 "kim_api_model_destroy", ier)
527 call kim_api_free(pkim, ier)
528 if (ier.lt.kim_status_ok)
then 529 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)
program ex_test_ar_fcc_cluster_fortran
integer(c_int) function, public get_neigh(pkim, mode, request, part, numnei, pnei1part, pRij)