kim-api 2.3.0+AppleClang.AppleClang.GNU
An Application Programming Interface (API) for the Knowledgebase of Interatomic Models (KIM).
Loading...
Searching...
No Matches
ex_model_Ar_SLJ_MultiCutoff.f90
Go to the documentation of this file.
1!
2! KIM-API: An API for interatomic models
3! Copyright (c) 2013--2022, Regents of the University of Minnesota.
4! All rights reserved.
5!
6! Contributors:
7! Ellad B. Tadmor
8!
9! SPDX-License-Identifier: LGPL-2.1-or-later
10!
11! This library is free software; you can redistribute it and/or
12! modify it under the terms of the GNU Lesser General Public
13! License as published by the Free Software Foundation; either
14! version 2.1 of the License, or (at your option) any later version.
15!
16! This library is distributed in the hope that it will be useful,
17! but WITHOUT ANY WARRANTY; without even the implied warranty of
18! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19! Lesser General Public License for more details.
20!
21! You should have received a copy of the GNU Lesser General Public License
22! along with this library; if not, write to the Free Software Foundation,
23! Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
24!
25!
26! Copyright (c) 2013--2022, Regents of the University of Minnesota.
27! All rights reserved.
28!
29
30!****************************************************************************
31!**
32!** MODULE ex_model_Ar_SLJ_MultiCutoff
33!**
34!** Spring-modified Lennard-Jones (SLJ) pair potential model for Ar
35!**
36!** V = 0.5 \sum_i \sum_j eps_i eps_j 4 [ (sig/r_ij)^12 - (sig/r_ij)^6 ] (1)
37!**
38!** where
39!** eps_i = 0.5 \sum_k spring (r_ik)^2 (2)
40!**
41!** See README for details.
42!**
43!** Language: Fortran 2003
44!**
45!** Author: Ellad B. Tadmor
46!**
47!** Date: August 23, 2018
48!**
49!****************************************************************************
50
52
53 use, intrinsic :: iso_c_binding
55 implicit none
56
57 save
58 private
59 public compute_energy_forces, &
65 speccode, &
67
68 ! Below are the definitions and values of all Model parameters
69 integer(c_int), parameter :: cd = c_double ! used for literal constants
70 integer(c_int), parameter :: DIM = 3 ! dimensionality of space
71 integer(c_int), parameter :: speccode = 1 ! internal species code
72
73 !-----------------------------------------------------------------------------
74 ! Below are the definitions and values of all additional model parameters
75 !
76 ! Recall that the Fortran 2003 format for declaring parameters is as follows:
77 !
78 ! integer(c_int), parameter :: parname = value ! This defines an integer
79 ! ! parameter called `parname'
80 ! ! with a value equal to
81 ! ! `value' (a number)
82 !
83 ! real(c_double), parameter :: parname = value ! This defines a real(c_double)
84 ! ! parameter called `parname'
85 ! ! with a value equal to
86 ! ! `value' (a number)
87 !-----------------------------------------------------------------------------
88 real(c_double), parameter :: lj_spring = 0.00051226_cd
89 real(c_double), parameter :: lj_sigma = 5.26_cd / 1.74724_cd
90 ! experimental fcc lattice constant is a0=5.26
91 real(c_double), parameter :: model_cutoff1 = 1.25 * lj_sigma
92 ! short-range nearest neighbor for fcc
93 real(c_double), parameter :: model_cutoff2 = 2.25 * lj_sigma
94 ! long-range third neighbor for fcc
95 real(c_double), parameter :: model_cutsq1 = model_cutoff1**2
96 real(c_double), parameter :: model_cutsq2 = model_cutoff2**2
97
98 type, bind(c) :: buffer_type
99 real(c_double) :: influence_distance
100 real(c_double) :: cutoff(2)
101 integer(c_int) :: &
102 model_will_not_request_neighbors_of_noncontributing_particles(2)
103 end type buffer_type
104
105contains
106
107 !-----------------------------------------------------------------------------
108 !
109 ! Calculate Lennard-Jones potential phi(r) and, if requested,
110 ! its derivative dphi(r)
111 !
112 !-----------------------------------------------------------------------------
113 recursive subroutine calc_phi(r, phi, dphi, calc_deriv)
114 implicit none
115
116 !-- Transferred variables
117 real(c_double), intent(in) :: r
118 real(c_double), intent(out) :: phi
119 real(c_double), intent(out) :: dphi
120 logical, intent(in) :: calc_deriv
121
122 !-- Local variables
123 real(c_double) rsq, sor, sor6, sor12
124
125 rsq = r * r ! r^2
126 sor = lj_sigma / r ! (sig/r)
127 sor6 = sor * sor * sor !
128 sor6 = sor6 * sor6 ! (sig/r)^6
129 sor12 = sor6 * sor6 ! (sig/r)^12
130 if (r > model_cutoff2) then
131 ! Argument exceeds cutoff radius
132 phi = 0.0_cd
133 if (calc_deriv) dphi = 0.0_cd
134 else
135 phi = 4.0_cd * (sor12 - sor6)
136 if (calc_deriv) dphi = 24.0_cd * (-2.0_cd * sor12 + sor6) / r
137 end if
138
139 end subroutine calc_phi
140
141 !-----------------------------------------------------------------------------
142 !
143 ! Calculate short-range linear spring-based energy amplitude for `atom`
144 !
145 !-----------------------------------------------------------------------------
146 recursive subroutine calc_spring_energyamp(model_compute_arguments_handle, &
147 atom, coor, eps, ierr)
148 implicit none
149
150 !-- Transferred variables
151 type(kim_model_compute_arguments_handle_type), &
152 intent(in) :: model_compute_arguments_handle
153 integer(c_int), intent(in) :: atom
154 real(c_double), intent(in) :: coor(:, :)
155 real(c_double), intent(out) :: eps
156 integer(c_int), intent(out) :: ierr
157
158 !-- Local variables
159 integer(c_int) k, kk
160 real(c_double) rrel(dim)
161 real(c_double) rsqrel
162 integer(c_int) numneishort
163 integer(c_int), pointer :: neishort(:)
164
165 ! Get short-range neighbors of `atom`
166 call kim_get_neighbor_list( &
167 model_compute_arguments_handle, 1, atom, numneishort, neishort, ierr)
168 if (ierr /= 0) return
169
170 eps = 0.0_cd
171 do kk = 1, numneishort
172 k = neishort(kk)
173 rrel(:) = coor(:, k) - coor(:, atom)
174 rsqrel = dot_product(rrel, rrel)
175 if (rsqrel < model_cutsq1) eps = eps + rsqrel
176 end do
177 eps = 0.5_cd * lj_spring * eps
178
179 end subroutine calc_spring_energyamp
180
181 !-----------------------------------------------------------------------------
182 !
183 ! Calculate short-range linear spring-based contribution to force
184 !
185 !-----------------------------------------------------------------------------
186 recursive subroutine calc_spring_force(model_compute_arguments_handle, atom, &
187 coor, eps, phi, force, ierr)
188 implicit none
189
190 !-- Transferred variables
191 type(kim_model_compute_arguments_handle_type), &
192 intent(in) :: model_compute_arguments_handle
193 integer(c_int), intent(in) :: atom
194 real(c_double), intent(in) :: coor(:, :)
195 real(c_double), intent(in) :: eps
196 real(c_double), intent(in) :: phi
197 real(c_double), intent(inout) :: force(:, :)
198 integer(c_int), intent(out) :: ierr
199
200 !-- Local variables
201 integer(c_int) k, kk
202 real(c_double) rrel(dim), dforce(dim)
203 real(c_double) rsqrel
204 integer(c_int) numneishort
205 integer(c_int), pointer :: neishort(:)
206
207 ! Get short-range neighbors of `atom`
208 call kim_get_neighbor_list( &
209 model_compute_arguments_handle, 1, atom, numneishort, neishort, ierr)
210 if (ierr /= 0) return
211
212 ! Add contribution to force on `atom` and its near neighbors that contribute
213 ! to the spring term
214 do kk = 1, numneishort
215 k = neishort(kk)
216 rrel(:) = coor(:, k) - coor(:, atom)
217 rsqrel = dot_product(rrel, rrel)
218 if (rsqrel < model_cutsq1) then
219 dforce(:) = 0.5_cd * eps * lj_spring * rrel(:) * phi
220 force(:, atom) = force(:, atom) + dforce(:) ! accumulate force on atom
221 force(:, k) = force(:, k) - dforce(:) ! accumulate force on k
222 end if
223 end do
224
225 end subroutine calc_spring_force
226
227 !-----------------------------------------------------------------------------
228 !
229 ! Compute energy and forces on particles from the positions.
230 !
231 !-----------------------------------------------------------------------------
232 recursive subroutine compute_energy_forces( &
233 model_compute_handle, model_compute_arguments_handle, ierr) bind(c)
234 implicit none
235
236 !-- Transferred variables
237 type(kim_model_compute_handle_type), intent(in) :: model_compute_handle
238 type(kim_model_compute_arguments_handle_type), intent(in) :: &
239 model_compute_arguments_handle
240 integer(c_int), intent(out) :: ierr
241
242 !-- Local variables
243 real(c_double) :: rij(dim), dforce(dim)
244 real(c_double) :: r, rsqij, phi, dphi, deidr, epsi, epsj
245 integer(c_int) :: i, j, jj, comp_force, comp_enepot, comp_energy
246 !integer(c_int) :: comp_virial
247 integer(c_int) :: numnei
248 integer(c_int) :: ierr2
249 logical :: calc_deriv
250
251 !-- KIM variables
252 integer(c_int), pointer :: n
253 real(c_double), pointer :: energy
254 real(c_double), pointer :: coor(:, :)
255 real(c_double), pointer :: force(:, :)
256 real(c_double), pointer :: enepot(:)
257 integer(c_int), pointer :: nei1part(:)
258 integer(c_int), pointer :: particlespeciescodes(:)
259 integer(c_int), pointer :: particlecontributing(:)
260 !real(c_double), pointer :: virial(:)
261
262 ! Unpack data from KIM object
263 !
264 ierr = 0
265 call kim_get_argument_pointer( &
266 model_compute_arguments_handle, &
267 kim_compute_argument_name_number_of_particles, n, ierr2)
268 ierr = ierr + ierr2
269 call kim_get_argument_pointer( &
270 model_compute_arguments_handle, &
271 kim_compute_argument_name_particle_species_codes, n, &
272 particlespeciescodes, ierr2)
273 ierr = ierr + ierr2
274 call kim_get_argument_pointer( &
275 model_compute_arguments_handle, &
276 kim_compute_argument_name_particle_contributing, n, &
277 particlecontributing, ierr2)
278 ierr = ierr + ierr2
279 call kim_get_argument_pointer( &
280 model_compute_arguments_handle, &
281 kim_compute_argument_name_coordinates, dim, n, coor, ierr2)
282 ierr = ierr + ierr2
283 call kim_get_argument_pointer( &
284 model_compute_arguments_handle, &
285 kim_compute_argument_name_partial_energy, energy, ierr2)
286 ierr = ierr + ierr2
287 call kim_get_argument_pointer( &
288 model_compute_arguments_handle, &
289 kim_compute_argument_name_partial_forces, dim, n, force, ierr2)
290 ierr = ierr + ierr2
291 call kim_get_argument_pointer( &
292 model_compute_arguments_handle, &
293 kim_compute_argument_name_partial_particle_energy, n, enepot, ierr2)
294 ierr = ierr + ierr2
295 !call kim_model_compute_arguments_get_argument_pointer( &
296 ! model_compute_arguments_handle, &
297 ! KIM_COMPUTE_ARGUMENT_NAME_PARTIAL_VIRIAL, 6, virial, ierr2)
298 !ierr = ierr + ierr2
299 if (ierr /= 0) then
300 call kim_log_entry(model_compute_arguments_handle, &
301 kim_log_verbosity_error, "get data")
302 return
303 end if
304
305 ! Check to see if we have been asked to compute the forces, energyperpart,
306 ! energy and virial
307 !
308 if (associated(energy)) then
309 comp_energy = 1
310 else
311 comp_energy = 0
312 end if
313 if (associated(force)) then
314 comp_force = 1
315 else
316 comp_force = 0
317 end if
318 if (associated(enepot)) then
319 comp_enepot = 1
320 else
321 comp_enepot = 0
322 end if
323 !if (associated(virial)) then
324 ! comp_virial = 1
325 !else
326 ! comp_virial = 0
327 !end if
328 calc_deriv = comp_force == 1 !.or.comp_virial.eq.1
329
330 ! Check to be sure that the species are correct
331 !
332 ierr = 1 ! assume an error
333 do i = 1, n
334 if (particlespeciescodes(i) /= speccode) then
335 call kim_log_entry( &
336 model_compute_handle, kim_log_verbosity_error, &
337 "Unexpected species code detected")
338 return
339 end if
340 end do
341 ierr = 0 ! everything is ok
342
343 ! Initialize potential energies, forces, virial term
344 !
345 if (comp_enepot == 1) enepot = 0.0_cd
346 if (comp_energy == 1) energy = 0.0_cd
347 if (comp_force == 1) force = 0.0_cd
348 !if (comp_virial.eq.1) virial = 0.0_cd
349 if (calc_deriv) deidr = 0.0_cd
350
351 !
352 ! Compute energy and forces
353 !
354
355 ! Loop over particles and compute energy and forces
356 !
357 do i = 1, n
358 if (particlecontributing(i) == 1) then
359 ! Set up neighbor list for next particle
360 call kim_get_neighbor_list( &
361 model_compute_arguments_handle, 2, i, numnei, nei1part, ierr)
362 if (ierr /= 0) then
363 ! some sort of problem, exit
364 call kim_log_entry( &
365 model_compute_arguments_handle, kim_log_verbosity_error, &
366 "GetNeighborList failed")
367 ierr = 1
368 return
369 end if
370
371 ! Get short range contribution for atom i to energy amplitude
372 call calc_spring_energyamp(model_compute_arguments_handle, &
373 i, coor, epsi, ierr)
374 if (ierr /= 0) then
375 ! some sort of problem, exit
376 call kim_log_entry( &
377 model_compute_handle, kim_log_verbosity_error, &
378 "GetNeighborList failed")
379 ierr = 1
380 return
381 end if
382
383 ! Loop over the neighbors of particle i
384 !
385 do jj = 1, numnei
386
387 j = nei1part(jj) ! get neighbor ID
388
389 ! Get short range contribution for atom j to energy amplitude
390 call calc_spring_energyamp(model_compute_arguments_handle, j, coor, &
391 epsj, ierr)
392 if (ierr /= 0) then
393 ! some sort of problem, exit
394 call kim_log_entry( &
395 model_compute_handle, kim_log_verbosity_error, &
396 "GetNeighborList failed")
397 ierr = 1
398 return
399 end if
400
401 ! compute relative position vector
402 !
403 rij(:) = coor(:, j) - coor(:, i) ! distance vector between i j
404
405 ! compute energy and forces
406 !
407 rsqij = dot_product(rij, rij) ! compute square distance
408 if (rsqij < model_cutsq2) then ! particles are interacting?
409
410 r = sqrt(rsqij) ! compute distance
411 call calc_phi(r, phi, dphi, calc_deriv) ! compute pair potential and deriv
412 if (calc_deriv) deidr = 0.5_cd * dphi
413
414 ! contribution to energy
415 !
416 if (comp_enepot == 1) then
417 enepot(i) = enepot(i) + 0.5_cd * epsi * epsj * phi ! accumulate energy
418 end if
419 if (comp_energy == 1) then
420 energy = energy + 0.5_cd * epsi * epsj * phi
421 end if
422
423 !!@@@@@@@@@@@@@@@@@@@@ NOT FIXED YET
424 ! ! contribution to virial tensor,
425 ! ! virial(i,j)=r(i)*r(j)*(dV/dr)/r
426 ! !
427 ! if (comp_virial.eq.1) then
428 ! virial(1) = virial(1) + Rij(1)*Rij(1)*dEidr/r
429 ! virial(2) = virial(2) + Rij(2)*Rij(2)*dEidr/r
430 ! virial(3) = virial(3) + Rij(3)*Rij(3)*dEidr/r
431 ! virial(4) = virial(4) + Rij(2)*Rij(3)*dEidr/r
432 ! virial(5) = virial(5) + Rij(1)*Rij(3)*dEidr/r
433 ! virial(6) = virial(6) + Rij(1)*Rij(2)*dEidr/r
434 ! endif
435 !!@@@@@@@@@@@@@@@@@@@@
436
437 ! contribution to forces
438 !
439 if (comp_force == 1) then
440 ! Contribution due to short range neighbors of i
441 call calc_spring_force(model_compute_arguments_handle, i, coor, &
442 epsj, phi, force, ierr)
443 if (ierr /= 0) then
444 ! some sort of problem, exit
445 call kim_log_entry( &
446 model_compute_handle, kim_log_verbosity_error, &
447 "GetNeighborList failed")
448 ierr = 1
449 return
450 end if
451 ! Contribution due to short range neighbors of j
452 call calc_spring_force(model_compute_arguments_handle, j, coor, &
453 epsi, phi, force, ierr)
454 if (ierr /= 0) then
455 ! some sort of problem, exit
456 call kim_log_entry( &
457 model_compute_handle, kim_log_verbosity_error, &
458 "GetNeighborList failed")
459 ierr = 1
460 return
461 end if
462 ! Contribution due to deriv of LJ term
463 dforce(:) = epsi * epsj * deidr * rij(:) / r
464 force(:, i) = force(:, i) + dforce(:) ! accumulate force on i
465 force(:, j) = force(:, j) - dforce(:) ! accumulate force on j
466 end if
467
468 end if
469
470 end do ! loop on jj
471
472 end if ! if particleContributing
473
474 end do ! do i
475
476 ! Everything is great
477 !
478 ierr = 0
479 return
480
481 end subroutine compute_energy_forces
482
483 !-----------------------------------------------------------------------------
484 !
485 ! Model destroy routine (REQUIRED)
486 !
487 !-----------------------------------------------------------------------------
488 recursive subroutine model_destroy_func(model_destroy_handle, ierr) bind(c)
489 use, intrinsic :: iso_c_binding
490 implicit none
491
492 !-- Transferred variables
493 type(kim_model_destroy_handle_type), intent(inout) :: model_destroy_handle
494 integer(c_int), intent(out) :: ierr
495
496 type(buffer_type), pointer :: buf; type(c_ptr) :: pbuf
497
498 call kim_get_model_buffer_pointer(model_destroy_handle, pbuf)
499 call c_f_pointer(pbuf, buf)
500 call kim_log_entry(model_destroy_handle, kim_log_verbosity_information, &
501 "deallocating model buffer")
502 deallocate (buf)
503 ierr = 0 ! everything is good
504 end subroutine model_destroy_func
505
506 !-----------------------------------------------------------------------------
507 !
508 ! Model compute arguments create routine (REQUIRED)
509 !
510 !-----------------------------------------------------------------------------
511 recursive subroutine model_compute_arguments_create( &
512 model_compute_handle, model_compute_arguments_create_handle, ierr) bind(c)
513 use, intrinsic :: iso_c_binding
514 implicit none
515
516 !-- Transferred variables
517 type(kim_model_compute_handle_type), intent(in) :: model_compute_handle
518 type(kim_model_compute_arguments_create_handle_type), intent(inout) :: &
519 model_compute_arguments_create_handle
520 integer(c_int), intent(out) :: ierr
521
522 integer(c_int) :: ierr2
523
524 ! avoid unsed dummy argument warnings
525 if (model_compute_handle == kim_model_compute_null_handle) continue
526
527 ierr = 0
528 ierr2 = 0
529
530 ! register arguments
531 call kim_set_argument_support_status( &
532 model_compute_arguments_create_handle, &
533 kim_compute_argument_name_partial_energy, &
534 kim_support_status_optional, ierr2)
535 ierr = ierr + ierr2
536 call kim_set_argument_support_status( &
537 model_compute_arguments_create_handle, &
538 kim_compute_argument_name_partial_forces, &
539 kim_support_status_optional, ierr2)
540 ierr = ierr + ierr2
541 call kim_set_argument_support_status( &
542 model_compute_arguments_create_handle, &
543 kim_compute_argument_name_partial_particle_energy, &
544 kim_support_status_optional, ierr2)
545 ierr = ierr + ierr2
546 ! call kim_set_argument_support_status( &
547 ! model_compute_arguments_create_handle, &
548 ! KIM_COMPUTE_ARGUMENT_NAME_PARTIAL_VIRIAL, &
549 ! KIM_SUPPORT_STATUS_OPTIONAL, ierr2)
550 ! ierr = ierr + ierr2
551
552 ! register call backs
553 ! NONE
554
555 if (ierr /= 0) then
556 ierr = 1
557 call kim_log_entry( &
558 model_compute_arguments_create_handle, &
559 kim_log_verbosity_error, &
560 "Unable to successfully create compute_arguments object")
561 end if
562
563 return
564 end subroutine model_compute_arguments_create
565
566 !-----------------------------------------------------------------------------
567 !
568 ! Model compute arguments destroy routine (REQUIRED)
569 !
570 !-----------------------------------------------------------------------------
571 recursive subroutine model_compute_arguments_destroy( &
572 model_compute_handle, model_compute_arguments_destroy_handle, ierr) bind(c)
573 use, intrinsic :: iso_c_binding
574 implicit none
575
576 !-- Transferred variables
577 type(kim_model_compute_handle_type), intent(in) :: model_compute_handle
578 type(kim_model_compute_arguments_destroy_handle_type), intent(inout) :: &
579 model_compute_arguments_destroy_handle
580 integer(c_int), intent(out) :: ierr
581
582 ! avoid unsed dummy argument warnings
583 if (model_compute_handle == kim_model_compute_null_handle) continue
584 if (model_compute_arguments_destroy_handle == &
585 kim_model_compute_arguments_destroy_null_handle) continue
586
587 ierr = 0
588 return
590
592
593!-------------------------------------------------------------------------------
594!
595! Model create routine (REQUIRED)
596!
597!-------------------------------------------------------------------------------
598recursive subroutine model_create_routine( &
599 model_create_handle, requested_length_unit, requested_energy_unit, &
600 requested_charge_unit, requested_temperature_unit, requested_time_unit, &
601 ierr) bind(c)
602 use, intrinsic :: iso_c_binding
605 implicit none
606
607 !-- Transferred variables
608 type(kim_model_create_handle_type), intent(inout) :: model_create_handle
609 type(kim_length_unit_type), intent(in), value :: requested_length_unit
610 type(kim_energy_unit_type), intent(in), value :: requested_energy_unit
611 type(kim_charge_unit_type), intent(in), value :: requested_charge_unit
612 type(kim_temperature_unit_type), intent(in), value :: &
613 requested_temperature_unit
614 type(kim_time_unit_type), intent(in), value :: requested_time_unit
615 integer(c_int), intent(out) :: ierr
616
617 !-- KIM variables
618 integer(c_int) :: ierr2
619 type(buffer_type), pointer :: buf
620
621 ierr = 0
622 ierr2 = 0
623
624 ! avoid unsed dummy argument warnings
625 if (requested_length_unit == kim_length_unit_unused) continue
626 if (requested_energy_unit == kim_energy_unit_unused) continue
627 if (requested_charge_unit == kim_charge_unit_unused) continue
628 if (requested_temperature_unit == kim_temperature_unit_unused) continue
629 if (requested_time_unit == kim_time_unit_unused) continue
630
631 ! set units
632 call kim_set_units(model_create_handle, &
633 kim_length_unit_a, &
634 kim_energy_unit_ev, &
635 kim_charge_unit_unused, &
636 kim_temperature_unit_unused, &
637 kim_time_unit_unused, &
638 ierr2)
639 ierr = ierr + ierr2
640
641 ! register species
642 call kim_set_species_code(model_create_handle, &
643 kim_species_name_ar, speccode, ierr2)
644 ierr = ierr + ierr2
645
646 ! register numbering
647 call kim_set_model_numbering(model_create_handle, &
648 kim_numbering_one_based, ierr2)
649 ierr = ierr + ierr2
650
651 ! register function pointers
652 call kim_set_routine_pointer( &
653 model_create_handle, &
654 kim_model_routine_name_compute, kim_language_name_fortran, &
655 1, c_funloc(compute_energy_forces), ierr2)
656 ierr = ierr + ierr2
657 call kim_set_routine_pointer( &
658 model_create_handle, &
659 kim_model_routine_name_compute_arguments_create, &
660 kim_language_name_fortran, &
661 1, c_funloc(model_compute_arguments_create), ierr2)
662 ierr = ierr + ierr2
663 call kim_set_routine_pointer( &
664 model_create_handle, &
665 kim_model_routine_name_compute_arguments_destroy, &
666 kim_language_name_fortran, &
667 1, c_funloc(model_compute_arguments_destroy), ierr2)
668 ierr = ierr + ierr2
669 call kim_set_routine_pointer( &
670 model_create_handle, &
671 kim_model_routine_name_destroy, kim_language_name_fortran, 1, &
672 c_funloc(model_destroy_func), ierr2)
673 ierr = ierr + ierr2
674
675 ! allocate buffer
676 allocate (buf)
677
678 ! store model buffer in KIM object
679 call kim_set_model_buffer_pointer(model_create_handle, &
680 c_loc(buf))
681
682 ! set buffer values
683 buf%influence_distance = model_cutoff1 + model_cutoff2
684 buf%cutoff(1) = model_cutoff1
685 buf%cutoff(2) = model_cutoff2
686 buf%model_will_not_request_neighbors_of_noncontributing_particles(1) = 0
687 buf%model_will_not_request_neighbors_of_noncontributing_particles(2) = 1
688
689 ! register influence distance
690 call kim_set_influence_distance_pointer( &
691 model_create_handle, buf%influence_distance)
692
693 ! register cutoff
694 call kim_set_neighbor_list_pointers( &
695 model_create_handle, 2, buf%cutoff, &
696 buf%model_will_not_request_neighbors_of_noncontributing_particles)
697
698 if (ierr /= 0) then
699 ierr = 1
700 deallocate (buf)
701 call kim_log_entry( &
702 model_create_handle, kim_log_verbosity_error, &
703 "Unable to successfully initialize model")
704 end if
705
706 return
707
708end subroutine model_create_routine
static void calc_phi(double const *epsilon, double const *C, double const *Rzero, double const *shift, double const cutoff, double const r, double *phi)
integer(c_int), parameter, public speccode
recursive subroutine, public model_destroy_func(model_destroy_handle, ierr)
recursive subroutine, public model_compute_arguments_create(model_compute_handle, model_compute_arguments_create_handle, ierr)
recursive subroutine, public compute_energy_forces(model_compute_handle, model_compute_arguments_handle, ierr)
real(c_double), parameter, public model_cutoff1
real(c_double), parameter, public model_cutoff2
recursive subroutine, public model_compute_arguments_destroy(model_compute_handle, model_compute_arguments_destroy_handle, ierr)