53 use,
intrinsic :: iso_c_binding
69 integer(c_int),
parameter :: cd = c_double
70 integer(c_int),
parameter :: dim = 3
88 real(c_double),
parameter :: lj_spring = 0.00051226_cd
89 real(c_double),
parameter :: lj_sigma = 5.26_cd/1.74724_cd
98 type, bind(c) :: buffer_type
99 real(c_double) :: influence_distance
100 real(c_double) :: cutoff(2)
102 model_will_not_request_neighbors_of_noncontributing_particles(2)
113 recursive subroutine calc_phi(r,phi,dphi,calc_deriv)
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
123 real(c_double) rsq,sor,sor6,sor12
133 if (calc_deriv) dphi = 0.0_cd
135 phi = 4.0_cd*(sor12-sor6)
136 if (calc_deriv) dphi = 24.0_cd*(-2.0_cd*sor12+sor6)/r
146 recursive subroutine calc_spring_energyamp(model_compute_arguments_handle, &
147 atom, coor, eps, ierr)
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
160 real(c_double) rrel(dim)
161 real(c_double) rsqrel
162 integer(c_int) numneishort
163 integer(c_int),
pointer :: neishort(:)
166 call kim_get_neighbor_list( &
167 model_compute_arguments_handle, 1, atom, numneishort, neishort, ierr)
168 if (ierr /= 0)
return 171 do kk = 1, numneishort
173 rrel(:) = coor(:,k) - coor(:,atom)
174 rsqrel = dot_product(rrel,rrel)
175 if (rsqrel .lt. model_cutsq1) eps = eps + rsqrel
177 eps = 0.5_cd*lj_spring*eps
179 end subroutine calc_spring_energyamp
186 recursive subroutine calc_spring_force(model_compute_arguments_handle, atom, &
187 coor, eps, phi, force, ierr)
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
202 real(c_double) rrel(dim), dforce(dim)
203 real(c_double) rsqrel
204 integer(c_int) numneishort
205 integer(c_int),
pointer :: neishort(:)
208 call kim_get_neighbor_list( &
209 model_compute_arguments_handle, 1, atom, numneishort, neishort, ierr)
210 if (ierr /= 0)
return 214 do kk = 1, numneishort
216 rrel(:) = coor(:,k) - coor(:,atom)
217 rsqrel = dot_product(rrel,rrel)
218 if (rsqrel .lt. model_cutsq1)
then 219 dforce(:) = 0.5_cd*eps*lj_spring*rrel(:)*phi
220 force(:,atom) = force(:,atom) + dforce(:)
221 force(:,k) = force(:,k) - dforce(:)
225 end subroutine calc_spring_force
234 model_compute_arguments_handle, ierr) bind(c)
238 type(kim_model_compute_handle_type),
intent(in) :: model_compute_handle
239 type(kim_model_compute_arguments_handle_type),
intent(in) :: &
240 model_compute_arguments_handle
241 integer(c_int),
intent(out) :: ierr
244 real(c_double) :: Rij(dim),dforce(dim)
245 real(c_double) :: r,Rsqij,phi,dphi,dEidr,epsi,epsj
246 integer(c_int) :: i,j,jj,comp_force,comp_enepot,comp_energy
248 integer(c_int) :: numnei
249 integer(c_int) :: ierr2
250 logical :: calc_deriv
253 integer(c_int),
pointer :: N
254 real(c_double),
pointer :: energy
255 real(c_double),
pointer :: coor(:,:)
256 real(c_double),
pointer :: force(:,:)
257 real(c_double),
pointer :: enepot(:)
258 integer(c_int),
pointer :: nei1part(:)
259 integer(c_int),
pointer :: particleSpeciesCodes(:)
260 integer(c_int),
pointer :: particleContributing(:)
266 call kim_get_argument_pointer( &
267 model_compute_arguments_handle, &
268 kim_compute_argument_name_number_of_particles, n, ierr2)
270 call kim_get_argument_pointer( &
271 model_compute_arguments_handle, &
272 kim_compute_argument_name_particle_species_codes, n, particlespeciescodes, &
275 call kim_get_argument_pointer( &
276 model_compute_arguments_handle, &
277 kim_compute_argument_name_particle_contributing, n, particlecontributing, &
280 call kim_get_argument_pointer( &
281 model_compute_arguments_handle, &
282 kim_compute_argument_name_coordinates, dim, n, coor, ierr2)
284 call kim_get_argument_pointer( &
285 model_compute_arguments_handle, &
286 kim_compute_argument_name_partial_energy, energy, ierr2)
288 call kim_get_argument_pointer( &
289 model_compute_arguments_handle, &
290 kim_compute_argument_name_partial_forces, dim, n, force, ierr2)
292 call kim_get_argument_pointer( &
293 model_compute_arguments_handle, &
294 kim_compute_argument_name_partial_particle_energy, n, enepot, ierr2)
301 call kim_log_entry(model_compute_arguments_handle, &
302 kim_log_verbosity_error,
"get data")
309 if (
associated(energy))
then 314 if (
associated(force))
then 319 if (
associated(enepot))
then 329 calc_deriv = comp_force.eq.1
335 if (particlespeciescodes(i).ne.
speccode)
then 336 call kim_log_entry(model_compute_handle, &
337 kim_log_verbosity_error,
"Unexpected species code detected")
345 if (comp_enepot.eq.1) enepot = 0.0_cd
346 if (comp_energy.eq.1) energy = 0.0_cd
347 if (comp_force.eq.1) force = 0.0_cd
349 if (calc_deriv) deidr = 0.0_cd
358 if (particlecontributing(i) == 1)
then 360 call kim_get_neighbor_list( &
361 model_compute_arguments_handle, 2, i, numnei, nei1part, ierr)
364 call kim_log_entry( &
365 model_compute_arguments_handle, kim_log_verbosity_error, &
366 "GetNeighborList failed")
372 call calc_spring_energyamp(model_compute_arguments_handle, i, coor, epsi, ierr)
375 call kim_log_entry( &
376 model_compute_handle, kim_log_verbosity_error, &
377 "GetNeighborList failed")
389 call calc_spring_energyamp(model_compute_arguments_handle, j, coor, epsj, ierr)
392 call kim_log_entry( &
393 model_compute_handle, kim_log_verbosity_error, &
394 "GetNeighborList failed")
401 rij(:) = coor(:,j) - coor(:,i)
405 rsqij = dot_product(rij,rij)
406 if ( rsqij .lt. model_cutsq2 )
then 409 call calc_phi(r,phi,dphi,calc_deriv)
410 if (calc_deriv) deidr = 0.5_cd*dphi
414 if (comp_enepot.eq.1)
then 415 enepot(i) = enepot(i) + 0.5_cd*epsi*epsj*phi
417 if (comp_energy.eq.1)
then 418 energy = energy + 0.5_cd*epsi*epsj*phi
437 if (comp_force.eq.1)
then 439 call calc_spring_force(model_compute_arguments_handle, i, coor, epsj, &
443 call kim_log_entry( &
444 model_compute_handle, kim_log_verbosity_error, &
445 "GetNeighborList failed")
450 call calc_spring_force(model_compute_arguments_handle, j, coor, epsi, &
454 call kim_log_entry( &
455 model_compute_handle, kim_log_verbosity_error, &
456 "GetNeighborList failed")
461 dforce(:) = epsi*epsj*deidr*rij(:)/r
462 force(:,i) = force(:,i) + dforce(:)
463 force(:,j) = force(:,j) - dforce(:)
487 use,
intrinsic :: iso_c_binding
491 type(kim_model_destroy_handle_type),
intent(inout) :: model_destroy_handle
492 integer(c_int),
intent(out) :: ierr
494 type(buffer_type),
pointer :: buf;
type(c_ptr) :: pbuf
496 call kim_get_model_buffer_pointer(model_destroy_handle, pbuf)
497 call c_f_pointer(pbuf, buf)
498 call kim_log_entry(model_destroy_handle, &
499 kim_log_verbosity_information,
"deallocating model buffer")
510 model_compute_arguments_create_handle, ierr) bind(c)
511 use,
intrinsic :: iso_c_binding
515 type(kim_model_compute_handle_type),
intent(in) :: model_compute_handle
516 type(kim_model_compute_arguments_create_handle_type),
intent(inout) :: &
517 model_compute_arguments_create_handle
518 integer(c_int),
intent(out) :: ierr
520 integer(c_int) :: ierr2
523 if (model_compute_handle .eq. kim_model_compute_null_handle)
continue 529 call kim_set_argument_support_status( &
530 model_compute_arguments_create_handle, &
531 kim_compute_argument_name_partial_energy, &
532 kim_support_status_optional, ierr2)
534 call kim_set_argument_support_status( &
535 model_compute_arguments_create_handle, &
536 kim_compute_argument_name_partial_forces, &
537 kim_support_status_optional, ierr2)
539 call kim_set_argument_support_status( &
540 model_compute_arguments_create_handle, &
541 kim_compute_argument_name_partial_particle_energy, &
542 kim_support_status_optional, ierr2)
555 call kim_log_entry( &
556 model_compute_arguments_create_handle, &
557 kim_log_verbosity_error, &
558 "Unable to successfully create compute_arguments object")
570 model_compute_arguments_destroy_handle, ierr) bind(c)
571 use,
intrinsic :: iso_c_binding
575 type(kim_model_compute_handle_type),
intent(in) :: model_compute_handle
576 type(kim_model_compute_arguments_destroy_handle_type),
intent(inout) :: &
577 model_compute_arguments_destroy_handle
578 integer(c_int),
intent(out) :: ierr
581 if (model_compute_handle .eq. kim_model_compute_null_handle)
continue 582 if (model_compute_arguments_destroy_handle .eq. &
583 kim_model_compute_arguments_destroy_null_handle)
continue 596 recursive subroutine model_create_routine(model_create_handle, &
597 requested_length_unit, requested_energy_unit, requested_charge_unit, &
598 requested_temperature_unit, requested_time_unit, ierr) bind(c)
599 use,
intrinsic :: iso_c_binding
605 type(kim_model_create_handle_type),
intent(inout) :: model_create_handle
606 type(kim_length_unit_type),
intent(in),
value :: requested_length_unit
607 type(kim_energy_unit_type),
intent(in),
value :: requested_energy_unit
608 type(kim_charge_unit_type),
intent(in),
value :: requested_charge_unit
609 type(kim_temperature_unit_type),
intent(in),
value :: requested_temperature_unit
610 type(kim_time_unit_type),
intent(in),
value :: requested_time_unit
611 integer(c_int),
intent(out) :: ierr
614 integer(c_int) :: ierr2
615 type(buffer_type),
pointer :: buf
621 if (requested_length_unit .eq. kim_length_unit_unused)
continue 622 if (requested_energy_unit .eq. kim_energy_unit_unused)
continue 623 if (requested_charge_unit .eq. kim_charge_unit_unused)
continue 624 if (requested_temperature_unit .eq. kim_temperature_unit_unused)
continue 625 if (requested_time_unit .eq. kim_time_unit_unused)
continue 628 call kim_set_units(model_create_handle, &
630 kim_energy_unit_ev, &
631 kim_charge_unit_unused, &
632 kim_temperature_unit_unused, &
633 kim_time_unit_unused, &
638 call kim_set_species_code(model_create_handle, &
639 kim_species_name_ar,
speccode, ierr2)
643 call kim_set_model_numbering(model_create_handle, &
644 kim_numbering_one_based, ierr2);
648 call kim_set_routine_pointer(model_create_handle, &
649 kim_model_routine_name_compute, kim_language_name_fortran, &
652 call kim_set_routine_pointer(model_create_handle, &
653 kim_model_routine_name_compute_arguments_create, kim_language_name_fortran, &
656 call kim_set_routine_pointer(model_create_handle, &
657 kim_model_routine_name_compute_arguments_destroy, kim_language_name_fortran, &
660 call kim_set_routine_pointer(model_create_handle, &
661 kim_model_routine_name_destroy, kim_language_name_fortran, 1, &
669 call kim_set_model_buffer_pointer(model_create_handle, &
676 buf%model_will_not_request_neighbors_of_noncontributing_particles(1) = 0
677 buf%model_will_not_request_neighbors_of_noncontributing_particles(2) = 1
680 call kim_set_influence_distance_pointer( &
681 model_create_handle, buf%influence_distance)
684 call kim_set_neighbor_list_pointers(model_create_handle, &
686 buf%model_will_not_request_neighbors_of_noncontributing_particles)
691 call kim_log_entry(model_create_handle, &
692 kim_log_verbosity_error,
"Unable to successfully initialize model")
697 end subroutine model_create_routine
recursive subroutine, public model_compute_arguments_create(model_compute_handle, model_compute_arguments_create_handle, ierr)
recursive subroutine, public model_compute_arguments_destroy(model_compute_handle, model_compute_arguments_destroy_handle, ierr)
real(c_double), parameter, public model_cutoff1
real(c_double), parameter, public model_cutoff2
integer(c_int), parameter, public speccode
static void calc_phi(double const *epsilon, double const *C, double const *Rzero, double const *shift, double const cutoff, double const r, double *phi)
recursive subroutine, public compute_energy_forces(model_compute_handle, model_compute_arguments_handle, ierr)
recursive subroutine, public model_destroy_func(model_destroy_handle, ierr)