47 use,
intrinsic :: iso_c_binding
63 integer(c_int),
parameter :: cd = c_double
64 integer(c_int),
parameter :: DIM = 3
65 integer(c_int),
parameter :: speccode = 1
66 real(c_double),
parameter :: model_cutoff = 8.15_cd
68 real(c_double),
parameter :: model_cutsq = model_cutoff**2
85 real(c_double),
parameter :: lj_epsilon = 0.0104_cd
86 real(c_double),
parameter :: lj_sigma = 3.40_cd
87 real(c_double),
parameter :: lj_cutnorm = model_cutoff/lj_sigma
88 real(c_double),
parameter :: lj_A = 12.0_cd*lj_epsilon*(-26.0_cd &
89 + 7.0_cd*lj_cutnorm**6)/(lj_cutnorm**14 &
91 real(c_double),
parameter :: lj_B = 96.0_cd*lj_epsilon*(7.0_cd &
92 - 2.0_cd*lj_cutnorm**6)/(lj_cutnorm**13*lj_sigma)
93 real(c_double),
parameter :: lj_C = 28.0_cd*lj_epsilon*(-13.0_cd &
94 + 4.0_cd*lj_cutnorm**6)/(lj_cutnorm**12)
96 type, bind(c) :: buffer_type
97 real(c_double) :: influence_distance
98 real(c_double) :: cutoff(1)
108 subroutine calc_phi(r,phi)
112 real(c_double),
intent(in) :: r
113 real(c_double),
intent(out) :: phi
116 real(c_double) rsq,sor,sor6,sor12
123 if (r .gt. model_cutoff)
then 127 phi = 4.0_cd*lj_epsilon*(sor12-sor6) + lj_a*rsq + lj_b*r + lj_c
130 end subroutine calc_phi
137 subroutine calc_phi_dphi(r,phi,dphi)
141 real(c_double),
intent(in) :: r
142 real(c_double),
intent(out) :: phi,dphi
145 real(c_double) rsq,sor,sor6,sor12
152 if (r .gt. model_cutoff)
then 157 phi = 4.0_cd*lj_epsilon*(sor12-sor6) + lj_a*rsq + lj_b*r + lj_c
158 dphi = 24.0_cd*lj_epsilon*(-2.0_cd*sor12+sor6)/r + 2.0_cd*lj_a*r + lj_b
161 end subroutine calc_phi_dphi
168 #include "kim_model_compute_log_macros.fd" 170 model_compute_arguments_handle, ierr) bind(c)
174 type(kim_model_compute_handle_type),
intent(in) :: model_compute_handle
175 type(kim_model_compute_arguments_handle_type),
intent(in) :: &
176 model_compute_arguments_handle
177 integer(c_int),
intent(out) :: ierr
180 real(c_double) :: rij(dim)
181 real(c_double) :: r,rsqij,phi,dphi,deidr = 0.0_cd
182 integer(c_int) :: i,j,jj,numnei,comp_force,comp_enepot,comp_virial, comp_energy
183 integer(c_int) :: ierr2
186 integer(c_int),
pointer :: n
187 real(c_double),
pointer :: energy
188 real(c_double),
pointer :: coor(:,:)
189 real(c_double),
pointer :: force(:,:)
190 real(c_double),
pointer :: enepot(:)
191 integer(c_int),
pointer :: nei1part(:)
192 integer(c_int),
pointer :: particlespeciescodes(:)
193 integer(c_int),
pointer :: particlecontributing(:)
194 real(c_double),
pointer :: virial(:)
196 kim_log_file = __file__
201 call kim_model_compute_arguments_get_argument_pointer( &
202 model_compute_arguments_handle, &
203 kim_compute_argument_name_number_of_particles, n, ierr2)
205 call kim_model_compute_arguments_get_argument_pointer( &
206 model_compute_arguments_handle, &
207 kim_compute_argument_name_particle_species_codes, n, particlespeciescodes, &
210 call kim_model_compute_arguments_get_argument_pointer( &
211 model_compute_arguments_handle, &
212 kim_compute_argument_name_particle_contributing, n, particlecontributing, &
215 call kim_model_compute_arguments_get_argument_pointer( &
216 model_compute_arguments_handle, &
217 kim_compute_argument_name_coordinates, dim, n, coor, ierr2)
219 call kim_model_compute_arguments_get_argument_pointer( &
220 model_compute_arguments_handle, &
221 kim_compute_argument_name_partial_energy, energy, ierr2)
223 call kim_model_compute_arguments_get_argument_pointer( &
224 model_compute_arguments_handle, &
225 kim_compute_argument_name_partial_forces, dim, n, force, ierr2)
227 call kim_model_compute_arguments_get_argument_pointer( &
228 model_compute_arguments_handle, &
229 kim_compute_argument_name_partial_particle_energy, n, enepot, ierr2)
231 call kim_model_compute_arguments_get_argument_pointer( &
232 model_compute_arguments_handle, &
233 kim_compute_argument_name_partial_virial, 6, virial, ierr2)
236 kim_log_message =
"get data" 244 if (
associated(energy))
then 249 if (
associated(force))
then 254 if (
associated(enepot))
then 259 if (
associated(virial))
then 269 if (particlespeciescodes(i).ne.speccode)
then 270 kim_log_message =
"Unexpected species code detected" 279 if (comp_enepot.eq.1) enepot = 0.0_cd
280 if (comp_energy.eq.1) energy = 0.0_cd
281 if (comp_force.eq.1) force = 0.0_cd
282 if (comp_virial.eq.1) virial = 0.0_cd
291 if (particlecontributing(i) == 1)
then 293 call kim_model_compute_arguments_get_neighbor_list( &
294 model_compute_arguments_handle, 1, i, numnei, nei1part, ierr)
297 kim_log_message =
"GetNeighborList failed" 311 rij(:) = coor(:,j) - coor(:,i)
315 rsqij = dot_product(rij,rij)
316 if ( rsqij .lt. model_cutsq )
then 319 if (comp_force.eq.1.or.comp_virial.eq.1)
then 320 call calc_phi_dphi(r,phi,dphi)
329 if (comp_enepot.eq.1)
then 330 enepot(i) = enepot(i) + 0.5_cd*phi
332 if (comp_energy.eq.1)
then 333 energy = energy + 0.5_cd*phi
338 if (comp_virial.eq.1)
then 339 virial(1) = virial(1) + rij(1)*rij(1)*deidr/r
340 virial(2) = virial(2) + rij(2)*rij(2)*deidr/r
341 virial(3) = virial(3) + rij(3)*rij(3)*deidr/r
342 virial(4) = virial(4) + rij(2)*rij(3)*deidr/r
343 virial(5) = virial(5) + rij(1)*rij(3)*deidr/r
344 virial(6) = virial(6) + rij(1)*rij(2)*deidr/r
349 if (comp_force.eq.1)
then 350 force(:,i) = force(:,i) + deidr*rij/r
351 force(:,j) = force(:,j) - deidr*rij/r
374 #include "kim_model_destroy_log_macros.fd" 376 use,
intrinsic :: iso_c_binding
380 type(kim_model_destroy_handle_type),
intent(inout) :: model_destroy_handle
381 integer(c_int),
intent(out) :: ierr
383 type(buffer_type),
pointer :: buf;
type(c_ptr) :: pbuf
385 kim_log_file = __file__
387 call kim_model_destroy_get_model_buffer_pointer(model_destroy_handle, pbuf)
388 call c_f_pointer(pbuf, buf)
389 kim_log_message =
"deallocating model buffer" 400 #include "kim_model_refresh_log_macros.fd" 402 use,
intrinsic :: iso_c_binding
406 type(kim_model_refresh_handle_type),
intent(inout) :: model_refresh_handle
407 integer(c_int),
intent(out) :: ierr
409 type(buffer_type),
pointer :: buf;
type(c_ptr) :: pbuf
411 kim_log_file = __file__
413 call kim_model_refresh_get_model_buffer_pointer(model_refresh_handle, pbuf)
414 call c_f_pointer(pbuf, buf)
416 kim_log_message =
"Resettings influence distance and cutoffs" 418 call kim_model_refresh_set_influence_distance_pointer( &
419 model_refresh_handle, buf%cutoff(1))
420 call kim_model_refresh_set_neighbor_list_cutoffs_pointer( &
421 model_refresh_handle, 1, buf%cutoff)
431 #include "kim_model_compute_arguments_create_log_macros.fd" 433 model_compute_arguments_create_handle, ierr) bind(c)
434 use,
intrinsic :: iso_c_binding
436 log_entry=>kim_model_compute_arguments_create_log_entry
440 type(kim_model_compute_handle_type),
intent(in) :: model_compute_handle
441 type(kim_model_compute_arguments_create_handle_type),
intent(inout) :: &
442 model_compute_arguments_create_handle
443 integer(c_int),
intent(out) :: ierr
445 integer(c_int) :: ierr2
451 call kim_model_compute_arguments_create_set_argument_support_status( &
452 model_compute_arguments_create_handle, &
453 kim_compute_argument_name_partial_energy, &
454 kim_support_status_optional, ierr2)
456 call kim_model_compute_arguments_create_set_argument_support_status( &
457 model_compute_arguments_create_handle, &
458 kim_compute_argument_name_partial_forces, &
459 kim_support_status_optional, ierr2)
461 call kim_model_compute_arguments_create_set_argument_support_status( &
462 model_compute_arguments_create_handle, &
463 kim_compute_argument_name_partial_particle_energy, &
464 kim_support_status_optional, ierr2)
466 call kim_model_compute_arguments_create_set_argument_support_status( &
467 model_compute_arguments_create_handle, &
468 kim_compute_argument_name_partial_virial, &
469 kim_support_status_optional, ierr2)
477 kim_log_message =
"Unable to successfully create compute_arguments object" 489 #include "kim_model_compute_arguments_destroy_log_macros.fd" 491 model_compute_arguments_destroy_handle, ierr) bind(c)
492 use,
intrinsic :: iso_c_binding
494 log_entry=>kim_model_compute_arguments_destroy_log_entry
498 type(kim_model_compute_handle_type),
intent(in) :: model_compute_handle
499 type(kim_model_compute_arguments_destroy_handle_type),
intent(inout) :: &
500 model_compute_arguments_destroy_handle
501 integer(c_int),
intent(out) :: ierr
503 integer(c_int) :: ierr2
520 #include "kim_model_create_log_macros.fd" 521 subroutine model_create_routine(model_create_handle, requested_length_unit, &
522 requested_energy_unit, requested_charge_unit, requested_temperature_unit, &
523 requested_time_unit, ierr) bind(c)
524 use,
intrinsic :: iso_c_binding
530 type(kim_model_create_handle_type),
intent(inout) :: model_create_handle
531 type(kim_length_unit_type),
intent(in) :: requested_length_unit
532 type(kim_energy_unit_type),
intent(in) :: requested_energy_unit
533 type(kim_charge_unit_type),
intent(in) :: requested_charge_unit
534 type(kim_temperature_unit_type),
intent(in) :: requested_temperature_unit
535 type(kim_time_unit_type),
intent(in) :: requested_time_unit
536 integer(c_int),
intent(out) :: ierr
539 integer(c_int) :: ierr2
540 type(buffer_type),
pointer :: buf
542 kim_log_file = __file__
548 call kim_model_create_set_units(model_create_handle, &
550 kim_energy_unit_ev, &
551 kim_charge_unit_unused, &
552 kim_temperature_unit_unused, &
553 kim_time_unit_unused, &
558 call kim_model_create_set_species_code(model_create_handle, &
559 kim_species_name_ar, speccode, ierr2)
563 call kim_model_create_set_model_numbering(model_create_handle, &
564 kim_numbering_one_based, ierr2);
568 call kim_model_create_set_compute_pointer(model_create_handle, &
571 call kim_model_create_set_compute_arguments_create_pointer( &
572 model_create_handle, kim_language_name_fortran, &
575 call kim_model_create_set_compute_arguments_destroy_pointer( &
576 model_create_handle, kim_language_name_fortran, &
579 call kim_model_create_set_destroy_pointer(model_create_handle, &
582 call kim_model_create_set_refresh_pointer( &
583 model_create_handle, kim_language_name_fortran, &
591 call kim_model_create_set_model_buffer_pointer(model_create_handle, &
595 buf%influence_distance = model_cutoff
596 buf%cutoff = model_cutoff
599 call kim_model_create_set_influence_distance_pointer( &
600 model_create_handle, buf%influence_distance)
603 call kim_model_create_set_neighbor_list_cutoffs_pointer(model_create_handle, &
609 kim_log_message =
"Unable to successfully initialize model" 615 end subroutine model_create_routine
subroutine, public model_destroy_func(model_destroy_handle, ierr)
subroutine, public model_refresh_func(model_refresh_handle, ierr)
subroutine, public model_compute_arguments_destroy(model_compute_handle, model_compute_arguments_destroy_handle, ierr)
subroutine, public model_compute_arguments_create(model_compute_handle, model_compute_arguments_create_handle, ierr)
integer(c_int) function, public compute_energy_forces(pkim)