KIM API V2
ex_test_Ar_fcc_cluster_fortran.F90
Go to the documentation of this file.
1 !
2 ! CDDL HEADER START
3 !
4 ! The contents of this file are subject to the terms of the Common Development
5 ! and Distribution License Version 1.0 (the "License").
6 !
7 ! You can obtain a copy of the license at
8 ! http://www.opensource.org/licenses/CDDL-1.0. See the License for the
9 ! specific language governing permissions and limitations under the License.
10 !
11 ! When distributing Covered Code, include this CDDL HEADER in each file and
12 ! include the License file in a prominent location with the name LICENSE.CDDL.
13 ! If applicable, add the following below this CDDL HEADER, with the fields
14 ! enclosed by brackets "[]" replaced with your own identifying information:
15 !
16 ! Portions Copyright (c) [yyyy] [name of copyright owner]. All rights reserved.
17 !
18 ! CDDL HEADER END
19 !
20 
21 !
22 ! Copyright (c) 2013--2018, Regents of the University of Minnesota.
23 ! All rights reserved.
24 !
25 ! Contributors:
26 ! Ellad B. Tadmor
27 ! Ryan S. Elliott
28 ! Stephen M. Whalen
29 !
30 
31 
32 #include "KIM_API_status.h"
33 #define THIS_FILE_NAME __FILE__
34 #define TRUEFALSE(TRUTH) merge(1,0,(TRUTH))
35 
36 !-------------------------------------------------------------------------------
37 !
38 ! module mod_neighborlist :
39 !
40 ! Module contains type and routines related to neighbor list calculation
41 !
42 !-------------------------------------------------------------------------------
43 
45 
46  use, intrinsic :: iso_c_binding
47  use kim_api_f03
48 
49  public get_neigh
50 
52  integer(c_int), pointer :: neighborlist(:,:)
53  real(c_double), pointer :: rijlist(:,:,:)
54  end type neighobject_type
55 contains
56 
57 !-------------------------------------------------------------------------------
58 !
59 ! get_neigh neighbor list access function
60 !
61 ! This function implements Locator and Iterator mode
62 !
63 !-------------------------------------------------------------------------------
64 integer(c_int) function get_neigh(pkim,mode,request,part,numnei,pnei1part, &
65  pRij) bind(c)
66  implicit none
67 
68  !-- Transferred variables
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
76 
77  !-- Local variables
78  integer(c_int), parameter :: dim = 3
79  integer(c_int), save :: iterval = 0
80  integer(c_int) n
81  integer(c_int) parttoreturn
82  integer(c_int), pointer :: numberofparticles; type(c_ptr) :: pnparts
83  type(neighobject_type), pointer :: neighobject; type(c_ptr) :: pneighobject
84  integer(c_int) ier, idum
85 
86  ! unpack number of particles
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)
91  stop
92  endif
93  call c_f_pointer(pnparts, numberofparticles)
94 
95  ! unpack neighbor list object
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)
100  stop
101  endif
102  call c_f_pointer(pneighobject, neighobject)
103 
104  n = size(neighobject%neighborList, 2)
105 
106  ! check mode and request
107  if (mode.eq.0) then ! iterator mode
108  if (request.eq.0) then ! reset iterator
109  iterval = 0
110  get_neigh = kim_status_neigh_iter_init_ok
111  return
112  elseif (request.eq.1) then ! increment iterator
113  iterval = iterval + 1
114  if (iterval.gt.n) then
115  get_neigh = kim_status_neigh_iter_past_end
116  return
117  else
118  parttoreturn = iterval
119  endif
120  else
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
125  return
126  endif
127  elseif (mode.eq.1) then ! locator mode
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
133  return
134  else
135  parttoreturn = request
136  endif
137  else ! not iterator or locator mode
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
142  return
143  endif
144 
145  ! set the returned part
146  part = parttoreturn
147 
148  ! set the returned number of neighbors for the returned part
149  numnei = neighobject%neighborList(1,part)
150 
151  ! set the location for the returned neighbor list
152  pnei1part = c_loc(neighobject%neighborList(2,part))
153 
154  ! set pointer to Rij to appropriate value
155  if (associated(neighobject%RijList)) then
156  prij = c_loc(neighobject%RijList(1,1,part))
157  else
158  prij = c_null_ptr
159  endif
160 
161  get_neigh = kim_status_ok
162  return
163 end function get_neigh
164 
165 end module mod_neighborlist
166 
167 
169  implicit none
170  public
171 
172 contains
173 
174 !-------------------------------------------------------------------------------
175 !
176 ! MI_OPBC_cluster_neighborlist : construct a half or full neighbor list using
177 ! the particle coordinates in coords()
178 !
179 !-------------------------------------------------------------------------------
180 subroutine mi_opbc_cluster_neighborlist(half, numberOfParticles, coords, rcut, &
181  boxSideLengths, neighObject)
182  use, intrinsic :: iso_c_binding
183  use kim_api_f03
184  use mod_neighborlist
185  implicit none
186  integer(c_int), parameter :: cd = c_double ! used for literal constants
187 
188  !-- Transferred variables
189  logical, intent(in) :: half
190  integer(c_int), intent(in) :: numberOfParticles
191  real(c_double), dimension(3,numberOfParticles), &
192  intent(in) :: coords
193  real(c_double), intent(in) :: rcut
194  real(c_double), dimension(3), intent(in) :: boxSideLengths
195  type(neighObject_type), intent(inout) :: neighObject
196 
197  !-- Local variables
198  integer(c_int) i, j, a
199  real(c_double) dx(3)
200  real(c_double) r2
201  real(c_double) rcut2
202 
203  rcut2 = rcut**2
204 
205  do i=1,numberofparticles
206  a = 1
207  do j=1,numberofparticles
208  dx(:) = coords(:, j) - coords(:, i)
209  where (abs(dx) > 0.5_cd*boxsidelengths) ! apply PBC
210  dx = dx - sign(boxsidelengths,dx)
211  endwhere
212  r2 = dot_product(dx, dx)
213  if (r2.le.rcut2) then
214  if (i.ne.j) then
215  if ( (j .gt. i) .or. ((.not. half) .AND. (i.ne.j)) ) then
216  ! part j is a neighbor of part i
217  a = a+1
218  neighobject%neighborList(a,i) = j
219  endif
220  endif
221  endif
222  enddo
223  ! part i has a-1 neighbors
224  neighobject%neighborList(1,i) = a-1
225  enddo
226 
227  return
228 
229 end subroutine mi_opbc_cluster_neighborlist
230 
231 !-------------------------------------------------------------------------------
232 !
233 ! NEIGH_PURE_cluster_neighborlist
234 !
235 !-------------------------------------------------------------------------------
236 subroutine neigh_pure_cluster_neighborlist(half, numberOfParticles, coords, &
237  cutoff, neighObject)
238  use, intrinsic :: iso_c_binding
239  use kim_api_f03
240  use mod_neighborlist
241  implicit none
242 
243  !-- Transferred variables
244  logical, intent(in) :: half
245  integer(c_int), intent(in) :: numberOfParticles
246  real(c_double), dimension(3,numberOfParticles), &
247  intent(in) :: coords
248  real(c_double), intent(in) :: cutoff
249  type(neighObject_type), intent(inout) :: neighObject
250 
251  !-- Local variables
252  integer(c_int) i, j, a
253  real(c_double) dx(3)
254  real(c_double) r2
255  real(c_double) cutoff2
256 
257  cutoff2 = cutoff**2
258 
259  do i=1,numberofparticles
260  a = 1
261  do j=1,numberofparticles
262  dx(:) = coords(:, j) - coords(:, i)
263  r2 = dot_product(dx, dx)
264  if (r2.le.cutoff2) then
265  ! part j is a neighbor of part i
266  if ( (j .gt. i) .OR. ((.not. half) .AND. (i.ne.j)) ) then
267  a = a+1
268  neighobject%neighborList(a,i) = j
269  endif
270  endif
271  enddo
272  ! part i has a-1 neighbors
273  neighobject%neighborList(1,i) = a-1
274  enddo
275 
276  return
277 
278 end subroutine neigh_pure_cluster_neighborlist
279 
280 !-------------------------------------------------------------------------------
281 !
282 ! NEIGH_RVEC_cluster_neighborlist
283 !
284 !-------------------------------------------------------------------------------
285 subroutine neigh_rvec_cluster_neighborlist(half, numberOfParticles, coords, &
286  cutoff, neighObject)
287  use, intrinsic :: iso_c_binding
288  use kim_api_f03
289  use mod_neighborlist
290  implicit none
291 
292  !-- Transferred variables
293  logical, intent(in) :: half
294  integer(c_int), intent(in) :: numberOfParticles
295  real(c_double), dimension(3,numberOfParticles), &
296  intent(in) :: coords
297  real(c_double), intent(in) :: cutoff
298  type(neighObject_type), intent(inout) :: neighObject
299 
300  !-- Local variables
301  integer(c_int) i, j, a
302  real(c_double) dx(3)
303  real(c_double) r2
304  real(c_double) cutoff2
305 
306  cutoff2 = cutoff**2
307 
308  do i=1,numberofparticles
309  a = 1
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
315  ! part j is a neighbor of part i
316  a = a+1
317  neighobject%neighborList(a,i) = j
318  neighobject%RijList(:,a-1,i) = dx
319  endif
320  endif
321  enddo
322  ! part i has a-1 neighbors
323  neighobject%neighborList(1,i) = a-1
324  enddo
325 
326  return
327 
328 end subroutine neigh_rvec_cluster_neighborlist
329 
330 !-------------------------------------------------------------------------------
331 !
332 ! create_FCC_configuration subroutine
333 !
334 ! creates a cubic configuration of FCC particles with lattice spacing
335 ! `FCCspacing' and `nCellsPerSide' cells along each direction.
336 !
337 ! With periodic==.true. this will result in a total number of particles equal
338 ! to 4*(nCellsPerSide)**3 + 6*(nCellsPerSide)**2 + 3*(nCellsPerSide) + 1
339 !
340 ! With periodic==.false. this will result in a total number of particles equal
341 ! to 4*(nCellsPerSide)**3
342 !
343 ! Returns the Id of the particle situated in the middle of the configuration
344 ! (this particle will have the most neighbors.)
345 !
346 !-------------------------------------------------------------------------------
347 subroutine create_fcc_configuration(FCCspacing, nCellsPerSide, periodic, &
348  coords, MiddlePartId)
349  use, intrinsic :: iso_c_binding
350  implicit none
351  integer(c_int), parameter :: cd = c_double ! used for literal constants
352 
353  !-- Transferred variables
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
359  !
360  ! cluster setup variables
361  !
362  real(c_double) FCCshifts(3,4)
363  real(c_double) latVec(3)
364  integer(c_int) a, i, j, k, m
365 
366  ! Create a cubic FCC cluster
367  !
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
380 
381  middlepartid = 1 ! Always put middle particle as #1
382  a = 1 ! leave space for middle particle as particle #1
383  do i=1,ncellsperside
384  latvec(1) = (i-1)*fccspacing
385  do j=1,ncellsperside
386  latvec(2) = (j-1)*fccspacing
387  do k=1,ncellsperside
388  latvec(3) = (k-1)*fccspacing
389  do m=1,4
390  a = a+1
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) ! put middle particle as #1
395  a = a - 1
396  endif
397  enddo
398  enddo
399  if (.not. periodic) then
400  ! Add in the remaining three faces
401  ! pos-x face
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)
407  ! pos-y face
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)
413  ! pos-z face
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)
419  endif
420  enddo
421  if (.not. periodic) then
422  ! Add in the remaining three edges
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
435  endif
436  enddo
437  if (.not. periodic) then
438  ! Add in the remaining corner
439  a = a+1; coords(:,a) = ncellsperside*fccspacing
440  endif
441 
442  return
443 
444 end subroutine create_fcc_configuration
445 
446 end module mod_utility
447 
448 !*******************************************************************************
449 !**
450 !** PROGRAM vc_forces_numer_deriv
451 !**
452 !** KIM compliant program to perform numerical derivative check on a model
453 !**
454 !*******************************************************************************
455 
456 !-------------------------------------------------------------------------------
457 !
458 ! Main program
459 !
460 !-------------------------------------------------------------------------------
462  use, intrinsic :: iso_c_binding
463  use kim_api_f03
464  use mod_neighborlist
465  use mod_utility
466  implicit none
467  integer(c_int), parameter :: cd = c_double ! used for literal constants
468 
469  integer(c_int), parameter :: ncellsperside = 2
470  integer(c_int), parameter :: dim = 3
471  integer(c_int), parameter :: aspecies = 1
472 
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
480 
481  integer(c_int), parameter :: &
482  n = 4*(ncellsperside)**3 + 6*(ncellsperside)**2 + 3*(ncellsperside) + 1
483  integer(c_int), parameter :: sizeone = 1
484 
485  type(neighobject_type), target :: neighobject
486 
487  !
488  ! KIM variables
489  !
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
493  ! 0- NEIGH_RVEC_H, 1- NEIGH_PURE_H, 2- NEIGH_RVEC_F, 3- NEIGH_PURE_F,
494  ! 4- MI_OPBC_H, 5- MI_OPBC_F, 6- CLUSTER
495  integer(c_int) nbc
496 
497  type(c_ptr) :: pkim
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
507 
508 
509  integer(c_int) ier, idum
510  integer :: middledum
511  integer :: i, j
512 
513  ! Initialize error flag
514  ier = 0
515 
516  ! Get KIM Model name to use
517  print '("Please enter a valid KIM model name: ")'
518  read(*,*) modelname
519 
520  ! Print output header
521  !
522  print *
523  print *,'This is Test : ex_test_Ar_fcc_cluster_fortran'
524  print *
525  print '(80(''-''))'
526  print '("Results for KIM Model : ",A)', trim(modelname)
527 
528  ! Initialize the KIM object
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)
533  stop
534  endif
535 
536  ! determine which NBC scenerio to use
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)
541  stop
542  endif
543  if (index(nbc_method,"NEIGH_RVEC_H").eq.1) then
544  nbc = 0
545  elseif (index(nbc_method,"NEIGH_PURE_H").eq.1) then
546  nbc = 1
547  elseif (index(nbc_method,"NEIGH_RVEC_F").eq.1) then
548  nbc = 2
549  elseif (index(nbc_method,"NEIGH_PURE_F").eq.1) then
550  nbc = 3
551  elseif (index(nbc_method,"MI_OPBC_H").eq.1) then
552  nbc = 4
553  elseif (index(nbc_method,"MI_OPBC_F").eq.1) then
554  nbc = 5
555  elseif (index(nbc_method,"CLUSTER").eq.1) then
556  nbc = 6
557  else
558  ier = kim_status_fail
559  idum = kim_api_report_error(__line__, this_file_name, &
560  "Unknown NBC method", ier)
561  stop
562  endif
563 
564  ! Allocate memory via the KIM system
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)
569  stop
570  endif
571 
572  ! Allocate storage for neighbor lists
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))
576  endif
577  !
578  if (nbc.ne.6) then
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)
583  stop
584  endif
585  endif
586 
587  if (nbc.ne.6) then
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)
592  stop
593  endif
594  endif
595 
596  ! call model's init routine
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)
601  stop
602  endif
603 
604 
605  ! Unpack data from KIM object
606  !
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)
620  stop
621  endif
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, &
624  numcontrib)
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])
633 
634 
635  ! Set values
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)
643  stop
644  endif
645 
646  ! set boxSideLengths large enough to make the cluster isolated
647  if (nbc.eq.4 .or. nbc.eq.5) boxsidelengths(:) = 600.0_cd
648 
649  ! print header
650  print '(3A20)', "Energy", "Force Norm", "Lattice Spacing"
651  ! do the computations
652  current_spacing = min_spacing
653  do while (current_spacing < max_spacing)
654 
655  call create_fcc_configuration(current_spacing, ncellsperside, .false., &
656  coords, middledum)
657  ! Compute neighbor lists
658  if (nbc.eq.0) then
659  call neigh_rvec_cluster_neighborlist(.true., n, coords, (cutoff+cutpad), &
660  neighobject)
661  elseif (nbc.eq.1) then
662  call neigh_pure_cluster_neighborlist(.true., n, coords, (cutoff+cutpad), &
663  neighobject)
664  elseif (nbc.eq.2) then
665  call neigh_rvec_cluster_neighborlist(.false., n, coords, (cutoff+cutpad), &
666  neighobject)
667  elseif (nbc.eq.3) then
668  call neigh_pure_cluster_neighborlist(.false., n, coords, (cutoff+cutpad), &
669  neighobject)
670  elseif (nbc.eq.4) then
671  call mi_opbc_cluster_neighborlist(.true., n, coords, (cutoff+cutpad), &
672  boxsidelengths, neighobject)
673  elseif (nbc.eq.5) then
674  call mi_opbc_cluster_neighborlist(.false., n, coords, (cutoff+cutpad), &
675  boxsidelengths, neighobject)
676  else
677  ! nothing to do for CLUSTER
678  endif
679 
680  ! Call model compute
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)
685  stop
686  endif
687 
688  ! compue force_norm
689  force_norm = 0.0;
690  do i = 1, n
691  do j = 1, dim
692  force_norm = force_norm + forces(j,i)*forces(j,i)
693  end do
694  end do
695  force_norm = sqrt(force_norm)
696 
697  ! Print results to screen
698  !
699  print '(3ES20.10)', energy, force_norm, current_spacing
700 
701  current_spacing = current_spacing + spacing_incr
702  enddo
703 
704  ! Deallocate neighbor list object
705  if (nbc.lt.6) deallocate( neighobject%neighborList )
706  if (nbc.eq.0 .or. nbc.eq.2) then
707  deallocate( neighobject%RijList )
708  endif
709 
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)
714  stop
715  endif
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, &
719  "kim_api_free", ier)
720  stop
721  endif
722 
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)