From 95841b0efd4bb2904a0486935591ad18e28dbc88 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Mon, 21 Nov 2022 22:38:10 -0600 Subject: [PATCH 01/19] Implementation (after several failures) of set_fix_external_callback --- fortran/lammps.f90 | 223 ++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 220 insertions(+), 3 deletions(-) diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index 3ab7a26d25..02664c54db 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -35,7 +35,7 @@ MODULE LIBLAMMPS USE, INTRINSIC :: ISO_C_BINDING, ONLY: c_ptr, c_null_ptr, C_ASSOCIATED, & C_LOC, c_int, c_int64_t, c_char, c_null_char, c_double, c_size_t, & - C_F_POINTER + C_F_POINTER, c_funptr, C_FUNLOC IMPLICIT NONE PRIVATE @@ -171,7 +171,7 @@ MODULE LIBLAMMPS PROCEDURE, PRIVATE :: lmp_decode_image_flags_bigbig GENERIC :: decode_image_flags => lmp_decode_image_flags, & lmp_decode_image_flags_bigbig -! + PROCEDURE :: set_fix_external_callback => lmp_set_fix_external_callback PROCEDURE :: flush_buffers => lmp_flush_buffers PROCEDURE :: is_running => lmp_is_running PROCEDURE :: force_timeout => lmp_force_timeout @@ -261,6 +261,50 @@ MODULE LIBLAMMPS assign_int64_to_lammps_image_data END INTERFACE + ! Interface templates for fix external callbacks + ABSTRACT INTERFACE + SUBROUTINE external_callback_smallsmall(caller, timestep, ids, x, fexternal) + IMPORT :: c_int, c_double + CLASS(*), INTENT(IN) :: caller + INTEGER(c_int), INTENT(IN) :: timestep + INTEGER(c_int), DIMENSION(:), INTENT(IN) :: ids + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: fexternal + END SUBROUTINE external_callback_smallsmall + SUBROUTINE external_callback_smallbig(caller, timestep, ids, x, fexternal) + IMPORT :: c_int, c_double, c_int64_t + CLASS(*), INTENT(IN) :: caller + INTEGER(c_int64_t), INTENT(IN) :: timestep + INTEGER(c_int), DIMENSION(:), INTENT(IN) :: ids + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: fexternal + END SUBROUTINE external_callback_smallbig + SUBROUTINE external_callback_bigbig(caller, timestep, ids, x, fexternal) + IMPORT :: c_double, c_int64_t + CLASS(*), INTENT(IN) :: caller + INTEGER(c_int64_t), INTENT(IN) :: timestep + INTEGER(c_int64_t), DIMENSION(:), INTENT(IN) :: ids + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: fexternal + END SUBROUTINE external_callback_bigbig + END INTERFACE + + ! Derived type for fix external callback data + TYPE fix_external_data + CHARACTER(LEN=:), ALLOCATABLE :: id + PROCEDURE(external_callback_smallsmall), NOPASS, POINTER :: & + callback_smallsmall => NULL() + PROCEDURE(external_callback_smallbig), NOPASS, POINTER :: & + callback_smallbig => NULL() + PROCEDURE(external_callback_bigbig), NOPASS, POINTER :: & + callback_bigbig => NULL() + CLASS(*), POINTER :: caller => NULL() + CLASS(lammps), POINTER :: lammps_instance => NULL() + END TYPE fix_external_data + + ! Array used to store Fortran-facing callback functions for fix external + TYPE(fix_external_data), DIMENSION(:), ALLOCATABLE, TARGET, SAVE :: ext_data + ! interface definitions for calling functions in library.cpp INTERFACE FUNCTION lammps_open(argc, argv, comm) BIND(C,name='lammps_open_fortran') @@ -705,7 +749,12 @@ MODULE LIBLAMMPS ! It is re-written in Fortran below. It was easier to do the same for ! lammps_decode_image_flags's equivalent. - !SUBROUTINE lammps_set_fix_external_callback ! may have trouble.... + SUBROUTINE lammps_set_fix_external_callback(handle, id, funcptr, ptr) & + BIND(C) + IMPORT :: c_ptr, c_funptr + TYPE(c_ptr), VALUE :: handle, id, ptr + TYPE(c_funptr), VALUE :: funcptr + END SUBROUTINE lammps_set_fix_external_callback !FUNCTION lammps_fix_external_get_force() ! returns real(c_double)(:) !SUBROUTINE lammps_fix_external_set_energy_global @@ -2326,6 +2375,174 @@ CONTAINS END IF END SUBROUTINE lmp_decode_image_flags_bigbig + ! equivalent function to lammps_set_fix_external_callback for -DSMALLSMALL + ! note that "caller" is wrapped into a fix_external_data derived type along + ! with the fix id and the Fortran calling function. + SUBROUTINE lmp_set_fix_external_callback(self, id, callback, caller) + CLASS(lammps), TARGET, INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + EXTERNAL :: callback + CLASS(*), INTENT(IN), TARGET, OPTIONAL :: caller + TYPE(c_ptr) :: c_id, c_caller + TYPE(c_funptr) :: c_callback + INTEGER :: i, this_fix, size_tagint, size_bigint + + size_tagint = lmp_extract_setting(self, 'tagint') + size_bigint = lmp_extract_setting(self, 'bigint') + + c_id = f2c_string(id) + IF (ALLOCATED(ext_data)) THEN + this_fix = SIZE(ext_data) + 1 + DO i = 1, SIZE(ext_data) + IF (ext_data(i)%id == id) THEN + this_fix = i + EXIT + END IF + END DO + IF (this_fix > SIZE(ext_data)) THEN + ext_data = [ext_data, fix_external_data()] ! extends ext_data by 1 + END IF + ELSE + ALLOCATE(ext_data(1)) + this_fix = 1 + END IF + ext_data(this_fix)%id = id + ext_data(this_fix)%lammps_instance => self + + IF (size_tagint == 4_c_int .AND. size_bigint == 4_c_int) THEN + ! -DSMALLSMALL + c_callback = C_FUNLOC(callback_wrapper_smallsmall) + CALL set_fix_external_callback_smallsmall(this_fix, callback) + ELSE IF (size_tagint == 8_c_int .AND. size_bigint == 8_c_int) THEN + ! -DBIGBIG + c_callback = C_FUNLOC(callback_wrapper_bigbig) + CALL set_fix_external_callback_bigbig(this_fix, callback) + ELSE + ! -DSMALLBIG + c_callback = C_FUNLOC(callback_wrapper_smallbig) + CALL set_fix_external_callback_smallbig(this_fix, callback) + END IF + + IF (PRESENT(caller)) THEN + ext_data(this_fix)%caller => caller + ELSE + NULLIFY(ext_data(this_fix)%caller) + END IF + c_caller = C_LOC(ext_data(this_fix)) + CALL lammps_set_fix_external_callback(self%handle, c_id, c_callback, & + c_caller) + END SUBROUTINE lmp_set_fix_external_callback + + ! Wrappers to assign callback pointers with explicit interfaces + SUBROUTINE set_fix_external_callback_smallsmall(id, callback) + INTEGER, INTENT(IN) :: id + PROCEDURE(external_callback_smallsmall) :: callback + + ext_data(id)%callback_smallsmall => callback + END SUBROUTINE set_fix_external_callback_smallsmall + + SUBROUTINE set_fix_external_callback_smallbig(id, callback) + INTEGER, INTENT(IN) :: id + PROCEDURE(external_callback_smallbig) :: callback + + ext_data(id)%callback_smallbig => callback + END SUBROUTINE set_fix_external_callback_smallbig + + SUBROUTINE set_fix_external_callback_bigbig(id, callback) + INTEGER, INTENT(IN) :: id + PROCEDURE(external_callback_bigbig) :: callback + + ext_data(id)%callback_bigbig => callback + END SUBROUTINE set_fix_external_callback_bigbig + + ! companions to lmp_set_fix_external_callback to change interface + SUBROUTINE callback_wrapper_smallsmall(caller, timestep, nlocal, ids, x, & + fexternal) BIND(C) + TYPE(c_ptr), INTENT(IN), VALUE :: caller + INTEGER(c_int), INTENT(IN), VALUE :: timestep + INTEGER(c_int), INTENT(IN), VALUE :: nlocal + TYPE(c_ptr), INTENT(IN), VALUE :: ids, x, fexternal + TYPE(c_ptr), DIMENSION(:), POINTER :: x0, f0 + INTEGER(c_int), DIMENSION(:), POINTER :: f_ids => NULL() + REAL(c_double), DIMENSION(:,:), POINTER :: f_x => NULL(), & + f_fexternal => NULL() + TYPE(fix_external_data), POINTER :: f_caller => NULL() + + CALL C_F_POINTER(ids, f_ids, [nlocal]) + CALL C_F_POINTER(x, x0, [nlocal]) + CALL C_F_POINTER(x0(1), f_x, [3, nlocal]) + CALL C_F_POINTER(fexternal, f0, [nlocal]) + CALL C_F_POINTER(f0(1), f_fexternal, [3, nlocal]) + IF (C_ASSOCIATED(caller)) THEN + CALL C_F_POINTER(caller, f_caller) + CALL f_caller%callback_smallsmall(f_caller%caller, timestep, f_ids, & + f_x, f_fexternal) + ELSE + CALL lmp_error(f_caller%lammps_instance, & + LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Got null pointer from "caller"; this should never happen;& + & please report a bug') + END IF + END SUBROUTINE callback_wrapper_smallsmall + + SUBROUTINE callback_wrapper_smallbig(caller, timestep, nlocal, ids, x, & + fexternal) BIND(C) + TYPE(c_ptr), INTENT(IN), VALUE :: caller + INTEGER(c_int64_t), INTENT(IN), VALUE :: timestep + INTEGER(c_int), INTENT(IN), VALUE :: nlocal + TYPE(c_ptr), INTENT(IN), VALUE :: ids, x, fexternal + TYPE(c_ptr), DIMENSION(:), POINTER :: x0, f0 + INTEGER(c_int), DIMENSION(:), POINTER :: f_ids => NULL() + REAL(c_double), DIMENSION(:,:), POINTER :: f_x => NULL(), & + f_fexternal => NULL() + TYPE(fix_external_data), POINTER :: f_caller => NULL() + + CALL C_F_POINTER(ids, f_ids, [nlocal]) + CALL C_F_POINTER(x, x0, [nlocal]) + CALL C_F_POINTER(x0(1), f_x, [3, nlocal]) + CALL C_F_POINTER(fexternal, f0, [nlocal]) + CALL C_F_POINTER(f0(1), f_fexternal, [3, nlocal]) + IF (C_ASSOCIATED(caller)) THEN + CALL C_F_POINTER(caller, f_caller) + CALL f_caller%callback_smallbig(f_caller%caller, timestep, f_ids, f_x, & + f_fexternal) + ELSE + CALL lmp_error(f_caller%lammps_instance, & + LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Got null pointer from "caller"; this should never happen;& + & please report a bug') + END IF + END SUBROUTINE callback_wrapper_smallbig + + SUBROUTINE callback_wrapper_bigbig(caller, timestep, nlocal, ids, x, & + fexternal) BIND(C) + TYPE(c_ptr), INTENT(IN), VALUE :: caller + INTEGER(c_int64_t), INTENT(IN), VALUE :: timestep + INTEGER(c_int), INTENT(IN), VALUE :: nlocal + TYPE(c_ptr), INTENT(IN), VALUE :: ids, x, fexternal + TYPE(c_ptr), DIMENSION(:), POINTER :: x0, f0 + INTEGER(c_int64_t), DIMENSION(:), POINTER :: f_ids => NULL() + REAL(c_double), DIMENSION(:,:), POINTER :: f_x => NULL(), & + f_fexternal => NULL() + TYPE(fix_external_data), POINTER :: f_caller => NULL() + + CALL C_F_POINTER(ids, f_ids, [nlocal]) + CALL C_F_POINTER(x, x0, [nlocal]) + CALL C_F_POINTER(x0(1), f_x, [3, nlocal]) + CALL C_F_POINTER(fexternal, f0, [nlocal]) + CALL C_F_POINTER(f0(1), f_fexternal, [3, nlocal]) + IF (C_ASSOCIATED(caller)) THEN + CALL C_F_POINTER(caller, f_caller) + CALL f_caller%callback_bigbig(f_caller%caller, timestep, f_ids, f_x, & + f_fexternal) + ELSE + CALL lmp_error(f_caller%lammps_instance, & + LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Got null pointer from "caller"; this should never happen;& + & please report a bug') + END IF + END SUBROUTINE callback_wrapper_bigbig + ! equivalent function to lammps_flush_buffers SUBROUTINE lmp_flush_buffers(self) CLASS(lammps), INTENT(IN) :: self From 170c312a0c5675ae53b208b64f36bdf8e8f75a8a Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Thu, 24 Nov 2022 21:07:46 -0600 Subject: [PATCH 02/19] Fixed oversight in set_fix_external_callback and wrote its documentation --- doc/src/Fortran.rst | 63 +++++++++++++++++++++++++++++++++++++++++++++ fortran/lammps.f90 | 54 ++++++++++++++++++++++++++++---------- 2 files changed, 103 insertions(+), 14 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index ec54a594d1..12d39e8af7 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -358,6 +358,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype encode_image_flags: function :f decode_image_flags: :f:subr:`decode_image_flags` :ftype decode_image_flags: subroutine + :f set_fix_external_callback: :f:subr:`set_fix_external_callback` + :ftype set_fix_external_callback: subroutine :f flush_buffers: :f:subr:`flush_buffers` :ftype flush_buffers: subroutine :f is_running: :f:func:`is_running` @@ -2096,6 +2098,67 @@ Procedures Bound to the :f:type:`lammps` Derived Type -------- +.. f:subroutine:: set_fix_external_callback(id, callback, caller) + + Set the callback function for a :doc:`fix external ` instance + with the given ID. + + .. versionadded:: TBD + + Fix :doc:`external ` allows programs that are running LAMMPS + through its library interface to modify certain LAMMPS properties on + specific time steps, similar to the way other fixes do. + + This subroutine sets the callback function for use with the "pf/callback" + mode. The function should have Fortran language bindings with the following + interface, which depends on how LAMMPS was compiled: + + .. code-block:: Fortran + + ABSTRACT INTERFACE + SUBROUTINE external_callback(caller, timestep, ids, x, fexternal) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double, c_int64_t + CLASS(*), INTENT(IN) :: caller + INTEGER(c_bigint), INTENT(IN) :: timestep + INTEGER(c_tagint), DIMENSION(:), INTENT(IN) :: ids + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: fexternal + END SUBROUTINE external_callback + END INTERFACE + + where ``c_bigint`` is ``c_int`` if ``-DLAMMPS_SMALLSMALL`` was used and + ``c_int64_t`` otherwise; and ``c_tagint`` is ``c_int64_t`` if + ``-DLAMMPS_BIGBIG`` was used and ``c_int`` otherwise. + + The argument *caller* to :f:subr:`set_fix_external_callback` is unlimited + polymorphic (i.e., it can be any Fortran object you want to pass to the + calling function) and will be available as the first argument to the + callback function. It can be your LAMMPS instance, which you might need if + the callback function needs access to the library interface. + + The array *ids* is an array of length *nlocal* (as accessed from the + :cpp:class:`Atom` class or through :f:func:`extract_global`). The arrays + *x* and *fexternal* are :math:`3 \times {}`\ *nlocal* arrays; these are + transposed from what they would look like in C (see note about array index + order at :f:func:`extract_atom`). + + The callback mechanism is one of two ways that forces can be applied to a + simulation with the help of :doc:`fix external `. The + alternative is *array* mode, where one calls + :f:func:`fix_external_get_force`. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and couple it with external + programs. + + :p character(len=*) id: ID of :doc:`fix external ` instance + :p callback: subroutine :doc:`fix external ` should call + :ptype callback: external + :p class(*) caller [optional]: object you wish to pass to the callback + procedure + +-------- + .. f:subroutine:: flush_buffers() This function calls :cpp:func:`lammps_flush_buffers`, which flushes buffered diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index 02664c54db..eb32dcd266 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -75,6 +75,9 @@ MODULE LIBLAMMPS LMP_VAR_VECTOR = 2, & ! vector variables LMP_VAR_STRING = 3 ! string variables (everything else) + ! Constants we set once (in the constructor) and never need to check again + INTEGER(c_int), SAVE :: SIZE_TAGINT, SIZE_BIGINT, SIZE_IMAGEINT + ! "Constants" to use with extract_compute and friends TYPE lammps_style INTEGER(c_int) :: global, atom, local @@ -865,6 +868,11 @@ CONTAINS lmp_open%type%scalar = LMP_TYPE_SCALAR lmp_open%type%vector = LMP_TYPE_VECTOR lmp_open%type%array = LMP_TYPE_ARRAY + + ! Assign constants for bigint and tagint for use elsewhere + SIZE_TAGINT = lmp_extract_setting(lmp_open, 'tagint') + SIZE_BIGINT = lmp_extract_setting(lmp_open, 'bigint') + SIZE_IMAGEINT = lmp_extract_setting(lmp_open, 'imageint') END FUNCTION lmp_open ! Combined Fortran wrapper around lammps_close() and lammps_mpi_finalize() @@ -1731,19 +1739,16 @@ CONTAINS SUBROUTINE lmp_gather_bonds_small(self, data) CLASS(lammps), INTENT(IN) :: self INTEGER(c_int), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data - INTEGER(c_int) :: size_tagint, size_bigint INTEGER(c_int), POINTER :: nbonds_small INTEGER(c_int64_t), POINTER :: nbonds_big TYPE(c_ptr) :: Cdata - size_tagint = lmp_extract_setting(self, 'tagint') - IF (size_tagint /= 4_c_int) THEN + IF (SIZE_TAGINT /= 4_c_int) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & 'Incompatible integer kind in gather_bonds [Fortran API]') END IF IF (ALLOCATED(data)) DEALLOCATE(data) - size_bigint = lmp_extract_setting(self, 'bigint') - IF (size_bigint == 4_c_int) THEN + IF (SIZE_BIGINT == 4_c_int) THEN nbonds_small = lmp_extract_global(self, 'nbonds') ALLOCATE(data(3*nbonds_small)) ELSE @@ -1758,12 +1763,10 @@ CONTAINS SUBROUTINE lmp_gather_bonds_big(self, data) CLASS(lammps), INTENT(IN) :: self INTEGER(c_int64_t), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data - INTEGER(c_int) :: size_tagint INTEGER(c_int64_t), POINTER :: nbonds TYPE(c_ptr) :: Cdata - size_tagint = lmp_extract_setting(self, 'tagint') - IF (size_tagint /= 8_c_int) THEN + IF (SIZE_TAGINT /= 8_c_int) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & 'Incompatible integer kind in gather_bonds [Fortran API]') END IF @@ -2385,10 +2388,7 @@ CONTAINS CLASS(*), INTENT(IN), TARGET, OPTIONAL :: caller TYPE(c_ptr) :: c_id, c_caller TYPE(c_funptr) :: c_callback - INTEGER :: i, this_fix, size_tagint, size_bigint - - size_tagint = lmp_extract_setting(self, 'tagint') - size_bigint = lmp_extract_setting(self, 'bigint') + INTEGER :: i, this_fix c_id = f2c_string(id) IF (ALLOCATED(ext_data)) THEN @@ -2400,7 +2400,10 @@ CONTAINS END IF END DO IF (this_fix > SIZE(ext_data)) THEN + ! reallocates ext_data; this requires us to re-bind "caller" on the C + ! side to the new data structure, which likely moved to a new address ext_data = [ext_data, fix_external_data()] ! extends ext_data by 1 + CALL rebind_external_callback_data() END IF ELSE ALLOCATE(ext_data(1)) @@ -2409,11 +2412,11 @@ CONTAINS ext_data(this_fix)%id = id ext_data(this_fix)%lammps_instance => self - IF (size_tagint == 4_c_int .AND. size_bigint == 4_c_int) THEN + IF (SIZE_TAGINT == 4_c_int .AND. SIZE_BIGINT == 4_c_int) THEN ! -DSMALLSMALL c_callback = C_FUNLOC(callback_wrapper_smallsmall) CALL set_fix_external_callback_smallsmall(this_fix, callback) - ELSE IF (size_tagint == 8_c_int .AND. size_bigint == 8_c_int) THEN + ELSE IF (SIZE_TAGINT == 8_c_int .AND. SIZE_BIGINT == 8_c_int) THEN ! -DBIGBIG c_callback = C_FUNLOC(callback_wrapper_bigbig) CALL set_fix_external_callback_bigbig(this_fix, callback) @@ -2431,6 +2434,7 @@ CONTAINS c_caller = C_LOC(ext_data(this_fix)) CALL lammps_set_fix_external_callback(self%handle, c_id, c_callback, & c_caller) + CALL lammps_free(c_id) END SUBROUTINE lmp_set_fix_external_callback ! Wrappers to assign callback pointers with explicit interfaces @@ -2455,6 +2459,28 @@ CONTAINS ext_data(id)%callback_bigbig => callback END SUBROUTINE set_fix_external_callback_bigbig + ! subroutine that re-binds all external callback data after a reallocation + SUBROUTINE rebind_external_callback_data() + INTEGER :: i + TYPE(c_ptr) :: c_id, c_caller + TYPE(c_funptr) :: c_callback + + DO i = 1, SIZE(ext_data) - 1 + c_id = f2c_string(ext_data(i)%id) + c_caller = C_LOC(ext_data(i)) + IF (SIZE_TAGINT == 4_c_int .AND. SIZE_BIGINT == 4_c_int) THEN + c_callback = C_FUNLOC(callback_wrapper_smallsmall) + ELSE IF (SIZE_TAGINT == 8_c_int .AND. SIZE_BIGINT == 8_c_int) THEN + c_callback = C_FUNLOC(callback_wrapper_bigbig) + ELSE + c_callback = C_FUNLOC(callback_wrapper_smallbig) + END IF + CALL lammps_set_fix_external_callback( & + ext_data(i)%lammps_instance%handle, c_id, c_callback, c_caller) + CALL lammps_free(c_id) + END DO + END SUBROUTINE rebind_external_callback_data + ! companions to lmp_set_fix_external_callback to change interface SUBROUTINE callback_wrapper_smallsmall(caller, timestep, nlocal, ids, x, & fexternal) BIND(C) From 5f9956405a1b86112be8db24402f9e9bfdabd147 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Tue, 29 Nov 2022 15:37:15 -0600 Subject: [PATCH 03/19] Updated docs and wrote unit tests for lmp_set_fix_external_callback; fixed typos --- doc/src/Fortran.rst | 46 +++- fortran/lammps.f90 | 46 ++-- unittest/fortran/CMakeLists.txt | 4 + unittest/fortran/test_fortran_fixexternal.f90 | 220 ++++++++++++++++++ unittest/fortran/wrap_fixexternal.cpp | 75 ++++++ 5 files changed, 363 insertions(+), 28 deletions(-) create mode 100644 unittest/fortran/test_fortran_fixexternal.f90 create mode 100644 unittest/fortran/wrap_fixexternal.cpp diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 12d39e8af7..813a52a44b 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -409,7 +409,7 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. USE LIBLAMMPS USE MPI_F08 TYPE(lammps) :: lmp - lmp = lammps(MPI_COMM_SELF%MPI_VAL) + lmp = lammps(comm=MPI_COMM_SELF%MPI_VAL) END PROGRAM testmpi .. f:type:: lammps_style @@ -773,8 +773,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type Note that this function actually does not return a pointer, but rather associates the pointer on the left side of the assignment to point to internal LAMMPS data. Pointers must be of the correct type, kind, and - rank (e.g., ``INTEGER(c_int), DIMENSION(:)`` for "type", "mask", or "tag"; - ``INTEGER(c_int64_t), DIMENSION(:)`` for "tag" if LAMMPS was compiled + rank (e.g., ``INTEGER(c_int), DIMENSION(:)`` for "type", "mask", or "id"; + ``INTEGER(c_int64_t), DIMENSION(:)`` for "id" if LAMMPS was compiled with the ``-DLAMMPS_BIGBIG`` flag; ``REAL(c_double), DIMENSION(:,:)`` for "x", "v", or "f"; and so forth). The pointer being associated with LAMMPS data is type-, kind-, and rank-checked at run-time. @@ -2118,7 +2118,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type ABSTRACT INTERFACE SUBROUTINE external_callback(caller, timestep, ids, x, fexternal) USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double, c_int64_t - CLASS(*), INTENT(IN) :: caller + CLASS(*), INTENT(INOUT) :: caller INTEGER(c_bigint), INTENT(IN) :: timestep INTEGER(c_tagint), DIMENSION(:), INTENT(IN) :: ids REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x @@ -2135,6 +2135,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type calling function) and will be available as the first argument to the callback function. It can be your LAMMPS instance, which you might need if the callback function needs access to the library interface. + The argument must be a scalar; to pass non-scalar data, wrap those data in + a derived type and pass an instance of the derived type to *caller*. The array *ids* is an array of length *nlocal* (as accessed from the :cpp:class:`Atom` class or through :f:func:`extract_global`). The arrays @@ -2155,7 +2157,41 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p callback: subroutine :doc:`fix external ` should call :ptype callback: external :p class(*) caller [optional]: object you wish to pass to the callback - procedure + procedure (must be a scalar; see note) + + .. note:: + + The interface for your callback function must match types precisely + with the abstract interface block given above. **The compiler probably + will not be able to check this for you.** In particular, the first + argument ("caller") must be of type ``CLASS(*)`` or you will probably + get a segmentation fault or at least a misinterpretation of whatever is + in memory there. You can resolve the object using the ``SELECT TYPE`` + construct. An example callback function (assuming LAMMPS was compiled + with ``-DLAMMPS_SMALLBIG``) that applies something akin to Hooke's Law + (with each atom having a different *k* value) is shown below. + + .. code-block:: Fortran + + TYPE shield + REAL(c_double), DIMENSION(:), ALLOCATABLE :: k + ! assume k gets allocated to dimension(3,nlocal) at some point + END TYPE shield + + SUBROUTINE my_callback(caller, timestep, ids, x, fexternal) + CLASS(*), INTENT(INOUT) :: caller + INTEGER(c_int), INTENT(IN) :: timestep + INTEGER(c_int64_t), INTENT(IN) :: ids + REAL(c_double), INTENT(IN) :: x(:,:) + REAL(c_double), INTENT(OUT) :: fexternal(:,:) + + SELECT TYPE (caller) + TYPE IS (shield) + fexternal = - caller%k * x + CLASS DEFAULT + WRITE(error_unit,*) 'UH OH...' + END SELECT + END SUBROUTINE my_callback -------- diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index eb32dcd266..b6e4d5f56b 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -55,7 +55,7 @@ MODULE LIBLAMMPS LAMMPS_DOUBLE_2D = 3, & ! two-dimensional 64-bit double array LAMMPS_INT64 = 4, & ! 64-bit integer (or array) LAMMPS_INT64_2D = 5, & ! two-dimensional 64-bit integer array - LAMMPS_STRING = 6, & ! C-String + LAMMPS_STRING = 6, & ! string LMP_STYLE_GLOBAL = 0, & ! request global compute/fix/etc. data LMP_STYLE_ATOM = 1, & ! request per-atom compute/fix/etc. data LMP_STYLE_LOCAL = 2, & ! request local compute/fix/etc. data @@ -64,7 +64,7 @@ MODULE LIBLAMMPS LMP_TYPE_ARRAY = 2, & ! request array LMP_SIZE_VECTOR = 3, & ! request size of vector LMP_SIZE_ROWS = 4, & ! request rows (actually columns) - LMP_SIZE_COLS = 5, & ! request colums (actually rows) + LMP_SIZE_COLS = 5, & ! request columns (actually rows) LMP_ERROR_WARNING = 0, & ! call Error::warning() LMP_ERROR_ONE = 1, & ! call Error::one() (from this MPI rank) LMP_ERROR_ALL = 2, & ! call Error::all() (from all MPI ranks) @@ -268,7 +268,7 @@ MODULE LIBLAMMPS ABSTRACT INTERFACE SUBROUTINE external_callback_smallsmall(caller, timestep, ids, x, fexternal) IMPORT :: c_int, c_double - CLASS(*), INTENT(IN) :: caller + CLASS(*), INTENT(INOUT) :: caller INTEGER(c_int), INTENT(IN) :: timestep INTEGER(c_int), DIMENSION(:), INTENT(IN) :: ids REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x @@ -276,7 +276,7 @@ MODULE LIBLAMMPS END SUBROUTINE external_callback_smallsmall SUBROUTINE external_callback_smallbig(caller, timestep, ids, x, fexternal) IMPORT :: c_int, c_double, c_int64_t - CLASS(*), INTENT(IN) :: caller + CLASS(*), INTENT(INOUT) :: caller INTEGER(c_int64_t), INTENT(IN) :: timestep INTEGER(c_int), DIMENSION(:), INTENT(IN) :: ids REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x @@ -284,7 +284,7 @@ MODULE LIBLAMMPS END SUBROUTINE external_callback_smallbig SUBROUTINE external_callback_bigbig(caller, timestep, ids, x, fexternal) IMPORT :: c_double, c_int64_t - CLASS(*), INTENT(IN) :: caller + CLASS(*), INTENT(INOUT) :: caller INTEGER(c_int64_t), INTENT(IN) :: timestep INTEGER(c_int64_t), DIMENSION(:), INTENT(IN) :: ids REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x @@ -836,28 +836,28 @@ CONTAINS INTEGER(c_int) :: i, c_comm, argc IF (PRESENT(args)) THEN - ! convert fortran argument list to c style - argc = SIZE(args) - ALLOCATE(argv(argc)) - DO i=1, argc - argv(i) = f2c_string(args(i)) - END DO + ! convert fortran argument list to c style + argc = SIZE(args) + ALLOCATE(argv(argc)) + DO i=1, argc + argv(i) = f2c_string(args(i)) + END DO ELSE - argc = 1 - ALLOCATE(argv(1)) - argv(1) = f2c_string("liblammps") + argc = 1 + ALLOCATE(argv(1)) + argv(1) = f2c_string("liblammps") ENDIF IF (PRESENT(comm)) THEN - c_comm = comm - lmp_open%handle = lammps_open(argc, argv, c_comm) + c_comm = comm + lmp_open%handle = lammps_open(argc, argv, c_comm) ELSE - lmp_open%handle = lammps_open_no_mpi(argc, argv, c_null_ptr) + lmp_open%handle = lammps_open_no_mpi(argc, argv, c_null_ptr) END IF ! Clean up allocated memory DO i=1, argc - CALL lammps_free(argv(i)) + CALL lammps_free(argv(i)) END DO DEALLOCATE(argv) @@ -883,10 +883,10 @@ CONTAINS CALL lammps_close(self%handle) IF (PRESENT(finalize)) THEN - IF (finalize) THEN - CALL lammps_kokkos_finalize() - CALL lammps_mpi_finalize() - END IF + IF (finalize) THEN + CALL lammps_kokkos_finalize() + CALL lammps_mpi_finalize() + END IF END IF END SUBROUTINE lmp_close @@ -1054,7 +1054,7 @@ CONTAINS length = 3 CASE DEFAULT length = 1 - ! string cases doesn't use "length" + ! string cases do not use "length" END SELECT Cname = f2c_string(name) diff --git a/unittest/fortran/CMakeLists.txt b/unittest/fortran/CMakeLists.txt index f66b333db1..d12b274765 100644 --- a/unittest/fortran/CMakeLists.txt +++ b/unittest/fortran/CMakeLists.txt @@ -104,6 +104,10 @@ if(CMAKE_Fortran_COMPILER) target_link_libraries(test_fortran_neighlist PRIVATE flammps lammps MPI::MPI_Fortran GTest::GMockMain) add_test(NAME FortranNeighlist COMMAND test_fortran_neighlist) + add_executable(test_fortran_fixexternal wrap_fixexternal.cpp test_fortran_fixexternal.f90) + target_link_libraries(test_fortran_fixexternal PRIVATE flammps lammps MPI::MPI_Fortran GTest::GMockMain) + add_test(NAME FortranFixExternal COMMAND test_fortran_fixexternal) + else() message(STATUS "Skipping Tests for the LAMMPS Fortran Module: no Fortran compiler") endif() diff --git a/unittest/fortran/test_fortran_fixexternal.f90 b/unittest/fortran/test_fortran_fixexternal.f90 new file mode 100644 index 0000000000..1da8909627 --- /dev/null +++ b/unittest/fortran/test_fortran_fixexternal.f90 @@ -0,0 +1,220 @@ +MODULE ext_stuff + USE, INTRINSIC :: ISO_Fortran_ENV, ONLY : error_unit + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int, c_int64_t, c_loc + USE LIBLAMMPS + IMPLICIT NONE + + REAL(c_double), SAVE :: direction = 1.0_c_double + +CONTAINS + + SUBROUTINE f_lammps_reverse_direction() BIND(C) + direction = -direction + END SUBROUTINE f_lammps_reverse_direction + + SUBROUTINE f_callback_ss(instance, timestep, id, x, f) + CLASS(*), INTENT(INOUT) :: instance + INTEGER(c_int) :: timestep + INTEGER(c_int), DIMENSION(:), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + + WHERE (id == 1) + f(1,:) = 1.0_c_double + f(2,:) = -1.0_c_double + f(3,:) = 1.25_c_double + ELSEWHERE + f(1,:) = -1.0_c_double + f(2,:) = +1.0_c_double + f(3,:) = -1.25_c_double + END WHERE + END SUBROUTINE f_callback_ss + + SUBROUTINE f_callback_sb(instance, timestep, id, x, f) + CLASS(*), INTENT(INOUT) :: instance + INTEGER(c_int64_t) :: timestep + INTEGER(c_int), DIMENSION(:), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + + WHERE (id == 1_c_int) + f(1,:) = 1.0_c_double + f(2,:) = -1.0_c_double + f(3,:) = 1.25_c_double + ELSEWHERE + f(1,:) = -1.0_c_double + f(2,:) = +1.0_c_double + f(3,:) = -1.25_c_double + END WHERE + END SUBROUTINE f_callback_sb + + SUBROUTINE f_callback_bb(instance, timestep, id, x, f) + CLASS(*), INTENT(INOUT) :: instance + INTEGER(c_int64_t) :: timestep + INTEGER(c_int64_t), DIMENSION(:), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + + WHERE (id == 1_c_int64_t) + f(1,:) = 1.0_c_double + f(2,:) = -1.0_c_double + f(3,:) = 1.25_c_double + ELSEWHERE + f(1,:) = -1.0_c_double + f(2,:) = +1.0_c_double + f(3,:) = -1.25_c_double + END WHERE + END SUBROUTINE f_callback_bb + + SUBROUTINE f_callback2_ss(entity, timestep, id, x, f) + CLASS(*), INTENT(INOUT), target :: entity + INTEGER(c_int) :: timestep + INTEGER(c_int), DIMENSION(:), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + + SELECT TYPE (entity) + TYPE IS (REAL(c_double)) + WHERE (id == 1_c_int) + f(1,:) = SIGN(1.0_c_double, entity) * 2.0_c_double + f(2,:) = SIGN(1.0_c_double, entity) * (-2.0_c_double) + f(3,:) = SIGN(1.0_c_double, entity) * 2.5_c_double + ELSEWHERE + f(1,:) = SIGN(1.0_c_double, entity) * (-2.0_c_double) + f(2,:) = SIGN(1.0_c_double, entity) * 2.0_c_double + f(3,:) = SIGN(1.0_c_double, entity) * (-2.5_c_double) + END WHERE + CLASS DEFAULT + WRITE(error_unit,'(A)') 'ERROR: Failed to resolve "entity" in& + & f_callback2_ss' + STOP 1 + END SELECT + END SUBROUTINE f_callback2_ss + + SUBROUTINE f_callback2_sb(entity, timestep, id, x, f) + CLASS(*), INTENT(INOUT), target :: entity + INTEGER(c_int64_t) :: timestep + INTEGER(c_int), DIMENSION(:), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + + SELECT TYPE (entity) + TYPE IS (REAL(c_double)) + WHERE (id == 1_c_int) + f(1,:) = SIGN(1.0_c_double, entity) * 2.0_c_double + f(2,:) = SIGN(1.0_c_double, entity) * (-2.0_c_double) + f(3,:) = SIGN(1.0_c_double, entity) * 2.5_c_double + ELSEWHERE + f(1,:) = SIGN(1.0_c_double, entity) * (-2.0_c_double) + f(2,:) = SIGN(1.0_c_double, entity) * 2.0_c_double + f(3,:) = SIGN(1.0_c_double, entity) * (-2.5_c_double) + END WHERE + CLASS DEFAULT + WRITE(error_unit,'(A)') 'ERROR: Failed to resolve "entity" in& + & f_callback2_sb' + STOP 1 + END SELECT + END SUBROUTINE f_callback2_sb + + SUBROUTINE f_callback2_bb(entity, timestep, id, x, f) + CLASS(*), INTENT(INOUT), target :: entity + INTEGER(c_int64_t) :: timestep + INTEGER(c_int64_t), DIMENSION(:), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x + REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + + SELECT TYPE (entity) + TYPE IS (REAL(c_double)) + WHERE (id == 1_c_int64_t) + f(1,:) = SIGN(1.0_c_double, entity) * 2.0_c_double + f(2,:) = SIGN(1.0_c_double, entity) * (-2.0_c_double) + f(3,:) = SIGN(1.0_c_double, entity) * 2.5_c_double + ELSEWHERE + f(1,:) = SIGN(1.0_c_double, entity) * (-2.0_c_double) + f(2,:) = SIGN(1.0_c_double, entity) * 2.0_c_double + f(3,:) = SIGN(1.0_c_double, entity) * (-2.5_c_double) + END WHERE + CLASS DEFAULT + WRITE(error_unit,'(A)') 'ERROR: Failed to resolve "entity" in& + & f_callback2_sb' + STOP 1 + END SELECT + END SUBROUTINE f_callback2_bb +END MODULE ext_stuff + +FUNCTION f_lammps_with_args() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY: c_ptr + USE liblammps + USE keepstuff, ONLY: lmp + IMPLICIT NONE + TYPE(c_ptr) :: f_lammps_with_args + + CHARACTER(len=12), DIMENSION(12), PARAMETER :: args = & + [ CHARACTER(len=12) :: 'liblammps', '-log', 'none', & + '-echo','screen','-nocite','-var','zpos','1.5','-var','x','2'] + + lmp = lammps(args) + f_lammps_with_args = lmp%handle +END FUNCTION f_lammps_with_args + +SUBROUTINE f_lammps_close() BIND(C) + USE ISO_C_BINDING, ONLY: c_null_ptr + USE liblammps + USE keepstuff, ONLY: lmp + IMPLICIT NONE + + CALL lmp%close() + lmp%handle = c_null_ptr +END SUBROUTINE f_lammps_close + +SUBROUTINE f_lammps_setup_fix_external() BIND(C) + USE LIBLAMMPS + USE keepstuff, ONLY : lmp, demo_input, cont_input, pair_input + IMPLICIT NONE + + CALL lmp%commands_list(demo_input) + CALL lmp%commands_list(cont_input) + CALL lmp%commands_list(pair_input) + CALL lmp%command('neigh_modify exclude group all all') + CALL lmp%command('fix ext1 all external pf/callback 1 1') + CALL lmp%command('fix ext2 all external pf/callback 1 1') +END SUBROUTINE f_lammps_setup_fix_external + +SUBROUTINE f_lammps_set_fix_external_callbacks() BIND(C) + USE ISO_C_BINDING, ONLY : c_int + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + USE ext_stuff + IMPLICIT NONE + INTEGER :: size_bigint, size_tagint, nlocal + + nlocal = lmp%extract_setting('nlocal') + + size_bigint = lmp%extract_setting('bigint') + size_tagint = lmp%extract_setting('tagint') + IF (size_bigint == 4_c_int .AND. size_tagint == 4_c_int) THEN + CALL lmp%set_fix_external_callback('ext1', f_callback_ss) + CALL lmp%set_fix_external_callback('ext2', f_callback2_ss, direction) + ELSE IF (size_bigint == 8_c_int .AND. size_tagint == 8_c_int) THEN + CALL lmp%set_fix_external_callback('ext1', f_callback_bb) + CALL lmp%set_fix_external_callback('ext2', f_callback2_bb, direction) + ELSE + CALL lmp%set_fix_external_callback('ext1', f_callback_sb) + CALL lmp%set_fix_external_callback('ext2', f_callback2_sb, direction) + END IF +END SUBROUTINE f_lammps_set_fix_external_callbacks + +SUBROUTINE f_lammps_get_force (i, ptr) BIND(C) + USE ISO_C_BINDING, ONLY : c_int, c_double, c_ptr, C_F_POINTER + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i + TYPE(c_ptr), INTENT(IN), VALUE :: ptr + REAL(c_double), DIMENSION(:,:), POINTER :: force => NULL() + REAL(c_double), DIMENSION(:), POINTER :: f => NULL() + + CALL C_F_POINTER(ptr, f, [3]) + force = lmp%extract_atom('f') + f = force(:,i) +END SUBROUTINE f_lammps_get_force diff --git a/unittest/fortran/wrap_fixexternal.cpp b/unittest/fortran/wrap_fixexternal.cpp new file mode 100644 index 0000000000..282b95a299 --- /dev/null +++ b/unittest/fortran/wrap_fixexternal.cpp @@ -0,0 +1,75 @@ + +// unit tests for gathering and scattering data from a LAMMPS instance through +// the Fortran wrapper + +#include "lammps.h" +#include "library.h" +#include +#include +#include +#include + +#include "gtest/gtest.h" + +// prototypes for Fortran reverse wrapper functions +extern "C" { +void *f_lammps_with_args(); +void f_lammps_close(); +void f_lammps_setup_fix_external(); +void f_lammps_set_fix_external_callbacks(); +void f_lammps_get_force(int, double*); +void f_lammps_reverse_direction(); +} + +using namespace LAMMPS_NS; + +class LAMMPS_fixexternal : public ::testing::Test { +protected: + LAMMPS_NS::LAMMPS *lmp; + LAMMPS_fixexternal() = default; + ~LAMMPS_fixexternal() override = default; + + void SetUp() override + { + ::testing::internal::CaptureStdout(); + lmp = (LAMMPS_NS::LAMMPS *)f_lammps_with_args(); + std::string output = ::testing::internal::GetCapturedStdout(); + EXPECT_STREQ(output.substr(0, 8).c_str(), "LAMMPS ("); + } + void TearDown() override + { + ::testing::internal::CaptureStdout(); + f_lammps_close(); + std::string output = ::testing::internal::GetCapturedStdout(); + EXPECT_STREQ(output.substr(0, 16).c_str(), "Total wall time:"); + lmp = nullptr; + } +}; + +TEST_F(LAMMPS_fixexternal, callback) +{ + f_lammps_setup_fix_external(); + f_lammps_set_fix_external_callbacks(); + lammps_command(lmp, "run 0"); + double f[3]; + f_lammps_get_force(1,f); + EXPECT_DOUBLE_EQ(f[0], 3.0); + EXPECT_DOUBLE_EQ(f[1], -3.0); + EXPECT_DOUBLE_EQ(f[2], 3.75); + f_lammps_get_force(2,f); + EXPECT_DOUBLE_EQ(f[0], -3.0); + EXPECT_DOUBLE_EQ(f[1], 3.0); + EXPECT_DOUBLE_EQ(f[2], -3.75); + + f_lammps_reverse_direction(); + f_lammps_set_fix_external_callbacks(); + lammps_command(lmp, "run 0"); + f_lammps_get_force(1,f); + EXPECT_DOUBLE_EQ(f[0], -1.0); + EXPECT_DOUBLE_EQ(f[1], 1.0); + EXPECT_DOUBLE_EQ(f[2], -1.25); + f_lammps_get_force(2,f); + EXPECT_DOUBLE_EQ(f[0], 1.0); + EXPECT_DOUBLE_EQ(f[1], -1.0); + EXPECT_DOUBLE_EQ(f[2], 1.25); +}; From c674f0864d0a551b0da804cacb0847119e092477 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Tue, 29 Nov 2022 15:56:26 -0600 Subject: [PATCH 04/19] Added lammps_fix_external_get_force to C library utility doc page --- doc/src/Library_utility.rst | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/doc/src/Library_utility.rst b/doc/src/Library_utility.rst index 9d16daae0b..e555b79c0b 100644 --- a/doc/src/Library_utility.rst +++ b/doc/src/Library_utility.rst @@ -7,6 +7,7 @@ functions. They do not directly call the LAMMPS library. - :cpp:func:`lammps_encode_image_flags` - :cpp:func:`lammps_decode_image_flags` - :cpp:func:`lammps_set_fix_external_callback` +- :cpp:func:`lammps_fix_external_get_force` - :cpp:func:`lammps_fix_external_set_energy_global` - :cpp:func:`lammps_fix_external_set_energy_peratom` - :cpp:func:`lammps_fix_external_set_virial_global` @@ -44,6 +45,11 @@ where such memory buffers were allocated that require the use of ----------------------- +.. doxygenfunction:: lammps_fix_external_get_force + :project: progguide + +----------------------- + .. doxygenfunction:: lammps_fix_external_set_energy_global :project: progguide From aecd3841bed9a02a3e705fd21dab143363e94fa9 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Tue, 29 Nov 2022 18:28:52 -0600 Subject: [PATCH 05/19] Initial implementation of fix_external_get_force --- fortran/lammps.f90 | 24 +++++++++++++++++++++++- 1 file changed, 23 insertions(+), 1 deletion(-) diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index b6e4d5f56b..f6818a21d9 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -758,7 +758,12 @@ MODULE LIBLAMMPS TYPE(c_ptr), VALUE :: handle, id, ptr TYPE(c_funptr), VALUE :: funcptr END SUBROUTINE lammps_set_fix_external_callback - !FUNCTION lammps_fix_external_get_force() ! returns real(c_double)(:) + + FUNCTION lammps_fix_external_get_force(handle, id) BIND(C) + IMPORT :: c_ptr + TYPE(c_ptr), VALUE :: handle, id + TYPE(c_ptr) :: lammps_fix_external_get_force + END FUNCTION lammps_fix_external_get_force !SUBROUTINE lammps_fix_external_set_energy_global !SUBROUTINE lammps_fix_external_set_energy_peratom @@ -2569,6 +2574,23 @@ CONTAINS END IF END SUBROUTINE callback_wrapper_bigbig + ! equivalent function to lammps_fix_external_get_force + FUNCTION lmp_fix_external_get_force(self, id) RESULT(fexternal) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), POINTER :: fexternal + TYPE(c_ptr) :: ptr, Cid + TYPE(c_ptr), DIMENSION(:), POINTER :: f + INTEGER(c_int) :: nmax + + Cid = f2c_string(id) + ptr = lammps_fix_external_get_force(self%handle, Cid) + nmax = lmp_extract_setting(self, 'nmax') + CALL C_F_POINTER(ptr, f, [nmax]) + CALL C_F_POINTER(f(1), fexternal, [3, nmax]) + CALL lammps_free(Cid) + END FUNCTION lmp_fix_external_get_force + ! equivalent function to lammps_flush_buffers SUBROUTINE lmp_flush_buffers(self) CLASS(lammps), INTENT(IN) :: self From 8579b117afdca1db7b7fc56fa63a5d2f1c79fcef Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Tue, 29 Nov 2022 23:23:14 -0600 Subject: [PATCH 06/19] Implemented remaining fix_external functions and documentation --- doc/src/Fortran.rst | 230 ++++++++++++++++++++++++++++++++++++++++++++ fortran/lammps.f90 | 140 +++++++++++++++++++++++++-- 2 files changed, 364 insertions(+), 6 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 813a52a44b..874b57cbd7 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -2195,6 +2195,236 @@ Procedures Bound to the :f:type:`lammps` Derived Type -------- +.. f:function:: fix_external_get_force(id) + + Get pointer to the force array storage in a fix external instance with the + given ID. + + .. versionadded:: TBD + + Fix :doc:`external ` allows programs that are running LAMMPS + through its library interfaces to add or modify certain LAMMPS properties on + specific time steps, similar to the way other fixes do. + + This function provides access to the per-atom force storage in a fix + external instance with the given fix-ID to be added to the individual + atoms when using the "pf/array" mode. The *fexternal* array can be + accessed like other "native" per-atom arrays accessible via the + :f:func:`extract_atom` function. Please note that the array + stores the forces for *local* atoms for each MPI rank, in the order + determined by the neighbor list build. Because the underlying data + structures can change as well as the order of atom as they migrate between + MPI processes because of the domain decomposition parallelization, this + function should be always called immediately before the forces are going to + be set to get an up-to-date pointer. You can use, for example, + :f:func:`extract_setting` to obtain the number of local atoms `nlocal` and + then assume the dimensions of the returned force array as + ``REAL(c_double) :: force(3,nlocal)``. + + This function is an alternative to the callback mechanism in fix external + set up by :f:subr:`set_fix_external_callback`. The main difference is that + this mechanism can be used when forces are to be pre-computed and the + control alternates between LAMMPS and the external driver, while the + callback mechanism can call an external subroutine to compute the force when + the fix is triggered and needs them. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and how to couple it with an + external program. + + :p character(len=*) id: ID of :doc:`fix external ` instance + :to: :cpp:func:`lammps_fix_external_get_force` + :r fexternal: pointer to the per-atom force array allocated by the fix + :rtype fexternal: real(c_double), dimension(3,nlocal) + +-------- + +.. f:subroutine:: fix_external_set_energy_global(id, eng) + + Set the global energy contribution for a :doc:`fix external ` + instance with the given ID. + + .. versionadded:: TBD + + This is a companion function to :f:func:`set_fix_external_callback` + and :f:func:`fix_external_get_force` that also sets the contribution to the + global energy from the external program. The value of the *eng* argument + will be stored in the fix and applied on the current and all following + time steps until changed by another call to this function. The energy is in + energy units as determined by the current :doc:`units ` settings and + is the **total** energy of the contribution. Thus, when running in + parallel, all MPI processes have to call this function with the **same** + value, and this will be returned as scalar property of the fix external + instance when accessed in LAMMPS input commands or from variables. + + Please see the documentation for :doc:`fix external ` for more + information about how to use the fix and how to couple it with an external + program. + + :p character(len=*) id: fix ID of fix external instance + :p real(c_double) eng: total energy to be added to the global energy + :to: :cpp:func:`lammps_fix_external_set_energy_global` + +-------- + +.. f:subroutine:: fix_external_set_virial_global(id, virial) + + Set the global virial contribution for a fix external instance with the + given ID. + + .. versionadded:: TBD + + This is a companion function to :f:subr:`set_fix_external_callback` + and :f:func:`fix_external_get_force` to set the contribution to the global + virial from an external program. + + The six values of the *virial* array will be stored in the fix and applied + on the current and all following time steps until changed by another call + to this function. The components of the virial need to be stored in the + following order: *xx*, *yy*, *zz*, *xy*, *xz*, *yz*. In LAMMPS, the virial + is stored internally as `stress*volume` in units of `pressure*volume` as + determined by the current :doc:`units ` settings and is the + **total** contribution. Thus, when running in parallel, all MPI processes + have to call this function with the **same** value, and this will then be + added by fix external. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and how to couple it with an + external code. + + :p character(len=*) id: fix ID of fix external instance + :p real(c_double) virial [dimension(6)]: the six global stress tensor + components to be added to the global virial + :to: :cpp:func:`lammps_fix_external_set_virial_global` + +-------- + +.. f:subroutine:: fix_external_set_energy_peratom(id, eng) + + Set the per-atom energy contribution for a fix external instance with the + given ID. + + .. versionadded:: TBD + + This is a companion function to :f:subr:`set_fix_external_callback` to set + the per-atom energy contribution due to the fix from the external code as + part of the callback function. For this to work, the handle to the LAMMPS + object must be passed as the *ptr* argument when registering the callback + function. + + .. note:: + + This function is fully independent from + :f:subr:`fix_external_set_energy_global` and will **NOT** add any + contributions to the global energy tally and will **NOT** check whether + the sum of the contributions added here are consistent with the global + added energy. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and how to couple it with an + external code. + + :p character(len=*) id: fix ID of the fix external instance + :p real(c_double) eng [dimension(:)]: array of length nlocal containing + the energy to add to the per-atom energy + :to: :cpp:func:`lammps_fix_external_set_energy_peratom` + +-------- + +.. f:subroutine:: set_fix_external_set_virial_peratom(id, virial) + + This is a companion function to :f:subr:`set_fix_external_callback` to set + the per-atom virial contribution due to the fix from the external program as + part of the callback function. For this to work, the LAMMPS object must be + passed as the *caller* argument when registering the callback function. + + .. versionadded:: TBD + + .. note:: + + This function is fully independent from + :f:subr:`fix_external_set_virial_global` and will **NOT** add any + contributions to the global virial tally and **NOT** check whether the + sum of the contributions added here are consistent with the global added + virial. + + The order and units of the per-atom stress tensor elements are the same + as for the global virial. The type and dimensions of the per-atom virial + array must be ``REAL(c_double), DIMENSION(6,nlocal)``. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and how to couple it with an + external program. + + :p character(len=*) id: fix ID of fix external instance + :p real(c_double) virial [dimension(:,:)]: an array of :math:`6 \times{}`\ + *nlocal* components to be added to the per-atom virial + :to: :cpp:func:`lammps_set_virial_peratom` + +-------- + +.. f:subroutine:: fix_external_set_vector_length(id, length) + + Set the vector length for a global vector stored with fix external for + analysis. + + .. versionadded:: TBD + + This is a companion function to :f:subr:`set_fix_external_callback` and + :f:func:`fix_external_get_force` to set the length of a global vector of + properties that will be stored with the fix via + :f:subr:`fix_external_set_vector`. + + This function needs to be called **before** a call to + :f:subr:`fix_external_set_vector` and **before** a run or minimize command. + When running in parallel, it must be called from **all** MPI + processes with the same length argument. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and how to couple it with an + external program. + + :p character(len=*) id: fix ID of fix external instance + :p integer(c_int) length: length of the global vector to be stored with the + fix + +-------- + +.. f:subroutine:: fix_external_set_vector(id, idx, val) + + Store a global vector value for a fix external instance with the given ID. + + .. versionadded:: TBD + + This is a companion function to :f:subr:`set_fix_external_callback` and + :f:func:`fix_external_get_force` to set the values of a global vector of + properties that will be stored with the fix and can be accessed from + within LAMMPS input commands (e.g., fix ave/time or variables) when used + in a vector context. + + This function needs to be called **after** a call to + :f:subr:`fix_external_set_vector_length` and **before** a run or minimize + command. When running in parallel, it must be called from **all** MPI + processes with the **same**\ *idx* and *val* parameters. The variable + *val* is assumed to be extensive. + + .. note:: + + The index in the *idx* parameter is 1-based (i.e., the first element + is set with *idx*\ :math:`{} = 1`, and the last element of the vector + with *idx*\ :math:`{} = N`, where :math:`N` is the value of the *length* + parameter of the call to :f:subr:`fix_external_set_vector_length`. + + Please see the documentation for :doc:`fix external ` for + more information about how to use the fix and how to couple it with an + external code. + + :p character(len=*) id: ID of fix external instance + :p integer(c_int) idx: 1-based index in global vector + :p integer(c_int) val: value to be stored in global vector at index *idx* + +-------- + .. f:subroutine:: flush_buffers() This function calls :cpp:func:`lammps_flush_buffers`, which flushes buffered diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index f6818a21d9..c5ccc17b4e 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -175,6 +175,18 @@ MODULE LIBLAMMPS GENERIC :: decode_image_flags => lmp_decode_image_flags, & lmp_decode_image_flags_bigbig PROCEDURE :: set_fix_external_callback => lmp_set_fix_external_callback + PROCEDURE :: fix_external_get_force => lmp_fix_external_get_force + PROCEDURE :: fix_external_set_energy_global & + => lmp_fix_external_set_energy_global + PROCEDURE :: fix_external_set_virial_global & + => lmp_fix_external_set_virial_global + PROCEDURE :: fix_external_set_energy_peratom & + => lmp_fix_external_set_energy_peratom + PROCEDURE :: fix_external_set_virial_peratom & + => lmp_fix_external_set_virial_peratom + PROCEDURE :: fix_external_set_vector_length & + => lmp_fix_external_set_vector_length + PROCEDURE :: fix_external_set_vector => lmp_fix_external_set_vector PROCEDURE :: flush_buffers => lmp_flush_buffers PROCEDURE :: is_running => lmp_is_running PROCEDURE :: force_timeout => lmp_force_timeout @@ -765,12 +777,42 @@ MODULE LIBLAMMPS TYPE(c_ptr) :: lammps_fix_external_get_force END FUNCTION lammps_fix_external_get_force - !SUBROUTINE lammps_fix_external_set_energy_global - !SUBROUTINE lammps_fix_external_set_energy_peratom - !SUBROUTINE lammps_fix_external_set_virial_global - !SUBROUTINE lammps_fix_external_set_virial_peratom - !SUBROUTINE lammps_fix_external_set_vector_length - !SUBROUTINE lammps_fix_external_set_vector + SUBROUTINE lammps_fix_external_set_energy_global(handle, id, eng) BIND(C) + IMPORT :: c_ptr, c_double + TYPE(c_ptr), VALUE :: handle, id + REAL(c_double), VALUE :: eng + END SUBROUTINE lammps_fix_external_set_energy_global + + SUBROUTINE lammps_fix_external_set_virial_global(handle, id, virial) & + BIND(C) + IMPORT :: c_ptr + TYPE(c_ptr), VALUE :: handle, id, virial + END SUBROUTINE lammps_fix_external_set_virial_global + + SUBROUTINE lammps_fix_external_set_energy_peratom(handle, id, eng) BIND(C) + IMPORT :: c_ptr + TYPE(c_ptr), VALUE :: handle, id, eng + END SUBROUTINE lammps_fix_external_set_energy_peratom + + SUBROUTINE lammps_fix_external_set_virial_peratom(handle, id, virial) & + BIND(C) + IMPORT :: c_ptr + TYPE(c_ptr), VALUE :: handle, id, virial + END SUBROUTINE lammps_fix_external_set_virial_peratom + + SUBROUTINE lammps_fix_external_set_vector_length(handle, id, length) & + BIND(C) + IMPORT :: c_ptr, c_int + TYPE(c_ptr), VALUE :: handle, id + INTEGER(c_int), VALUE :: length + END SUBROUTINE lammps_fix_external_set_vector_length + + SUBROUTINE lammps_fix_external_set_vector(handle, id, idx, val) BIND(C) + IMPORT :: c_ptr, c_int, c_double + TYPE(c_ptr), VALUE :: handle, id + INTEGER(c_int), VALUE :: idx + REAL(c_double), VALUE :: val + END SUBROUTINE lammps_fix_external_set_vector SUBROUTINE lammps_flush_buffers(handle) BIND(C) IMPORT :: c_ptr @@ -2591,6 +2633,92 @@ CONTAINS CALL lammps_free(Cid) END FUNCTION lmp_fix_external_get_force + SUBROUTINE lmp_fix_external_set_energy_global(self, id, eng) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + REAL(c_double), INTENT(OUT) :: eng + TYPE(c_ptr) :: Cid + + Cid = f2c_string(id) + CALL lammps_fix_external_set_energy_global(self%handle, Cid, eng) + CALL lammps_free(Cid) + END SUBROUTINE lmp_fix_external_set_energy_global + + SUBROUTINE lmp_fix_external_set_virial_global(self, id, virial) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + REAL(c_double), DIMENSION(6), TARGET, INTENT(IN) :: virial + TYPE(c_ptr) :: Cid, Cvirial + + Cid = f2c_string(id) + Cvirial = C_LOC(virial(1)) + CALL lammps_fix_external_set_virial_global(self%handle, Cid, Cvirial) + CALL lammps_free(Cid) + END SUBROUTINE lmp_fix_external_set_virial_global + + SUBROUTINE lmp_fix_external_set_energy_peratom(self, id, eng) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + REAL(c_double), DIMENSION(:), TARGET, INTENT(IN) :: eng + TYPE(c_ptr) :: Cid, Ceng + INTEGER(c_int) :: nlocal + + nlocal = lmp_extract_setting(self, 'nlocal') + IF (SIZE(eng) < nlocal) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Array "eng" should be length nlocal or greater & + &[Fortran/fix_external_set_energy_peratom]') + END IF + Cid = f2c_string(id) + Ceng = C_LOC(eng) + CALL lammps_fix_external_set_energy_peratom(self%handle, Cid, Ceng) + CALL lammps_free(Cid) + END SUBROUTINE lmp_fix_external_set_energy_peratom + + SUBROUTINE lmp_fix_external_set_virial_peratom(self, id, virial) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + REAL(c_double), DIMENSION(:,:), TARGET, INTENT(IN) :: virial + TYPE(c_ptr) :: Cid, Cvirial + TYPE(c_ptr), TARGET :: Cptr + INTEGER(c_int) :: nlocal + + nlocal = lmp_extract_setting(self, 'nlocal') + IF (SIZE(virial,2) < nlocal .OR. SIZE(virial,1) /= 6) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Array "virial" should be size 6 x nlocal or greater & + &[Fortran/fix_external_set_energy_peratom]') + END IF + Cid = f2c_string(id) + Cptr = C_LOC(virial(1,1)) + Cvirial = C_LOC(Cptr) + CALL lammps_fix_external_set_virial_peratom(self%handle, Cid, Cvirial) + CALL lammps_free(Cid) + END SUBROUTINE lmp_fix_external_set_virial_peratom + + SUBROUTINE lmp_fix_external_set_vector_length(self, id, length) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + INTEGER(c_int), INTENT(IN) :: length + TYPE(c_ptr) :: Cid + + Cid = f2c_string(id) + CALL lammps_fix_external_set_vector_length(self%handle, Cid, length) + CALL lammps_free(Cid) + END SUBROUTINE lmp_fix_external_set_vector_length + + SUBROUTINE lmp_fix_external_set_vector(self, id, idx, val) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: id + INTEGER(c_int), INTENT(IN) :: idx + REAL(c_double), INTENT(IN) :: val + TYPE(c_ptr) :: Cid + + Cid = f2c_string(id) + CALL lammps_fix_external_set_vector(self%handle, Cid, idx, val) + CALL lammps_free(Cid) + END SUBROUTINE lmp_fix_external_set_vector + ! equivalent function to lammps_flush_buffers SUBROUTINE lmp_flush_buffers(self) CLASS(lammps), INTENT(IN) :: self From a87aff7b878d8170e023292c5016be2563307e71 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Wed, 30 Nov 2022 22:48:29 -0600 Subject: [PATCH 07/19] Fixed bug and wrote unit tests for fix_external_array functions --- doc/src/Fortran.rst | 13 +- fortran/lammps.f90 | 10 +- unittest/fortran/test_fortran_fixexternal.f90 | 132 ++++++++++++++++-- unittest/fortran/wrap_fixexternal.cpp | 66 ++++++++- 4 files changed, 201 insertions(+), 20 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 874b57cbd7..dc5001f4c9 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -2254,7 +2254,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type energy units as determined by the current :doc:`units ` settings and is the **total** energy of the contribution. Thus, when running in parallel, all MPI processes have to call this function with the **same** - value, and this will be returned as scalar property of the fix external + value, and this will be returned as a scalar property of the fix external instance when accessed in LAMMPS input commands or from variables. Please see the documentation for :doc:`fix external ` for more @@ -2307,10 +2307,11 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. versionadded:: TBD This is a companion function to :f:subr:`set_fix_external_callback` to set - the per-atom energy contribution due to the fix from the external code as - part of the callback function. For this to work, the handle to the LAMMPS - object must be passed as the *ptr* argument when registering the callback - function. + the per-atom energy contribution due to the fix from the external program as + part of the callback function. For this to work, the LAMMPS object must be + passed as part of the *caller* argument when registering the callback + function, or the callback function must otherwise have access to the + LAMMPS object, such as through a module-based pointer. .. note:: @@ -2325,7 +2326,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type external code. :p character(len=*) id: fix ID of the fix external instance - :p real(c_double) eng [dimension(:)]: array of length nlocal containing + :p real(c_double) eng [dimension(:)]: array of length *nlocal* containing the energy to add to the per-atom energy :to: :cpp:func:`lammps_fix_external_set_energy_peratom` diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index c5ccc17b4e..fc318091bf 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -2618,9 +2618,9 @@ CONTAINS ! equivalent function to lammps_fix_external_get_force FUNCTION lmp_fix_external_get_force(self, id) RESULT(fexternal) - CLASS(lammps), INTENT(IN) :: self + CLASS(lammps), TARGET, INTENT(IN) :: self CHARACTER(LEN=*), INTENT(IN) :: id - REAL(c_double), DIMENSION(:,:), POINTER :: fexternal + TYPE(lammps_fix_data) :: fexternal TYPE(c_ptr) :: ptr, Cid TYPE(c_ptr), DIMENSION(:), POINTER :: f INTEGER(c_int) :: nmax @@ -2629,14 +2629,16 @@ CONTAINS ptr = lammps_fix_external_get_force(self%handle, Cid) nmax = lmp_extract_setting(self, 'nmax') CALL C_F_POINTER(ptr, f, [nmax]) - CALL C_F_POINTER(f(1), fexternal, [3, nmax]) + fexternal%datatype = DATA_DOUBLE_2D + fexternal%lammps_instance => self + CALL C_F_POINTER(f(1), fexternal%r64_mat, [3, nmax]) CALL lammps_free(Cid) END FUNCTION lmp_fix_external_get_force SUBROUTINE lmp_fix_external_set_energy_global(self, id, eng) CLASS(lammps), INTENT(IN) :: self CHARACTER(LEN=*), INTENT(IN) :: id - REAL(c_double), INTENT(OUT) :: eng + REAL(c_double), INTENT(IN) :: eng TYPE(c_ptr) :: Cid Cid = f2c_string(id) diff --git a/unittest/fortran/test_fortran_fixexternal.f90 b/unittest/fortran/test_fortran_fixexternal.f90 index 1da8909627..ca08e3c4d9 100644 --- a/unittest/fortran/test_fortran_fixexternal.f90 +++ b/unittest/fortran/test_fortran_fixexternal.f90 @@ -5,6 +5,7 @@ MODULE ext_stuff IMPLICIT NONE REAL(c_double), SAVE :: direction = 1.0_c_double + REAL(c_double), DIMENSION(:,:), POINTER, SAVE :: f3 => NULL(), f4 => NULL() CONTAINS @@ -18,16 +19,26 @@ CONTAINS INTEGER(c_int), DIMENSION(:), INTENT(IN) :: id REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + REAL(c_double), DIMENSION(SIZE(id)) :: e WHERE (id == 1) f(1,:) = 1.0_c_double f(2,:) = -1.0_c_double f(3,:) = 1.25_c_double + e = 1.0_c_double ELSEWHERE f(1,:) = -1.0_c_double f(2,:) = +1.0_c_double f(3,:) = -1.25_c_double + e = 10.0_c_double END WHERE + SELECT TYPE (instance) + CLASS IS (lammps) + CALL instance%fix_external_set_energy_peratom('ext1', e) + CLASS DEFAULT + WRITE(0,*) 'UMM...this should never happen.' + STOP 1 + END SELECT END SUBROUTINE f_callback_ss SUBROUTINE f_callback_sb(instance, timestep, id, x, f) @@ -36,16 +47,26 @@ CONTAINS INTEGER(c_int), DIMENSION(:), INTENT(IN) :: id REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + REAL(c_double), DIMENSION(SIZE(id)) :: e WHERE (id == 1_c_int) f(1,:) = 1.0_c_double f(2,:) = -1.0_c_double f(3,:) = 1.25_c_double + e = 1.0_c_double ELSEWHERE f(1,:) = -1.0_c_double f(2,:) = +1.0_c_double f(3,:) = -1.25_c_double + e = 10.0_c_double END WHERE + SELECT TYPE (instance) + CLASS IS (lammps) + CALL instance%fix_external_set_energy_peratom('ext1', e) + CLASS DEFAULT + WRITE(0,*) 'UMM...this should never happen.' + STOP 1 + END SELECT END SUBROUTINE f_callback_sb SUBROUTINE f_callback_bb(instance, timestep, id, x, f) @@ -54,16 +75,26 @@ CONTAINS INTEGER(c_int64_t), DIMENSION(:), INTENT(IN) :: id REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f + REAL(c_double), DIMENSION(SIZE(id)) :: e WHERE (id == 1_c_int64_t) f(1,:) = 1.0_c_double f(2,:) = -1.0_c_double f(3,:) = 1.25_c_double + e = 1.0_c_double ELSEWHERE f(1,:) = -1.0_c_double f(2,:) = +1.0_c_double f(3,:) = -1.25_c_double + e = 10.0_c_double END WHERE + SELECT TYPE (instance) + CLASS IS (lammps) + CALL instance%fix_external_set_energy_peratom('ext1', e) + CLASS DEFAULT + WRITE(0,*) 'UMM...this should never happen.' + STOP 1 + END SELECT END SUBROUTINE f_callback_bb SUBROUTINE f_callback2_ss(entity, timestep, id, x, f) @@ -158,7 +189,7 @@ FUNCTION f_lammps_with_args() BIND(C) END FUNCTION f_lammps_with_args SUBROUTINE f_lammps_close() BIND(C) - USE ISO_C_BINDING, ONLY: c_null_ptr + USE, INTRINSIC :: ISO_C_BINDING, ONLY: c_null_ptr USE liblammps USE keepstuff, ONLY: lmp IMPLICIT NONE @@ -167,7 +198,7 @@ SUBROUTINE f_lammps_close() BIND(C) lmp%handle = c_null_ptr END SUBROUTINE f_lammps_close -SUBROUTINE f_lammps_setup_fix_external() BIND(C) +SUBROUTINE f_lammps_setup_fix_external_callback() BIND(C) USE LIBLAMMPS USE keepstuff, ONLY : lmp, demo_input, cont_input, pair_input IMPLICIT NONE @@ -178,10 +209,29 @@ SUBROUTINE f_lammps_setup_fix_external() BIND(C) CALL lmp%command('neigh_modify exclude group all all') CALL lmp%command('fix ext1 all external pf/callback 1 1') CALL lmp%command('fix ext2 all external pf/callback 1 1') -END SUBROUTINE f_lammps_setup_fix_external +END SUBROUTINE f_lammps_setup_fix_external_callback + +SUBROUTINE f_lammps_setup_fix_external_array() BIND(C) + USE LIBLAMMPS + USE keepstuff, ONLY : lmp, demo_input, cont_input, pair_input + USE ext_stuff, ONLY : f3, f4 + IMPLICIT NONE + + CALL lmp%commands_list(demo_input) + CALL lmp%commands_list(cont_input) + CALL lmp%commands_list(pair_input) + CALL lmp%command('neigh_modify exclude group all all') + CALL lmp%command('fix ext3 all external pf/array 1') + CALL lmp%command('fix ext4 all external pf/array 1') + CALL lmp%command('thermo_style custom step pxx pe etotal') + CALL lmp%command('thermo_modify norm no') + CALL lmp%command('thermo 100') + f3 = lmp%fix_external_get_force('ext3') + f4 = lmp%fix_external_get_force('ext4') +END SUBROUTINE f_lammps_setup_fix_external_array SUBROUTINE f_lammps_set_fix_external_callbacks() BIND(C) - USE ISO_C_BINDING, ONLY : c_int + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int USE LIBLAMMPS USE keepstuff, ONLY : lmp USE ext_stuff @@ -193,19 +243,19 @@ SUBROUTINE f_lammps_set_fix_external_callbacks() BIND(C) size_bigint = lmp%extract_setting('bigint') size_tagint = lmp%extract_setting('tagint') IF (size_bigint == 4_c_int .AND. size_tagint == 4_c_int) THEN - CALL lmp%set_fix_external_callback('ext1', f_callback_ss) + CALL lmp%set_fix_external_callback('ext1', f_callback_ss, lmp) CALL lmp%set_fix_external_callback('ext2', f_callback2_ss, direction) ELSE IF (size_bigint == 8_c_int .AND. size_tagint == 8_c_int) THEN - CALL lmp%set_fix_external_callback('ext1', f_callback_bb) + CALL lmp%set_fix_external_callback('ext1', f_callback_bb, lmp) CALL lmp%set_fix_external_callback('ext2', f_callback2_bb, direction) ELSE - CALL lmp%set_fix_external_callback('ext1', f_callback_sb) + CALL lmp%set_fix_external_callback('ext1', f_callback_sb, lmp) CALL lmp%set_fix_external_callback('ext2', f_callback2_sb, direction) END IF END SUBROUTINE f_lammps_set_fix_external_callbacks SUBROUTINE f_lammps_get_force (i, ptr) BIND(C) - USE ISO_C_BINDING, ONLY : c_int, c_double, c_ptr, C_F_POINTER + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double, c_ptr, C_F_POINTER USE LIBLAMMPS USE keepstuff, ONLY : lmp IMPLICIT NONE @@ -218,3 +268,69 @@ SUBROUTINE f_lammps_get_force (i, ptr) BIND(C) force = lmp%extract_atom('f') f = force(:,i) END SUBROUTINE f_lammps_get_force + +SUBROUTINE f_lammps_find_forces() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int, c_int64_t + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + USE ext_stuff, ONLY : f3, f4 + IMPLICIT NONE + INTEGER(c_int) :: size_tagint + INTEGER(c_int), DIMENSION(:), POINTER :: id + INTEGER(c_int64_t), DIMENSION(:), POINTER :: tag + + f3(:,:) = 0.0_c_double + f4(:,:) = 0.0_c_double + size_tagint = lmp%extract_setting('tagint') + IF (size_tagint == 4_c_int) THEN + id = lmp%extract_atom('id') + WHERE (id == 1_c_int) + f3(1,:) = 4.0_c_double + f3(2,:) = -4.0_c_double + f3(3,:) = 6.0_c_double + f4(1,:) = 10.0_c_double + f4(2,:) = -10.0_c_double + f4(3,:) = 12.0_c_double + ELSEWHERE + f3(1,:) = 5.0_c_double + f3(2,:) = -5.0_c_double + f3(3,:) = 7.0_c_double + f4(1,:) = 11.0_c_double + f4(2,:) = -11.0_c_double + f4(3,:) = 13.0_c_double + END WHERE + ELSE + tag = lmp%extract_atom('id') + WHERE (tag == 1_c_int64_t) + f3(1,:) = 4.0_c_double + f3(2,:) = -4.0_c_double + f3(3,:) = 6.0_c_double + ELSEWHERE + f3(1,:) = 5.0_c_double + f3(2,:) = -5.0_c_double + f3(3,:) = 7.0_c_double + f4(1,:) = 11.0_c_double + f4(2,:) = -11.0_c_double + f4(3,:) = 13.0_c_double + END WHERE + END IF +END SUBROUTINE f_lammps_find_forces + +SUBROUTINE f_lammps_add_energy() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + + CALL lmp%fix_external_set_energy_global('ext3', -20.2_c_double); +END SUBROUTINE f_lammps_add_energy + +SUBROUTINE f_lammps_set_virial() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + + CALL lmp%fix_external_set_virial_global('ext4', [1.0_c_double, & + 2.0_c_double, 2.5_c_double, -1.0_c_double, -2.25_c_double, -3.02_c_double]) +END SUBROUTINE f_lammps_set_virial diff --git a/unittest/fortran/wrap_fixexternal.cpp b/unittest/fortran/wrap_fixexternal.cpp index 282b95a299..aa6cb2b921 100644 --- a/unittest/fortran/wrap_fixexternal.cpp +++ b/unittest/fortran/wrap_fixexternal.cpp @@ -15,10 +15,14 @@ extern "C" { void *f_lammps_with_args(); void f_lammps_close(); -void f_lammps_setup_fix_external(); +void f_lammps_setup_fix_external_callback(); +void f_lammps_setup_fix_external_array(); void f_lammps_set_fix_external_callbacks(); void f_lammps_get_force(int, double*); void f_lammps_reverse_direction(); +void f_lammps_find_forces(); +void f_lammps_add_energy(); +void f_lammps_set_virial(); } using namespace LAMMPS_NS; @@ -48,7 +52,7 @@ protected: TEST_F(LAMMPS_fixexternal, callback) { - f_lammps_setup_fix_external(); + f_lammps_setup_fix_external_callback(); f_lammps_set_fix_external_callbacks(); lammps_command(lmp, "run 0"); double f[3]; @@ -73,3 +77,61 @@ TEST_F(LAMMPS_fixexternal, callback) EXPECT_DOUBLE_EQ(f[1], -1.0); EXPECT_DOUBLE_EQ(f[2], 1.25); }; + +TEST_F(LAMMPS_fixexternal, array) +{ + f_lammps_setup_fix_external_array(); + double **f; + f = (double**) lammps_extract_atom(lmp, "f"); + f_lammps_find_forces(); + lammps_command(lmp, "run 0"); + EXPECT_DOUBLE_EQ(f[0][0], 14.0); + EXPECT_DOUBLE_EQ(f[0][1], -14.0); + EXPECT_DOUBLE_EQ(f[0][2], 18.0); + EXPECT_DOUBLE_EQ(f[1][0], 16.0); + EXPECT_DOUBLE_EQ(f[1][1], -16.0); + EXPECT_DOUBLE_EQ(f[1][2], 20.0); +}; + +TEST_F(LAMMPS_fixexternal, energy_global) +{ + f_lammps_setup_fix_external_array(); + double energy; + f_lammps_add_energy(); + lammps_command(lmp, "run 0"); + energy = lammps_get_thermo(lmp, "etotal"); + EXPECT_DOUBLE_EQ(energy, -20.2); +}; + +TEST_F(LAMMPS_fixexternal, virial_global) +{ + f_lammps_setup_fix_external_array(); + double virial[6], volume; + f_lammps_set_virial(); + lammps_command(lmp, "run 0"); + volume = lammps_get_thermo(lmp, "vol"); + virial[0] = lammps_get_thermo(lmp, "pxx"); + virial[1] = lammps_get_thermo(lmp, "pyy"); + virial[2] = lammps_get_thermo(lmp, "pzz"); + virial[3] = lammps_get_thermo(lmp, "pxy"); + virial[4] = lammps_get_thermo(lmp, "pxz"); + virial[5] = lammps_get_thermo(lmp, "pyz"); + EXPECT_DOUBLE_EQ(virial[0], 1.0/volume); + EXPECT_DOUBLE_EQ(virial[1], 2.0/volume); + EXPECT_DOUBLE_EQ(virial[2], 2.5/volume); + EXPECT_DOUBLE_EQ(virial[3], -1.0/volume); + EXPECT_DOUBLE_EQ(virial[4], -2.25/volume); + EXPECT_DOUBLE_EQ(virial[5], -3.02/volume); +}; + +TEST_F(LAMMPS_fixexternal, energy_peratom) +{ + f_lammps_setup_fix_external_callback(); + f_lammps_set_fix_external_callbacks(); + double energy; + lammps_command(lmp, "run 0"); +/* FIXME: the per-atom energy is NOT summed up by this function! We need + another test. */ + energy = lammps_get_thermo(lmp, "pe"); + EXPECT_DOUBLE_EQ(energy, 11.0); +}; From 713c7d3508df9029085cf97c1f7c6d7ff7decdb1 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Thu, 1 Dec 2022 16:49:18 -0600 Subject: [PATCH 08/19] Cleaned up documentation --- doc/src/Fortran.rst | 185 ++++++++++++++++++++++++++------------------ 1 file changed, 111 insertions(+), 74 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index dc5001f4c9..6161ccaebe 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -30,7 +30,7 @@ as a static or dynamic library `. If the LAMMPS library itself has been compiled with MPI support, the resulting executable will still be able to run LAMMPS in parallel with -``mpirun``, ``mpiexec`` or equivalent. Please also note that the order +``mpirun``, ``mpiexec``, or equivalent. Please also note that the order of the source files matters: the ``lammps.f90`` file needs to be compiled first, since it provides the :f:mod:`LIBLAMMPS` module that is imported by the Fortran code that uses the interface. A working example @@ -42,14 +42,14 @@ can be found together with equivalent examples in C and C++ in the .. admonition:: Work in Progress :class: note - This Fortran module is work in progress and only the documented + This Fortran module is a work in progress, and only the documented functionality is currently available. The final implementation should cover the entire range of functionality available in the C and Python library interfaces. .. note:: - A contributed (and more complete!) Fortran interface that more + A contributed Fortran interface that more closely resembles the C library interface is available in the ``examples/COUPLE/fortran2`` folder. Please see the ``README`` file in that folder for more information about it and how to contact its @@ -79,7 +79,7 @@ function and triggered with the optional logical argument set to USE LIBLAMMPS ! include the LAMMPS library interface IMPLICIT NONE TYPE(lammps) :: lmp ! derived type to hold LAMMPS instance - CHARACTER(LEN=*), PARAMETER :: args(3) = & + CHARACTER(LEN=12), PARAMETER :: args(3) = & [ CHARACTER(LEN=12) :: 'liblammps', '-log', 'none' ] ! create a LAMMPS instance (and initialize MPI) @@ -136,11 +136,11 @@ is done similarly to how it is implemented in the :doc:`C library interface `. Before handing off the calls to the C library interface, the corresponding Fortran versions of the calls (:f:func:`file`, :f:func:`command`, :f:func:`commands_list`, and -:f:func:`commands_string`) have to make a copy of the strings passed as +:f:func:`commands_string`) have to make copies of the strings passed as arguments so that they can be modified to be compatible with the requirements of strings in C without affecting the original strings. Those copies are automatically deleted after the functions return. -Below is a small demonstration of the uses of the different functions: +Below is a small demonstration of the uses of the different functions. .. code-block:: fortran @@ -360,6 +360,20 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype decode_image_flags: subroutine :f set_fix_external_callback: :f:subr:`set_fix_external_callback` :ftype set_fix_external_callback: subroutine + :f fix_external_get_force: :f:func:`fix_external_get_force` + :ftype fix_external_get_force: function + :f fix_external_set_energy_global: :f:subr:`fix_external_set_energy_global` + :ftype fix_external_set_energy_global: subroutine + :f fix_external_set_virial_global: :f:subr:`fix_external_set_virial_global` + :ftype fix_external_set_virial_global: subroutine + :f fix_external_set_energy_peratom: :f:subr:`fix_external_set_energy_peratom` + :ftype fix_external_set_energy_peratom: subroutine + :f fix_external_set_virial_peratom: :f:subr:`fix_external_set_virial_peratom` + :ftype fix_external_set_virial_peratom: subroutine + :f fix_external_set_vector_length: :f:subr:`fix_external_set_vector_length` + :ftype fix_external_set_vector_length: subroutine + :f fix_external_set_vector: :f:subr:`fix_external_set_vector` + :ftype fix_external_set_vector: subroutine :f flush_buffers: :f:subr:`flush_buffers` :ftype flush_buffers: subroutine :f is_running: :f:func:`is_running` @@ -577,9 +591,9 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. code-block:: fortran ! code to start up - logical :: periodic(3) + LOGICAL :: periodic(3) ! code to initialize LAMMPS / run things / etc. - call lmp%extract_box(pflags = periodic) + CALL lmp%extract_box(pflags = periodic) -------- @@ -645,8 +659,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type USE MPI_F08 USE LIBLAMMPS - TYPE (lammps) :: lmp - TYPE (MPI_Comm) :: comm + TYPE(lammps) :: lmp + TYPE(MPI_Comm) :: comm ! ... [commands to set up LAMMPS/etc.] comm%MPI_VAL = lmp%get_mpi_comm() @@ -686,7 +700,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type LAMMPS data is type-, kind-, and rank-checked at run-time via an overloaded assignment operator. The pointers returned by this function are generally persistent; therefore it is not necessary to - call the function again, unless a :doc:`clear` command has been + call the function again unless a :doc:`clear` command has been issued, which wipes out and recreates the contents of the :cpp:class:`LAMMPS ` class. @@ -735,7 +749,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type CHARACTER(LEN=:), ALLOCATABLE :: str lmp = lammps() CALL lmp%command('units metal') - ALLOCATE ( CHARACTER(LEN=80) :: str ) + ALLOCATE(CHARACTER(LEN=80) :: str) str = lmp%extract_global('units') str = TRIM(str) ! re-allocates to length len_trim(str) here PRINT*, LEN(str), LEN_TRIM(str) @@ -899,7 +913,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type TYPE(lammps) :: lmp REAL(c_double), DIMENSION(:), POINTER :: COM ! code to setup, create atoms, etc. - CALL lmp%compute('compute COM all com') + CALL lmp%command('compute COM all com') COM = lmp%extract_compute('COM', lmp%style%global, lmp%style%type) will bind the variable *COM* to the center of mass of the atoms created in @@ -1073,10 +1087,10 @@ Procedures Bound to the :f:type:`lammps` Derived Type - Per-atom array In the case of global data, this function returns a value of type - ``real(c_double)``. For per-atom or local data, this function does not + ``REAL(c_double)``. For per-atom or local data, this function does not return a value but instead associates the pointer on the left side of the assignment to point to internal LAMMPS data. Pointers must be of the correct - data type to point to said data (i.e., ``REAL(c_double)``) and have + type and kind to point to said data (i.e., ``REAL(c_double)``) and have compatible rank. The pointer being associated with LAMMPS data is type-, kind-, and rank-checked at run-time via an overloaded assignment operator. @@ -1256,10 +1270,10 @@ Procedures Bound to the :f:type:`lammps` Derived Type of atoms, see :f:func:`gather_atoms_subset`. The *data* array will be ordered in groups of *count* values, sorted by atom - ID (e.g., if *name* is *x* and *count* = 3, then *data* = x[1][1], x[2][1], - x[3][1], x[1][2], x[2][2], x[3][2], x[1][3], :math:`\dots`); - *data* must be ``ALLOCATABLE`` and will be allocated to length - (*count* :math:`\times` *natoms*), as queried by + ID (e.g., if *name* is *x* and *count* = 3, then *data* = *x*\ (1,1), + *x*\ (2,1), *x*\ (3,1), *x*\ (1,2), *x*\ (2,2), *x*\ (3,2), *x*\ (1,3), + :math:`\dots`); *data* must be ``ALLOCATABLE`` and will be allocated to + length (*count* :math:`\times` *natoms*), as queried by :f:func:`extract_setting`. :p character(len=\*) name: desired quantity (e.g., *x* or *mask*) @@ -1340,9 +1354,10 @@ Procedures Bound to the :f:type:`lammps` Derived Type The *data* array will be in groups of *count* values, sorted by atom ID in the same order as the array *ids* (e.g., if *name* is *x*, *count* = 3, and *ids* is [100, 57, 210], then *data* might look like - [x(1,100), x(2,100), x(3,100), x(1,57), x(2,57), x(3,57), x(1,210), - :math:`\dots`]; *ids* must be provided by the user, and *data* must be - of rank 1 (i.e., ``DIMENSION(:)``) and have the ``ALLOCATABLE`` attribute. + [*x*\ (1,100), *x*\ (2,100), *x*\ (3,100), *x*\ (1,57), *x*\ (2,57), + *x*\ (3,57), *x*\ (1,210), :math:`\dots`]; *ids* must be provided by the + user, and *data* must be of rank 1 (i.e., ``DIMENSION(:)``) and have the + ``ALLOCATABLE`` attribute. :p character(len=\*) name: desired quantity (e.g., *x* or *mask*) :p integer(c_int) count: number of per-atom values you expect per atom @@ -1374,13 +1389,13 @@ Procedures Bound to the :f:type:`lammps` Derived Type The *data* array needs to be ordered in groups of *count* values, sorted by atom ID (e.g., if *name* is *x* and *count* = 3, then - *data* = [x(1,1) x(2,1) x(3,1) x(1,2) x(2,2) x(3,2) x(1,3) :math:`\dots`]; - *data* must be of length (*count* :math:`\times` *natoms*). + *data* = [x(1,1), x(2,1), x(3,1), x(1,2), x(2,2), x(3,2), x(1,3), + :math:`\dots`]; *data* must be of length (*count* :math:`\times` *natoms*). :p character(len=\*) name: quantity to be scattered (e.g., *x* or *charge*) :p data: per-atom values packed in a one-dimensional array containing the data to be scattered. This array must have length *natoms* - (e.g., for *type* or *charge*) or length *natoms*\ :math:`\times 3` + (e.g., for *type* or *charge*) or length *natoms*\ :math:`{}\times 3` (e.g., for *x* or *f*). The array *data* must be rank 1 (i.e., ``DIMENSION(:)``) and be of type ``INTEGER(c_int)`` (e.g., for *mask* or *type*) or of type ``REAL(c_double)`` (e.g., for *x* or *charge* or *f*). @@ -1405,8 +1420,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type The *data* array needs to be organized in groups of 1 or 3 values, depending on which quantity is being scattered, with the groups in the same order as the array *ids*. For example, if you want *data* to be the array - [x(1,1) x(2,1) x(3,1) x(1,100) x(2,100) x(3,100) x(1,57) x(2,57) x(3,57)], - then *ids* would be [1 100 57] and *name* would be *x*. + [x(1,1), x(2,1), x(3,1), x(1,100), x(2,100), x(3,100), x(1,57), x(2,57), + x(3,57)], then *ids* would be [1, 100, 57] and *name* would be *x*. :p character(len=\*) name: quantity to be scattered (e.g., *x* or *charge*) :p integer(c_int) ids [dimension(:)]: atom IDs corresponding to the atoms @@ -1452,10 +1467,10 @@ Procedures Bound to the :f:type:`lammps` Derived Type INTEGER :: i ! other commands to initialize LAMMPS, create bonds, etc. CALL lmp%gather_bonds(bonds) - bonds_array(1:3,1:SIZE(bonds)/3) => bonds + bonds_array(1:3, 1:SIZE(bonds)/3) => bonds DO i = 1, SIZE(bonds)/3 - WRITE(OUTPUT_UNIT,'(A,1X,I4,A,I4,1X,I4)') 'bond', bonds(1,i), & - '; type = ', bonds(2,i), bonds(3,i) + WRITE(OUTPUT_UNIT,'(A,1X,I4,A,I4,1X,I4)') 'bond', bonds_array(1,i), & + '; type = ', bonds_array(2,i), bonds_array(3,i) END DO END PROGRAM bonds @@ -1480,21 +1495,22 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p real(c_double) x [dimension(3N)]: vector of :math:`3N` x/y/z positions of the new atoms, arranged as :math:`[x_1,y_1,z_1,x_2,y_2,\dotsc]` (required/see note below) - :o integer(kind=\*) id [dimension(N)]: vector of :math:`N` atom IDs; if - absent, LAMMPS will generate them for you. \*The ``KIND`` parameter should - be ``c_int`` unless LAMMPS was compiled with ``-DLAMMPS_BIGBIG``, in which - case it should be ``c_int64_t``. - :o real(c_double) v [dimension(3N)]: vector of :math:`3N` x/y/z velocities - of the new atoms, arranged as :math:`[v_{1,x},v_{1,y},v_{1,z},v_{2,x}, - \dotsc]`; if absent, they will be set to zero - :o integer(kind=\*) image [dimension(N)]: vector of :math:`N` image flags; - if absent, they are set to zero. \*The ``KIND`` parameter should be + :o integer(kind=\*) id [dimension(N),optional]: vector of :math:`N` atom + IDs; if absent, LAMMPS will generate them for you. \*The ``KIND`` parameter + should be ``c_int`` unless LAMMPS was compiled with ``-DLAMMPS_BIGBIG``, in + which case it should be ``c_int64_t``. + :o real(c_double) v [dimension(3N),optional]: vector of :math:`3N` + *x*\ /*y*\ /*z* velocities of the new atoms, arranged as + :math:`[v_{1,x},v_{1,y},v_{1,z},v_{2,x}, \dotsc]`; if absent, they will be + set to zero + :o integer(kind=\*) image [dimension(N),optional]: vector of :math:`N` image + flags; if absent, they are set to zero. \*The ``KIND`` parameter should be ``c_int`` unless LAMMPS was compiled with ``-DLAMMPS_BIGBIG``, in which case it should be ``c_int64_t``. See note below. :o bexpand: if ``.TRUE.``, atoms outside of shrink-wrap boundaries will be created, not dropped, and the box dimensions will be extended. Default is ``.FALSE.`` - :otype bexpand: logical + :otype bexpand: logical,optional :to: :cpp:func:`lammps_create_atoms` .. note:: @@ -1543,11 +1559,11 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p character(len=\*) style: String used to search for pair style instance. :o exact: Flag to control whether style should match exactly or only a regular expression/sub-string match is applied. Default: ``.TRUE.``. - :otype exact: logical - :o integer(c_int) nsub: Match *nsub*\ th hybrid sub-style instance of - the same style. Default: 0. - :o integer(c_int) reqid: Request ID to identify the neighbor list in - case there are multiple requests from the same pair style instance. + :otype exact: logical,optional + :o integer(c_int) nsub [optional]: Match *nsub*\ th hybrid sub-style + instance of the same style. Default: 0. + :o integer(c_int) reqid [optional]: Request ID to identify the neighbor list + in case there are multiple requests from the same pair style instance. Default: 0. :to: :cpp:func:`lammps_find_pair_neighlist` :r integer(c_int) index: Neighbor list index if found, otherwise @@ -1566,8 +1582,9 @@ Procedures Bound to the :f:type:`lammps` Derived Type fixes with multiple neighbor list requests. :p character(len=\*) id: Identifier of fix instance - :o integer(c_int) reqid: request ID to identify the neighbor list in cases - in which there are multiple requests from the same fix. Default: 0. + :o integer(c_int) reqid [optional]: request ID to identify the neighbor list + in cases in which there are multiple requests from the same fix. + Default: 0. :to: :cpp:func:`lammps_find_fix_neighlist` :r index: neighbor list index if found, otherwise :math:`-1` :rtype index: integer(c_int) @@ -1585,8 +1602,9 @@ Procedures Bound to the :f:type:`lammps` Derived Type in case a compute has multiple neighbor list requests. :p character(len=\*) id: Identifier of compute instance. - :o integer(c_int) reqid: request ID to identify the neighbor list in cases - in which there are multiple requests from the same compute. Default: 0. + :o integer(c_int) reqid [optional]: request ID to identify the neighbor list + in cases in which there are multiple requests from the same compute. + Default: 0. :to: :cpp:func:`lammps_find_compute_neighlist` :r index: neighbor list index if found, otherwise :math:`-1`. :rtype index: integer(c_int) @@ -1808,9 +1826,10 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. versionadded:: 3Nov2022 - This function is analogous to the :py:func`installed_packages` function in - the Python API. The optional argument *length* sets the length of each - string in the vector *package* (default: 31). + This function is analogous to the :py:meth:`installed_packages + ` function in the Python API. + The optional argument *length* sets the length of each string in the vector + *package* (default: 31). :p character(len=:) package [dimension(:),allocatable]: list of packages; *must* have the ``ALLOCATABLE`` attribute and be of rank 1 @@ -2046,11 +2065,12 @@ Procedures Bound to the :f:type:`lammps` Derived Type necessary to combine the values of three integers representing the image flags in the :math:`x`-, :math:`y`-, and :math:`z`-directions. Unless LAMMPS is compiled with ``-DLAMMPS_BIGBIG``, those integers are limited to 10-bit - signed integers :math:`[-512,512]`. If ``-DLAMMPS_BIGBIG`` was used when + signed integers :math:`[-512,512)`. If ``-DLAMMPS_BIGBIG`` was used when compiling, then the return value is of kind ``c_int64_t`` instead of kind ``c_int``, and the valid range for the individual image flags becomes - :math:`[-1048576,1048575]` (i.e., the range of a 21-bit signed integer). - There is no check on whether the arguments conform to these requirements. + :math:`[-1048576,1048575)` (i.e., the range of a 21-bit signed integer). + There is no check on whether the arguments conform to these requirements; + values out of range will simply be wrapped back into the interval. :p integer(c_int) ix: image flag in :math:`x`-direction :p integer(c_int) iy: image flag in :math:`y`-direction @@ -2173,25 +2193,42 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. code-block:: Fortran - TYPE shield - REAL(c_double), DIMENSION(:), ALLOCATABLE :: k - ! assume k gets allocated to dimension(3,nlocal) at some point - END TYPE shield + MODULE stuff + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double, c_int64_t + IMPLICIT NONE - SUBROUTINE my_callback(caller, timestep, ids, x, fexternal) - CLASS(*), INTENT(INOUT) :: caller - INTEGER(c_int), INTENT(IN) :: timestep - INTEGER(c_int64_t), INTENT(IN) :: ids - REAL(c_double), INTENT(IN) :: x(:,:) - REAL(c_double), INTENT(OUT) :: fexternal(:,:) + TYPE shield + REAL(c_double), DIMENSION(:), ALLOCATABLE :: k + ! assume k gets allocated to dimension(3,nlocal) at some point + ! and assigned values + END TYPE shield - SELECT TYPE (caller) - TYPE IS (shield) - fexternal = - caller%k * x - CLASS DEFAULT - WRITE(error_unit,*) 'UH OH...' - END SELECT - END SUBROUTINE my_callback + SUBROUTINE my_callback(caller, timestep, ids, x, fexternal) + CLASS(*), INTENT(INOUT) :: caller + INTEGER(c_int), INTENT(IN) :: timestep + INTEGER(c_int64_t), INTENT(IN) :: ids + REAL(c_double), INTENT(IN) :: x(:,:) + REAL(c_double), INTENT(OUT) :: fexternal(:,:) + + SELECT TYPE (caller) + TYPE IS (shield) + fexternal = - caller%k * x + CLASS DEFAULT + WRITE(error_unit,*) 'UH OH...' + END SELECT + END SUBROUTINE my_callback + END MODULE stuff + + ! then, when assigning the callback function, do this: + PROGRAM example + USE LIBLAMMPS + USE stuff + TYPE(lammps) :: lmp + TYPE(shield) :: my_shield + lmp = lammps() + CALL lmp%command('fix ext all external pf/callback 1 1') + CALL lmp%set_fix_external_callback('ext', my_callback, my_shield) + END PROGRAM example -------- @@ -2235,7 +2272,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p character(len=*) id: ID of :doc:`fix external ` instance :to: :cpp:func:`lammps_fix_external_get_force` :r fexternal: pointer to the per-atom force array allocated by the fix - :rtype fexternal: real(c_double), dimension(3,nlocal) + :rtype fexternal: real(c_double),dimension(3,nlocal) -------- @@ -2414,7 +2451,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type The index in the *idx* parameter is 1-based (i.e., the first element is set with *idx*\ :math:`{} = 1`, and the last element of the vector with *idx*\ :math:`{} = N`, where :math:`N` is the value of the *length* - parameter of the call to :f:subr:`fix_external_set_vector_length`. + parameter of the call to :f:subr:`fix_external_set_vector_length`). Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an From c2a066011243bc6968f0bd18831dde464e3dae41 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Thu, 1 Dec 2022 23:49:17 -0600 Subject: [PATCH 09/19] Bug fix and unit tests for fix external-related commands --- fortran/lammps.f90 | 12 ++- unittest/fortran/test_fortran_fixexternal.f90 | 88 +++++++++++++++++++ unittest/fortran/wrap_fixexternal.cpp | 65 +++++++++++++- 3 files changed, 157 insertions(+), 8 deletions(-) diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index fc318091bf..18b8fa89a3 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -2682,8 +2682,8 @@ CONTAINS CHARACTER(LEN=*), INTENT(IN) :: id REAL(c_double), DIMENSION(:,:), TARGET, INTENT(IN) :: virial TYPE(c_ptr) :: Cid, Cvirial - TYPE(c_ptr), TARGET :: Cptr - INTEGER(c_int) :: nlocal + TYPE(c_ptr), DIMENSION(:), ALLOCATABLE, TARGET :: Cptr + INTEGER(c_int) :: nlocal, i nlocal = lmp_extract_setting(self, 'nlocal') IF (SIZE(virial,2) < nlocal .OR. SIZE(virial,1) /= 6) THEN @@ -2692,10 +2692,14 @@ CONTAINS &[Fortran/fix_external_set_energy_peratom]') END IF Cid = f2c_string(id) - Cptr = C_LOC(virial(1,1)) - Cvirial = C_LOC(Cptr) + ALLOCATE(Cptr(nlocal)) + DO i = 1, nlocal + Cptr(i) = C_LOC(virial(1,i)) + END DO + Cvirial = C_LOC(Cptr(1)) CALL lammps_fix_external_set_virial_peratom(self%handle, Cid, Cvirial) CALL lammps_free(Cid) + DEALLOCATE(Cptr) END SUBROUTINE lmp_fix_external_set_virial_peratom SUBROUTINE lmp_fix_external_set_vector_length(self, id, length) diff --git a/unittest/fortran/test_fortran_fixexternal.f90 b/unittest/fortran/test_fortran_fixexternal.f90 index ca08e3c4d9..61fea1be9d 100644 --- a/unittest/fortran/test_fortran_fixexternal.f90 +++ b/unittest/fortran/test_fortran_fixexternal.f90 @@ -4,6 +4,7 @@ MODULE ext_stuff USE LIBLAMMPS IMPLICIT NONE + INTEGER, PARAMETER :: vec_length = 8 REAL(c_double), SAVE :: direction = 1.0_c_double REAL(c_double), DIMENSION(:,:), POINTER, SAVE :: f3 => NULL(), f4 => NULL() @@ -20,21 +21,35 @@ CONTAINS REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f REAL(c_double), DIMENSION(SIZE(id)) :: e + REAL(c_double), DIMENSION(6,SIZE(id)) :: v WHERE (id == 1) f(1,:) = 1.0_c_double f(2,:) = -1.0_c_double f(3,:) = 1.25_c_double e = 1.0_c_double + v(1,:) = 1.0_c_double + v(2,:) = 2.0_c_double + v(3,:) = -1.0_c_double + v(4,:) = -2.0_c_double + v(5,:) = 3.0_c_double + v(6,:) = -3.0_c_double ELSEWHERE f(1,:) = -1.0_c_double f(2,:) = +1.0_c_double f(3,:) = -1.25_c_double e = 10.0_c_double + v(1,:) = 10.0_c_double + v(2,:) = 20.0_c_double + v(3,:) = -10.0_c_double + v(4,:) = -20.0_c_double + v(5,:) = 30.0_c_double + v(6,:) = -30.0_c_double END WHERE SELECT TYPE (instance) CLASS IS (lammps) CALL instance%fix_external_set_energy_peratom('ext1', e) + CALL instance%fix_external_set_virial_peratom('ext1', v) CLASS DEFAULT WRITE(0,*) 'UMM...this should never happen.' STOP 1 @@ -48,21 +63,35 @@ CONTAINS REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f REAL(c_double), DIMENSION(SIZE(id)) :: e + REAL(c_double), DIMENSION(6,SIZE(id)) :: v WHERE (id == 1_c_int) f(1,:) = 1.0_c_double f(2,:) = -1.0_c_double f(3,:) = 1.25_c_double e = 1.0_c_double + v(1,:) = 1.0_c_double + v(2,:) = 2.0_c_double + v(3,:) = -1.0_c_double + v(4,:) = -2.0_c_double + v(5,:) = 3.0_c_double + v(6,:) = -3.0_c_double ELSEWHERE f(1,:) = -1.0_c_double f(2,:) = +1.0_c_double f(3,:) = -1.25_c_double e = 10.0_c_double + v(1,:) = 10.0_c_double + v(2,:) = 20.0_c_double + v(3,:) = -10.0_c_double + v(4,:) = -20.0_c_double + v(5,:) = 30.0_c_double + v(6,:) = -30.0_c_double END WHERE SELECT TYPE (instance) CLASS IS (lammps) CALL instance%fix_external_set_energy_peratom('ext1', e) + CALL instance%fix_external_set_virial_peratom('ext1', v) CLASS DEFAULT WRITE(0,*) 'UMM...this should never happen.' STOP 1 @@ -76,21 +105,35 @@ CONTAINS REAL(c_double), DIMENSION(:,:), INTENT(IN) :: x REAL(c_double), DIMENSION(:,:), INTENT(OUT) :: f REAL(c_double), DIMENSION(SIZE(id)) :: e + REAL(c_double), DIMENSION(6,SIZE(id)) :: v WHERE (id == 1_c_int64_t) f(1,:) = 1.0_c_double f(2,:) = -1.0_c_double f(3,:) = 1.25_c_double e = 1.0_c_double + v(1,:) = 1.0_c_double + v(2,:) = 2.0_c_double + v(3,:) = -1.0_c_double + v(4,:) = -2.0_c_double + v(5,:) = 3.0_c_double + v(6,:) = -3.0_c_double ELSEWHERE f(1,:) = -1.0_c_double f(2,:) = +1.0_c_double f(3,:) = -1.25_c_double e = 10.0_c_double + v(1,:) = 10.0_c_double + v(2,:) = 20.0_c_double + v(3,:) = -10.0_c_double + v(4,:) = -20.0_c_double + v(5,:) = 30.0_c_double + v(6,:) = -30.0_c_double END WHERE SELECT TYPE (instance) CLASS IS (lammps) CALL instance%fix_external_set_energy_peratom('ext1', e) + CALL instance%fix_external_set_virial_peratom('ext1', v) CLASS DEFAULT WRITE(0,*) 'UMM...this should never happen.' STOP 1 @@ -201,6 +244,7 @@ END SUBROUTINE f_lammps_close SUBROUTINE f_lammps_setup_fix_external_callback() BIND(C) USE LIBLAMMPS USE keepstuff, ONLY : lmp, demo_input, cont_input, pair_input + USE ext_stuff, ONLY : vec_length IMPLICIT NONE CALL lmp%commands_list(demo_input) @@ -209,6 +253,7 @@ SUBROUTINE f_lammps_setup_fix_external_callback() BIND(C) CALL lmp%command('neigh_modify exclude group all all') CALL lmp%command('fix ext1 all external pf/callback 1 1') CALL lmp%command('fix ext2 all external pf/callback 1 1') + CALL lmp%fix_external_set_vector_length('ext2', vec_length) END SUBROUTINE f_lammps_setup_fix_external_callback SUBROUTINE f_lammps_setup_fix_external_array() BIND(C) @@ -305,6 +350,9 @@ SUBROUTINE f_lammps_find_forces() BIND(C) f3(1,:) = 4.0_c_double f3(2,:) = -4.0_c_double f3(3,:) = 6.0_c_double + f4(1,:) = 10.0_c_double + f4(2,:) = -10.0_c_double + f4(3,:) = 12.0_c_double ELSEWHERE f3(1,:) = 5.0_c_double f3(2,:) = -5.0_c_double @@ -334,3 +382,43 @@ SUBROUTINE f_lammps_set_virial() BIND(C) CALL lmp%fix_external_set_virial_global('ext4', [1.0_c_double, & 2.0_c_double, 2.5_c_double, -1.0_c_double, -2.25_c_double, -3.02_c_double]) END SUBROUTINE f_lammps_set_virial + +FUNCTION f_lammps_find_peratom_energy(i) RESULT(energy) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i + REAL(c_double) :: energy + REAL(c_double), DIMENSION(:), POINTER :: e + + e = lmp%extract_compute('peratom', lmp%style%atom, lmp%type%vector) + energy = e(i) +END FUNCTION f_lammps_Find_peratom_energy + +SUBROUTINE f_lammps_find_peratom_virial(v, i) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + REAL(c_double), DIMENSION(6) :: v + INTEGER(c_int), INTENT(IN), VALUE :: i + REAL(c_double), DIMENSION(:,:), POINTER :: virial + + virial = lmp%extract_compute('vperatom', lmp%style%atom, lmp%type%array) + v = virial(:,i) +END SUBROUTINE f_lammps_find_peratom_virial + +SUBROUTINE f_lammps_fixexternal_set_vector() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + USE ext_stuff, ONLY : vec_length + IMPLICIT NONE + REAL(c_double), DIMENSION(vec_length) :: v + INTEGER :: i + DO i = 1, vec_length + v(i) = REAL(i, c_double) + CALL lmp%fix_external_set_vector('ext2', i, v(i)) + END DO +END SUBROUTINE f_lammps_fixexternal_set_vector diff --git a/unittest/fortran/wrap_fixexternal.cpp b/unittest/fortran/wrap_fixexternal.cpp index aa6cb2b921..694aceb009 100644 --- a/unittest/fortran/wrap_fixexternal.cpp +++ b/unittest/fortran/wrap_fixexternal.cpp @@ -23,6 +23,9 @@ void f_lammps_reverse_direction(); void f_lammps_find_forces(); void f_lammps_add_energy(); void f_lammps_set_virial(); +double f_lammps_find_peratom_energy(int); +void f_lammps_find_peratom_virial(double[6], int); +void f_lammps_fixexternal_set_vector(); } using namespace LAMMPS_NS; @@ -128,10 +131,64 @@ TEST_F(LAMMPS_fixexternal, energy_peratom) { f_lammps_setup_fix_external_callback(); f_lammps_set_fix_external_callbacks(); + lammps_command(lmp, "compute peratom all pe/atom"); double energy; lammps_command(lmp, "run 0"); -/* FIXME: the per-atom energy is NOT summed up by this function! We need - another test. */ - energy = lammps_get_thermo(lmp, "pe"); - EXPECT_DOUBLE_EQ(energy, 11.0); + int nlocal = lammps_extract_setting(lmp, "nlocal"); + for (int i = 1; i <= nlocal; i++) + { + energy = f_lammps_find_peratom_energy(i); + if (i == 1) + EXPECT_DOUBLE_EQ(energy, 1.0); + else + EXPECT_DOUBLE_EQ(energy, 10.0); + } +}; + +TEST_F(LAMMPS_fixexternal, virial_peratom) +{ + f_lammps_setup_fix_external_callback(); + f_lammps_set_fix_external_callbacks(); + lammps_command(lmp, "compute vperatom all stress/atom NULL"); + double virial[6]; + lammps_command(lmp, "run 0"); + int nlocal = lammps_extract_setting(lmp, "nlocal"); + for (int i = 1; i <= nlocal; i++) + { + f_lammps_find_peratom_virial(virial, i); + if (i == 1) + { + EXPECT_DOUBLE_EQ(virial[0], -1.0); + EXPECT_DOUBLE_EQ(virial[1], -2.0); + EXPECT_DOUBLE_EQ(virial[2], 1.0); + EXPECT_DOUBLE_EQ(virial[3], 2.0); + EXPECT_DOUBLE_EQ(virial[4], -3.0); + EXPECT_DOUBLE_EQ(virial[5], 3.0); + } + else + { + EXPECT_DOUBLE_EQ(virial[0], -10.0); + EXPECT_DOUBLE_EQ(virial[1], -20.0); + EXPECT_DOUBLE_EQ(virial[2], 10.0); + EXPECT_DOUBLE_EQ(virial[3], 20.0); + EXPECT_DOUBLE_EQ(virial[4], -30.0); + EXPECT_DOUBLE_EQ(virial[5], 30.0); + } + } +}; + +TEST_F(LAMMPS_fixexternal, vector) +{ + f_lammps_setup_fix_external_callback(); + f_lammps_set_fix_external_callbacks(); + f_lammps_fixexternal_set_vector(); + lammps_command(lmp, "run 0"); + double *v; + for (int i = 0; i < 8; i++) + { + v = (double*) lammps_extract_fix(lmp, "ext2", LMP_STYLE_GLOBAL, + LMP_TYPE_VECTOR, i, 1); + EXPECT_DOUBLE_EQ(i+1, *v); + std::free(v); + } }; From b1664ce8eaf54e1e462469d8ccfb99618fc7483b Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Fri, 2 Dec 2022 00:00:57 -0600 Subject: [PATCH 10/19] replaced unit 0 with error_unit --- unittest/fortran/test_fortran_fixexternal.f90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/unittest/fortran/test_fortran_fixexternal.f90 b/unittest/fortran/test_fortran_fixexternal.f90 index 61fea1be9d..d46e1ec673 100644 --- a/unittest/fortran/test_fortran_fixexternal.f90 +++ b/unittest/fortran/test_fortran_fixexternal.f90 @@ -51,7 +51,7 @@ CONTAINS CALL instance%fix_external_set_energy_peratom('ext1', e) CALL instance%fix_external_set_virial_peratom('ext1', v) CLASS DEFAULT - WRITE(0,*) 'UMM...this should never happen.' + WRITE(error_unit,*) 'UMM...this should never happen.' STOP 1 END SELECT END SUBROUTINE f_callback_ss @@ -93,7 +93,7 @@ CONTAINS CALL instance%fix_external_set_energy_peratom('ext1', e) CALL instance%fix_external_set_virial_peratom('ext1', v) CLASS DEFAULT - WRITE(0,*) 'UMM...this should never happen.' + WRITE(error_unit,*) 'UMM...this should never happen.' STOP 1 END SELECT END SUBROUTINE f_callback_sb @@ -135,7 +135,7 @@ CONTAINS CALL instance%fix_external_set_energy_peratom('ext1', e) CALL instance%fix_external_set_virial_peratom('ext1', v) CLASS DEFAULT - WRITE(0,*) 'UMM...this should never happen.' + WRITE(error_unit,*) 'UMM...this should never happen.' STOP 1 END SELECT END SUBROUTINE f_callback_bb From 71f086e159ed22e69cd20024d76dea28e5578acd Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Fri, 2 Dec 2022 17:19:42 -0600 Subject: [PATCH 11/19] implemented scatter, gather, and friends; wrote and updated documentation --- doc/src/Fortran.rst | 209 ++++++++++++++++++++++-- fortran/lammps.f90 | 378 ++++++++++++++++++++++++++++++++++++++++++-- src/library.cpp | 209 ++++++++++++++++++++++-- 3 files changed, 764 insertions(+), 32 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 6161ccaebe..80d3447437 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -298,6 +298,16 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype scatter_atoms_subset: subroutine :f gather_bonds: :f:subr:`gather_bonds` :ftype gather_bonds: subroutine + :f gather: :f:subr:`gather` + :ftype gather: subroutine + :f gather_concat: :f:subr:`gather_concat` + :ftype gather_concat: subroutine + :f gather_subset: :f:subr:`gather_subset` + :ftype gather_subset: subroutine + :f scatter: :f:subr:`scatter` + :ftype scatter: subroutine + :f scatter_subset: :f:subr:`scatter_subset` + :ftype scatter_subset: subroutine :f create_atoms: :f:subr:`create_atoms` :ftype create_atoms: subroutine :f find_pair_neighlist: :f:func:`find_pair_neighlist` @@ -1194,8 +1204,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. f:function:: extract_variable(name[,group]) - This function calls :cpp:func:`lammps_extract_variable` and returns a scalar, - vector, or string containing the value of the variable identified by + This function calls :cpp:func:`lammps_extract_variable` and returns a + scalar, vector, or string containing the value of the variable identified by *name*. When the variable is an *equal*-style variable (or one compatible with that style such as *internal*), the variable is evaluated and the corresponding value returned. When the variable is an *atom*-style variable, @@ -1276,15 +1286,19 @@ Procedures Bound to the :f:type:`lammps` Derived Type length (*count* :math:`\times` *natoms*), as queried by :f:func:`extract_setting`. + This function is not compatible with ``-DLAMMPS_BIGBIG``. + :p character(len=\*) name: desired quantity (e.g., *x* or *mask*) :p integer(c_int) count: number of per-atom values you expect per atom (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use *count* = 3 with *image* if you want a single image flag unpacked into *x*/*y*/*z* components. - :p real(c_double) data [dimension(:),allocatable]: array into which to store + :p polymorphic data [dimension(:),allocatable]: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be - reallocated to fit the length of the incoming data. + reallocated to fit the length of the incoming data. It should have type + ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if + expecting floating-point data. :to: :cpp:func:`lammps_gather_atoms` .. note:: @@ -1324,15 +1338,19 @@ Procedures Bound to the :f:type:`lammps` Derived Type :f:func:`gather_atoms`; for a similar array but for a subset of atoms, see :f:func:`gather_atoms_subset`. + This function is not compatible with ``-DLAMMPS_BIGBIG``. + :p character(len=\*) name: desired quantity (e.g., *x* or *mask*) :p integer(c_int) count: number of per-atom values you expect per atom (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use *count* = 3 with *image* if you want a single image flag unpacked into *x*/*y*/*z* components. - :p real(c_double) data [dimension(:),allocatable]: array into which to store + :p polymorphic data [dimension(:),allocatable]: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be - reallocated to fit the length of the incoming data. + reallocated to fit the length of the incoming data. It should have type + ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if + expecting floating-point data. :to: :cpp:func:`lammps_gather_atoms_concat` -------- @@ -1346,7 +1364,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. versionadded:: 3Nov2022 This subroutine gathers data for the requested atom IDs and stores them in a - one-dimensional array allocated by the user. The data will be ordered by + one-dimensional allocatable array. The data will be ordered by atom ID, but there is no requirement that the IDs be consecutive. If you wish to return a similar array for *all* the atoms, use :f:func:`gather_atoms` or :f:func:`gather_atoms_concat`. @@ -1359,6 +1377,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type user, and *data* must be of rank 1 (i.e., ``DIMENSION(:)``) and have the ``ALLOCATABLE`` attribute. + This function is not compatible with ``-DLAMMPS_BIGBIG``. + :p character(len=\*) name: desired quantity (e.g., *x* or *mask*) :p integer(c_int) count: number of per-atom values you expect per atom (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use @@ -1366,10 +1386,12 @@ Procedures Bound to the :f:type:`lammps` Derived Type *x*/*y*/*z* components. :p integer(c_int) ids [dimension(:)]: atom IDs corresponding to the atoms to be gathered - :p real(c_double) data [dimension(:),allocatable]: array into which to store + :p polymorphic data [dimension(:),allocatable]: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be - reallocated to fit the length of the incoming data. + reallocated to fit the length of the incoming data. It should have type + ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if + expecting floating-point data. :to: :cpp:func:`lammps_gather_atoms_subset` -------- @@ -1443,7 +1465,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. versionadded:: 3Nov2022 - This function copies the list of all bonds into an allocated array. + This function copies the list of all bonds into an allocatable array. The array will be filled with (bond type, bond atom 1, bond atom 2) for each bond. The array is allocated to the right length (i.e., three times the number of bonds). The array *data* must be of the same type as the LAMMPS @@ -1482,6 +1504,173 @@ Procedures Bound to the :f:type:`lammps` Derived Type -------- +.. f:subroutine:: gather(self, name, count, data) + + Gather the named per-atom, per-atom fix, per-atom compute, or fix + property/atom-based entities from all processes, in order by atom ID. + + .. versionadded:: TBD + + This subroutine gathers data from all processes and stores them in a + one-dimensional allocatable array. The array *data* will be + ordered by atom ID, which requires consecutive IDs (1 to *natoms*\ ). If you + need a similar array but for non-consecutive atom IDs, see + :cpp:func:`lammps_gather_concat`; for a similar array but for a subset of + atoms, see :cpp:func:`lammps_gather_subset`. + + The *data* array will be ordered in groups of *count* values, sorted by atom + ID (e.g., if *name* is *x*, then *data* is [x(1,1), x(2,1), x(3,1), x(1,2), + x(2,2), x(3,2), x(1,3), :math:`\dots`]); *data* must be ``ALLOCATABLE`` and + will be allocated to length (*count*\ :math:`{}\times{}`\ *natoms*), as + queried by :f:func:`extract_setting`. + + This function will return an error if fix or compute data are requested and + the fix or compute ID given does not have per-atom data. See the note about + re-interpreting the vector as a matrix at :f:subr:`gather_atoms`. + + This function is not compatible with ``-DLAMMPS_BIGBIG``. + + :p character(len=\*) name: desired quantity (e.g., "x" or "mask" for atom + properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, + "d_name" or "i_name" for fix property/atom vectors with *count* = 1, + "d2_name" or "i2_name" for fix propert/atom vectors with + *count*\ :math:`{}> 1`) + :p integer(c_int) count: number of per-atom values (e.g., 1 for *type* or + *charge*, 3 for *x* or *f*); use *count* = 3 with *image* if you want the + image flags unpacked into (*x*,\ *y*,\ *z*) components. + :p real(c_double) data [dimension(:),allocatable]: array into which to store + the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 + (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be + reallocated to fit the length of the incoming data. + :to: :cpp:func:`lammps_gather` + +-------- + +.. f:subroutine:: gather_concat(self, name, count, data) + + Gather the named per-atom, per-atom fix, per-atom compute, or fix + property/atom-based entities from all processes, unordered. + + This subroutine gathers data for all atoms and stores them in a + one-dimensional allocatable array. The data will be a + concatenation of chunks from each processor's owned atoms, in whatever order + the atoms are in on each processor. This process has no requirement that the + atom IDs be consecutive. If you need the ID of each atom, you can do another + call to either :f:subr:`gather_atoms_concat` or :f:subr:`gather_concat` with + *name* set to ``id``. If you have consecutive IDs and want the data to be in + order, use :f:subr:`gather`; for a similar array but for a subset of + atoms, use :f:subr:`gather_subset`. + + The *data* array will be in groups of *count* values, with *natoms* groups + total, but not in order by atom ID (e.g., if *name* is *x* and *count* is 3, + then *data* might be something like [x(1,11), x(2,11), x(3,11), x(1,3), + x(2,3), x(3,3), x(1,5), :math:`\dots`]); *data* must be ``ALLOCATABLE`` and + will be allocated to length (*count* :math:`\times` *natoms*), as queried by + :f:func:`extract_setting`. + + This function is not compatible with ``-DLAMMPS_BIGBIG``. + + :p character(len=\*) name: desired quantity (e.g., "x" or "mask" for atom + properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, + "d_name" or "i_name" for fix property/atom vectors with *count* = 1, + "d2_name" or "i2_name" for fix propert/atom vectors with + *count*\ :math:`{}> 1`) + :p integer(c_int) count: number of per-atom values you expect per atom + (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use + *count* = 3 with *image* if you want a single image flag unpacked into + *x*/*y*/*z* components. + :p polymorphic data [dimension(:),allocatable]: array into which to store + the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 + (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be + reallocated to fit the length of the incoming data. It should have type + ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if + expecting floating-point data. + :to: :cpp:func:`lammps_gather_concat` + +-------- + +.. f:subroutine:: gather_subset(name, count, ids, data) + + Gather the named per-atom, per-atom fix, per-atom compute, or fix + property/atom-based entities from all processes for a subset of atoms. + + This subroutine gathers data for the requested atom IDs and stores them in a + one-dimensional allocatable array. The data will be ordered by atom ID, but + there is no requirement that the IDs be consecutive. If you wish to return a + similar array for *all* the atoms, use :f:subr:`gather` or + :f:subr:`gather_concat`. + + The *data* array will be in groups of *count* values, sorted by atom ID in + the same order as the array *ids* (e.g., if *name* is *x*, *count* = 3, and + *ids* is [100, 57, 210], then *data* might look like [x(1,100), x(2,100), + x(3,100), x(1,57), x(2,57), x(3,57), x(1,210), :math:`\dots`]); *ids* must + be provided by the user, and *data* must have the ``ALLOCATABLE`` attribute + and be of rank 1 (i.e., ``DIMENSION(:)``). If *data* is already allocated, + it will be reallocated to fit the length of the incoming data. + + This function is not compatible with ``-DLAMMPS_BIGBIG``. + + :p character(len=\*) name: quantity to be scattered + + :p integer(c_int) ids [dimension(:)]: atom IDs corresponding to the atoms + being scattered (e.g., "x" or "f" for atom properties, "f_id" for per-atom + fix data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix + property/atom vectors with *count* = 1, "d2_name" or "i2_name" for fix + property/atom vectors with *count*\ :math:`{} > 1`) + :p data: per-atom values packed into a one-dimensional array containing the + data to be scattered. This array must have either the same length as *ids* + (for *mask*, *type*, etc.) or three times its length (for *x*, *f*, etc.); + the array must be rank 1 and be of type ``INTEGER(c_int)`` (e.g., for + *mask* or *type*) or of type ``REAL(c_double)`` (e.g., for *charge*, *x*, + or *f*). + :ptype data: polymorphic,dimension(:) + :to: :cpp:func:`lammps_scatter_subset` + +-------- + +.. f:subroutine:: scatter(name, data) + + This function calls :cpp:func:`lammps_scatter` to scatter the named + atom-based entities in *data* to all processes. + + *This function is not yet documented, as the underlying C routine has not + been thoroughly tested yet.* + +-------- + +.. f:subroutine:: scatter_subset(name, ids, data) + + This function calls :cpp:func:`lammps_scatter_subset` to scatter the named + per-atom, per-atom fix, per-atom compute, or fix property/atom-based + entities in *data* from a subset of atoms to all processes. + + .. versionadded:: TBD + + This subroutine takes data stored in a one-dimensional array supplied by the + user and scatters them to a subset of atoms on all processes. The array + *data* contains data associated with atom IDs, but there is no requirement + that the IDs be consecutive, as they are provided in a separate array. + Use :f:subr:`scatter` to scatter data for all atoms, in order. + + The *data* array needs to be organized in groups of *count* values, with the + groups in the same order as the array *ids*. For example, if you want *data* + to be the array [x(1,1), x(2,1), x(3,1), x(1,100), x(2,100), x(3,100), + x(1,57), x(2,57), x(3,57)], then *count* = 3 and *ids* = [1, 100, 57]. + + This function is not compatible with ``-DLAMMPS_BIGBIG``. + + :p character(len=\*) name: desired quantity (e.g., "x" or "mask" for atom + properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, + "d_name" or "i_name" for fix property/atom vectors with *count* = 1, + "d2_name" or "i2_name" for fix propert/atom vectors with + *count*\ :math:`{}> 1`) + :p integer(c_int) ids: list of atom IDs to scatter data for + :p polymorphic data [dimension(:)]: per-atom values packed in a + one-dimensional array of length *size(ids)* \* *count*. + :to: :cpp:func:`lammps_scatter_subset` + +-------- + .. f:subroutine:: create_atoms([id,] type, x, [v,] [image,] [bexpand]) This method calls :cpp:func:`lammps_create_atoms` to create additional atoms diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index 18b8fa89a3..4cf90e22a2 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -135,7 +135,24 @@ MODULE LIBLAMMPS PROCEDURE, PRIVATE :: lmp_gather_bonds_big GENERIC :: gather_bonds => lmp_gather_bonds_small, & lmp_gather_bonds_big -! + PROCEDURE, PRIVATE :: lmp_gather_int + PROCEDURE, PRIVATE :: lmp_gather_double + GENERIC :: gather => lmp_gather_int, lmp_gather_double + PROCEDURE, PRIVATE :: lmp_gather_concat_int + PROCEDURE, PRIVATE :: lmp_gather_concat_double + GENERIC :: gather_concat => lmp_gather_concat_int, & + lmp_gather_concat_double + PROCEDURE, PRIVATE :: lmp_gather_subset_int + PROCEDURE, PRIVATE :: lmp_gather_subset_double + GENERIC :: gather_subset => lmp_gather_subset_int, & + lmp_gather_subset_double + PROCEDURE, PRIVATE :: lmp_scatter_int + PROCEDURE, PRIVATE :: lmp_scatter_double + GENERIC :: scatter => lmp_scatter_int, lmp_scatter_double + PROCEDURE, PRIVATE :: lmp_scatter_subset_int + PROCEDURE, PRIVATE :: lmp_scatter_subset_double + GENERIC :: scatter_subset => lmp_scatter_subset_int, & + lmp_scatter_subset_double PROCEDURE, PRIVATE :: lmp_create_atoms_int PROCEDURE, PRIVATE :: lmp_create_atoms_bigbig GENERIC :: create_atoms => lmp_create_atoms_int, & @@ -552,13 +569,42 @@ MODULE LIBLAMMPS TYPE(c_ptr), VALUE :: handle, data END SUBROUTINE lammps_gather_bonds - !SUBROUTINE lammps_gather + SUBROUTINE lammps_gather(handle, name, type, count, data) BIND(C) + IMPORT :: c_ptr, c_int + IMPLICIT NONE + TYPE(c_ptr), VALUE :: handle, name, data + INTEGER(c_int), VALUE :: type, count + END SUBROUTINE lammps_gather - !SUBROUTINE lammps_gather_concat + SUBROUTINE lammps_gather_concat(handle, name, type, count, data) BIND(C) + IMPORT :: c_ptr, c_int + IMPLICIT NONE + TYPE(c_ptr), VALUE :: handle, name, data + INTEGER(c_int), VALUE :: type, count + END SUBROUTINE lammps_gather_concat - !SUBROUTINE lammps_gather_subset + SUBROUTINE lammps_gather_subset(handle, name, type, count, ndata, ids, & + data) BIND(C) + IMPORT :: c_ptr, c_int + IMPLICIT NONE + TYPE(c_ptr), VALUE :: handle, name, ids, data + INTEGER(c_int), VALUE :: type, count, ndata + END SUBROUTINE lammps_gather_subset - !SUBROUTINE lammps_scatter_subset + SUBROUTINE lammps_scatter(handle, name, type, count, data) BIND(C) + IMPORT :: c_ptr, c_int + IMPLICIT NONE + TYPE(c_ptr), VALUE :: handle, name, data + INTEGER(c_int), VALUE :: type, count + END SUBROUTINE lammps_scatter + + SUBROUTINE lammps_scatter_subset(handle, name, type, count, ndata, ids, & + data) BIND(C) + IMPORT :: c_ptr, c_int + IMPLICIT NONE + TYPE(c_ptr), VALUE :: handle, name, ids, data + INTEGER(c_int), VALUE :: count, ndata, type + END SUBROUTINE lammps_scatter_subset FUNCTION lammps_create_atoms(handle, n, id, type, x, v, image, bexpand) & BIND(C) @@ -1622,7 +1668,7 @@ CONTAINS IF (count /= 1 .AND. count /= 3) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & 'gather_atoms_subset requires "count" to be 1 or 3 & - &[Fortran/gather_atoms]') + &[Fortran/gather_atoms_subset]') END IF ndata = SIZE(ids, KIND=c_int) @@ -1652,7 +1698,7 @@ CONTAINS IF (count /= 1 .AND. count /= 3) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & 'gather_atoms_subset requires "count" to be 1 or 3 & - &[Fortran/gather_atoms]') + &[Fortran/gather_atoms_subset]') END IF ndata = SIZE(ids, KIND=c_int) @@ -1746,7 +1792,8 @@ CONTAINS Ccount = SIZE(data, KIND=c_int) / Cndata IF (Ccount /= 1 .AND. Ccount /= 3) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & - 'scatter_atoms_subset requires either 1 or 3 data per atom') + 'scatter_atoms_subset requires either 1 or 3 data per atom & + &[Fortran/scatter_atoms_subset]') END IF Cname = f2c_string(name) @@ -1771,7 +1818,8 @@ CONTAINS Ccount = SIZE(data, KIND=c_int) / Cndata IF (Ccount /= 1 .AND. Ccount /= 3) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & - 'scatter_atoms_subset requires either 1 or 3 data per atom') + 'scatter_atoms_subset requires either 1 or 3 data per atom & + &[Fortran/scatter_atoms_subset]') END IF Cname = f2c_string(name) @@ -1792,7 +1840,7 @@ CONTAINS IF (SIZE_TAGINT /= 4_c_int) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & - 'Incompatible integer kind in gather_bonds [Fortran API]') + 'Incompatible integer kind in gather_bonds [Fortran/gather_bonds]') END IF IF (ALLOCATED(data)) DEALLOCATE(data) IF (SIZE_BIGINT == 4_c_int) THEN @@ -1815,7 +1863,7 @@ CONTAINS IF (SIZE_TAGINT /= 8_c_int) THEN CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & - 'Incompatible integer kind in gather_bonds [Fortran API]') + 'Incompatible integer kind in gather_bonds [Fortran/gather_bonds]') END IF nbonds = lmp_extract_global(self, 'nbonds') IF (ALLOCATED(data)) DEALLOCATE(data) @@ -1824,6 +1872,314 @@ CONTAINS CALL lammps_gather_bonds(self%handle, Cdata) END SUBROUTINE lmp_gather_bonds_big + ! equivalent function to lammps_gather (for int data) + SUBROUTINE lmp_gather_int(self, name, count, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), INTENT(IN) :: count + INTEGER(c_int), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data + TYPE(c_ptr) :: Cdata, Cname + INTEGER(c_int) :: natoms + INTEGER(c_int), PARAMETER :: Ctype = 0_c_int + REAL(c_double) :: dnatoms + CHARACTER(LEN=100) :: error_msg + + IF (count /= 1 .AND. count /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'gather requires "count" to be 1 or 3 [Fortran/gather]') + END IF + + dnatoms = lmp_get_natoms(self) + IF (dnatoms > HUGE(1_c_int)) THEN + WRITE(error_msg,'(A,1X,I0,1X,A)') & + 'Cannot use library function gather with more than', & + HUGE(0_c_int), 'atoms [Fortran/gather]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) + END IF + natoms = NINT(dnatoms, c_int) + + Cname = f2c_string(name) + IF (ALLOCATED(data)) DEALLOCATE(data) + ALLOCATE(data(natoms*count)) + Cdata = C_LOC(data(1)) + CALL lammps_gather(self%handle, Cname, Ctype, count, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_gather_int + + ! equivalent function to lammps_gather_atoms (for doubles) + SUBROUTINE lmp_gather_double(self, name, count, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), INTENT(IN) :: count + REAL(c_double), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data + TYPE(c_ptr) :: Cdata, Cname + INTEGER(c_int) :: natoms + INTEGER(c_int), PARAMETER :: Ctype = 1_c_int + REAL(c_double) :: dnatoms + CHARACTER(LEN=100) :: error_msg + + IF (count /= 1 .AND. count /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'gather requires "count" to be 1 or 3 [Fortran/gather]') + END IF + + dnatoms = lmp_get_natoms(self) + IF (dnatoms > HUGE(1_c_int)) THEN + WRITE(error_msg,'(A,1X,I0,1X,A)') & + 'Cannot use library function gather with more than', & + HUGE(0_c_int), 'atoms [Fortran/gather]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) + END IF + natoms = NINT(dnatoms, c_int) + + Cname = f2c_string(name) + IF (ALLOCATED(data)) DEALLOCATE(data) + ALLOCATE(data(natoms*count)) + Cdata = C_LOC(data(1)) + CALL lammps_gather(self%handle, Cname, Ctype, count, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_gather_double + + ! equivalent function to lammps_gather_concat (for ints) + SUBROUTINE lmp_gather_concat_int(self, name, count, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), INTENT(IN) :: count + INTEGER(c_int), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data + TYPE(c_ptr) :: Cdata, Cname + INTEGER(c_int) :: natoms + INTEGER(c_int), PARAMETER :: Ctype = 0_c_int + REAL(c_double) :: dnatoms + CHARACTER(LEN=100) :: error_msg + + IF (count /= 1 .AND. count /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'gather_concat requires "count" to be 1 or 3 [Fortran/gather_concat]') + END IF + + dnatoms = lmp_get_natoms(self) + IF (dnatoms > HUGE(1_c_int)) THEN + WRITE(error_msg,'(A,1X,I0,1X,A)') & + 'Cannot use library function gather_concat with more than', & + HUGE(0_c_int), 'atoms [Fortran/gather_concat]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) + END IF + natoms = NINT(dnatoms, c_int) + + Cname = f2c_string(name) + IF (ALLOCATED(data)) DEALLOCATE(data) + ALLOCATE(data(natoms*count)) + Cdata = C_LOC(data(1)) + CALL lammps_gather_concat(self%handle, Cname, Ctype, count, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_gather_concat_int + + ! equivalent function to lammps_gather_concat (for doubles) + SUBROUTINE lmp_gather_concat_double(self, name, count, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), INTENT(IN) :: count + REAL(c_double), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data + TYPE(c_ptr) :: Cdata, Cname + INTEGER(c_int) :: natoms + INTEGER(c_int), PARAMETER :: Ctype = 1_c_int + REAL(c_double) :: dnatoms + CHARACTER(LEN=100) :: error_msg + + IF (count /= 1 .AND. count /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'gather_concat requires "count" to be 1 or 3 [Fortran/gather_concat]') + END IF + + dnatoms = lmp_get_natoms(self) + IF (dnatoms > HUGE(1_c_int)) THEN + WRITE(error_msg,'(A,1X,I0,1X,A)') & + 'Cannot use library function gather_concat with more than', & + HUGE(0_c_int), 'atoms [Fortran/gather_concat]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) + END IF + natoms = NINT(dnatoms, c_int) + + Cname = f2c_string(name) + IF (ALLOCATED(data)) DEALLOCATE(data) + ALLOCATE(data(natoms*count)) + Cdata = C_LOC(data(1)) + CALL lammps_gather_concat(self%handle, Cname, Ctype, count, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_gather_concat_double + + ! equivalent function to lammps_gather_subset (for integers) + SUBROUTINE lmp_gather_subset_int(self, name, count, ids, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), INTENT(IN) :: count + INTEGER(c_int), DIMENSION(:), TARGET, INTENT(IN) :: ids + INTEGER(c_int), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data + INTEGER(c_int) :: ndata + TYPE(c_ptr) :: Cdata, Cname, Cids + INTEGER(c_int), PARAMETER :: Ctype = 0_c_int + + IF (count /= 1 .AND. count /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'gather_subset requires "count" to be 1 or 3 [Fortran/gather_subset]') + END IF + + ndata = SIZE(ids, KIND=c_int) + + Cname = f2c_string(name) + IF (ALLOCATED(data)) DEALLOCATE(data) + ALLOCATE(data(ndata*count)) + data = -1_c_int + Cdata = C_LOC(data(1)) + Cids = C_LOC(ids(1)) + CALL lammps_gather_subset(self%handle, Cname, Ctype, count, & + ndata, Cids, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_gather_subset_int + + ! equivalent function to lammps_gather_subset (for doubles) + SUBROUTINE lmp_gather_subset_double(self, name, count, ids, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), INTENT(IN) :: count + INTEGER(c_int), DIMENSION(:), TARGET, INTENT(IN) :: ids + REAL(c_double), DIMENSION(:), ALLOCATABLE, TARGET, INTENT(OUT) :: data + INTEGER(c_int) :: ndata + TYPE(c_ptr) :: Cdata, Cname, Cids + INTEGER(c_int), PARAMETER :: Ctype = 1_c_int + + IF (count /= 1 .AND. count /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'gather_subset requires "count" to be 1 or 3 [Fortran/gather_subset]') + END IF + + ndata = SIZE(ids, KIND=c_int) + + Cname = f2c_string(name) + IF (ALLOCATED(data)) DEALLOCATE(data) + ALLOCATE(data(ndata*count)) + Cdata = C_LOC(data(1)) + Cids = C_LOC(ids(1)) + CALL lammps_gather_subset(self%handle, Cname, Ctype, count, & + ndata, Cids, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_gather_subset_double + + ! equivalent function to lammps_scatter (for integers) + SUBROUTINE lmp_scatter_int(self, name, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), DIMENSION(:), TARGET :: data + INTEGER(c_int) :: natoms, Ccount + INTEGER(c_int), PARAMETER :: Ctype = 0_c_int + TYPE(c_ptr) :: Cname, Cdata + REAL(c_double) :: dnatoms + CHARACTER(LEN=100) :: error_msg + + dnatoms = lmp_get_natoms(self) + IF (dnatoms > HUGE(1_c_int)) THEN + WRITE(error_msg,'(A,1X,I0,1X,A)') & + 'Cannot use library function scatter with more than', & + HUGE(0_c_int), 'atoms [Fortran/scatter]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) + END IF + natoms = NINT(dnatoms, c_int) + + Cname = f2c_string(name) + Cdata = C_LOC(data(1)) + Ccount = SIZE(data) / natoms + + IF (Ccount /= 1 .AND. Ccount /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'lammps_scatter requires either 1 or 3 data per atom') + END IF + CALL lammps_scatter(self%handle, Cname, Ctype, Ccount, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_scatter_int + + ! equivalent function to lammps_scatter (for doubles) + SUBROUTINE lmp_scatter_atoms_double(self, name, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + REAL(c_double), DIMENSION(:), TARGET :: data + INTEGER(c_int) :: natoms, Ccount + INTEGER(c_int), PARAMETER :: Ctype = 1_c_int + TYPE(c_ptr) :: Cname, Cdata + REAL(c_double) :: dnatoms + CHARACTER(LEN=100) :: error_msg + + dnatoms = lmp_get_natoms(self) + IF (dnatoms > HUGE(1_c_int)) THEN + WRITE(error_msg,'(A,1X,I0,1X,A)') & + 'Cannot use library function scatter with more than', & + HUGE(0_c_int), 'atoms [Fortran/scatter]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) + END IF + natoms = NINT(dnatoms, c_int) + + Cname = f2c_string(name) + Cdata = C_LOC(data(1)) + Ccount = SIZE(data) / natoms + + IF (Ccount /= 1 .AND. Ccount /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'scatter requires either 1 or 3 data per atom [Fortran/scatter]') + END IF + CALL lammps_scatter(self%handle, Cname, Ctype, Ccount, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_scatter_double + + ! equivalent function to lammps_scatter_subset (for integers) + SUBROUTINE lmp_scatter_subset_int(self, name, ids, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), DIMENSION(:), TARGET :: ids + INTEGER(c_int), DIMENSION(:), TARGET :: data + INTEGER(c_int), PARAMETER :: Ctype = 0_c_int + INTEGER(c_int) :: Cndata, Ccount + TYPE(c_ptr) :: Cdata, Cname, Cids + + Cndata = SIZE(ids, KIND=c_int) + Ccount = SIZE(data, KIND=c_int) / Cndata + IF (Ccount /= 1 .AND. Ccount /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'scatter_subset requires either 1 or 3 data per atom & + &[Fortran/scatter_subset]') + END IF + + Cname = f2c_string(name) + Cdata = C_LOC(data(1)) + Cids = C_LOC(ids(1)) + CALL lammps_scatter_subset(self%handle, Cname, Ctype, Ccount, & + Cndata, Cids, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_scatter_subset_int + + ! equivalent function to lammps_scatter_subset (for doubles) + SUBROUTINE lmp_scatter_subset_double(self, name, ids, data) + CLASS(lammps), INTENT(IN) :: self + CHARACTER(LEN=*), INTENT(IN) :: name + INTEGER(c_int), DIMENSION(:), TARGET :: ids + REAL(c_double), DIMENSION(:), TARGET :: data + INTEGER(c_int), PARAMETER :: Ctype = 1_c_int + INTEGER(c_int) :: Cndata, Ccount + TYPE(c_ptr) :: Cdata, Cname, Cids + + Cndata = SIZE(ids, KIND=c_int) + Ccount = SIZE(data, KIND=c_int) / Cndata + IF (Ccount /= 1 .AND. Ccount /= 3) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'scatter_subset requires either 1 or 3 data per atom') + END IF + + Cname = f2c_string(name) + Cdata = C_LOC(data(1)) + Cids = C_LOC(ids(1)) + CALL lammps_scatter_subset(self%handle, Cname, Ctype, Ccount, & + Cndata, Cids, Cdata) + CALL lammps_free(Cname) + END SUBROUTINE lmp_scatter_subset_double + ! equivalent function to lammps_create_atoms (int ids or id absent) SUBROUTINE lmp_create_atoms_int(self, id, type, x, v, image, bexpand) CLASS(lammps), INTENT(IN) :: self diff --git a/src/library.cpp b/src/library.cpp index 1e3c96cdcf..9dc9a701f3 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -685,7 +685,7 @@ void lammps_commands_string(void *handle, const char *str) \verbatim embed:rst This number may be very large when running large simulations across -multiple processors. Depending on compile time choices, LAMMPS may be +multiple processes. Depending on compile time choices, LAMMPS may be using either 32-bit or a 64-bit integer to store this number. For portability this function returns thus a double precision floating point number, which can represent up to a 53-bit signed @@ -2263,7 +2263,7 @@ int lammps_set_variable(void *handle, char *name, char *str) // Library functions for scatter/gather operations of data // ---------------------------------------------------------------------- -/** Gather the named atom-based entity for all atoms across all processors, +/** Gather the named atom-based entity for all atoms across all processes, * in order. * \verbatim embed:rst @@ -2282,6 +2282,8 @@ x[0][2], x[1][0], x[1][1], x[1][2], x[2][0], :math:`\dots`); *natoms*), as queried by :cpp:func:`lammps_get_natoms`, :cpp:func:`lammps_extract_global`, or :cpp:func:`lammps_extract_setting`. +This function is not compatible with ``-DLAMMPS_BIGBIG``. + \endverbatim * * \param handle pointer to a previously created LAMMPS instance @@ -2415,7 +2417,7 @@ void lammps_gather_atoms(void *handle, char *name, int type, int count, void *da END_CAPTURE } -/** Gather the named atom-based entity for all atoms across all processors, +/** Gather the named atom-based entity for all atoms across all processes, * unordered. * \verbatim embed:rst @@ -2432,12 +2434,14 @@ of atoms, use :cpp:func:`lammps_gather_atoms_subset`. The *data* array will be in groups of *count* values, with *natoms* groups total, but not in order by atom ID (e.g., if *name* is *x* and *count* -is 3, then *data* might be something like = x[10][0], x[10][1], x[10][2], +is 3, then *data* might be something like x[10][0], x[10][1], x[10][2], x[2][0], x[2][1], x[2][2], x[4][0], :math:`\dots`); *data* must be pre-allocated by the caller to length (*count* :math:`\times` *natoms*), as queried by :cpp:func:`lammps_get_natoms`, :cpp:func:`lammps_extract_global`, or :cpp:func:`lammps_extract_setting`. +This function is not compatible with ``-DLAMMPS_BIGBIG``. + \endverbatim * * \param handle: pointer to a previously created LAMMPS instance @@ -2599,6 +2603,8 @@ x[100][2], x[57][0], x[57][1], x[57][2], x[210][0], :math:`\dots`); *data* must be pre-allocated by the caller to length (*count* :math:`\times` *ndata*). +This function is not compatible with ``-DLAMMPS_BIGBIG``. + \endverbatim * * \param handle: pointer to a previously created LAMMPS instance @@ -2745,20 +2751,22 @@ void lammps_gather_atoms_subset(void *handle, char *name, int type, int count, END_CAPTURE } -/** Scatter the named atom-based entities in *data* to all processors. +/** Scatter the named atom-based entities in *data* to all processes. * \verbatim embed:rst This subroutine takes data stored in a one-dimensional array supplied by the -user and scatters them to all atoms on all processors. The data must be +user and scatters them to all atoms on all processes. The data must be ordered by atom ID, with the requirement that the IDs be consecutive. Use :cpp:func:`lammps_scatter_atoms_subset` to scatter data for some (or all) atoms, unordered. The *data* array needs to be ordered in groups of *count* values, sorted by atom ID (e.g., if *name* is *x* and *count* = 3, then -*data* = x[0][0], x[0][1], x[0][2], x[1][0], x[1][1], x[1][2], x[2][0], -:math:`\dots`); *data* must be of length (*count* :math:`\times` *natoms*). +*data* = {x[0][0], x[0][1], x[0][2], x[1][0], x[1][1], x[1][2], x[2][0], +:math:`\dots`}); *data* must be of length (*count* :math:`\times` *natoms*). + +This function is not compatible with ``-DLAMMPS_BIGBIG``. \endverbatim * @@ -2768,7 +2776,7 @@ atom ID (e.g., if *name* is *x* and *count* = 3, then * \param count number of per-atom values (e.g., 1 for *type* or *charge*, * 3 for *x* or *f*); use *count* = 3 with *image* if you have * a single image flag packed into (*x*,*y*,*z*) components. - * \param data per-atom values packed in a 1-dimensional array of length + * \param data per-atom values packed in a one-dimensional array of length * *natoms* \* *count*. * */ @@ -2879,12 +2887,12 @@ void lammps_scatter_atoms(void *handle, char *name, int type, int count, void *d } /** Scatter the named atom-based entities in *data* from a subset of atoms - * to all processors. + * to all processes. * \verbatim embed:rst This subroutine takes data stored in a one-dimensional array supplied by the -user and scatters them to a subset of atoms on all processors. The array +user and scatters them to a subset of atoms on all processes. The array *data* contains data associated with atom IDs, but there is no requirement that the IDs be consecutive, as they are provided in a separate array. Use :cpp:func:`lammps_scatter_atoms` to scatter data for all atoms, in order. @@ -2895,6 +2903,8 @@ to be the array {x[1][0], x[1][1], x[1][2], x[100][0], x[100][1], x[100][2], x[57][0], x[57][1], x[57][2]}, then *count* = 3, *ndata* = 3, and *ids* would be {1, 100, 57}. +This function is not compatible with ``-DLAMMPS_BIGBIG``. + \endverbatim * * \param handle: pointer to a previously created LAMMPS instance @@ -3146,6 +3156,43 @@ void lammps_gather_bonds(void *handle, void *data) END_CAPTURE } +/** Gather the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based entities + * from all processes, in order by atom ID. + * +\verbatim embed:rst + +This subroutine gathers data from all processes and stores them in a one-dimensional array +allocated by the user. The array *data* will be ordered by atom ID, which requires consecutive IDs +(1 to *natoms*\ ). If you need a similar array but for non-consecutive atom IDs, see +:cpp:func:`lammps_gather_concat`; for a similar array but for a subset of atoms, see +:cpp:func:`lammps_gather_subset`. + +The *data* array will be ordered in groups of *count* values, sorted by atom ID (e.g., if *name* is +*x*, then *data* is {x[0][0], x[0][1], x[0][2], x[1][0], x[1][1], x[1][2], x[2][0], +:math:`\dots`}); *data* must be pre-allocated by the caller to the correct length +(*count*\ :math:`{}\times{}`\ *natoms*), as queried by :cpp:func:`lammps_get_natoms`, +:cpp:func:`lammps_extract_global`, or :cpp:func:`lammps_extract_setting`. + +This function will return an error if fix or compute data are requested and the fix or compute ID +given does not have per-atom data. + +This function is not compatible with ``-DLAMMPS_BIGBIG``. + +\endverbatim + * + * \param handle pointer to a previously created LAMMPS instance + * \param name desired quantity (e.g., "x" or "f" for atom properties, "f_id" for per-atom fix + * data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix + * property/atom vectors with *count* = 1, "d2_name" or "i2_name" for fix + * property/atom vectors with *count* > 1) + * \param type 0 for ``int`` values, 1 for ``double`` values + * \param count number of per-atom values (e.g., 1 for *type* or *charge*, 3 for *x* or *f*); + * use *count* = 3 with *image* if you want the image flags unpacked into + * (*x*,*y*,*z*) components. + * \param data per-atom values packed into a one-dimensional array of length + * *natoms* \* *count*. + * + */ /* ---------------------------------------------------------------------- Contributing author: Thomas Swinburne (CNRS & CINaM, Marseille, France) gather the named atom-based entity for all atoms @@ -3381,6 +3428,44 @@ void lammps_gather(void *handle, char *name, int type, int count, void *data) END_CAPTURE } +/** Gather the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based entities + * from all processes, unordered. + * +\verbatim embed:rst + +This subroutine gathers data for all atoms and stores them in a one-dimensional array allocated by +the user. The data will be a concatenation of chunks from each processor's owned atoms, in +whatever order the atoms are in on each processor. This process has no requirement that the atom +IDs be consecutive. If you need the ID of each atom, you can do another call to either +:cpp:func:`lammps_gather_atoms_concat` or :cpp:func:`lammps_gather_concat` with *name* set to +``id``. If you have consecutive IDs and want the data to be in order, use +:cpp:func:`lammps_gather`; for a similar array but for a subset of atoms, use +:cpp:func:`lammps_gather_subset`. + +The *data* array will be in groups of *count* values, with *natoms* groups total, but not in order +by atom ID (e.g., if *name* is *x* and *count* is 3, then *data* might be something like +{x[10][0], x[10][1], x[10][2], x[2][0], x[2][1], x[2][2], x[4][0], :math:`\dots`}); *data* must be +pre-allocated by the caller to length (*count* :math:`\times` *natoms*), as queried by +:cpp:func:`lammps_get_natoms`, :cpp:func:`lammps_extract_global`, or +:cpp:func:`lammps_extract_setting`. + +This function is not compatible with ``-DLAMMPS_BIGBIG``. + +\endverbatim + * + * \param handle: pointer to a previously created LAMMPS instance + * \param name: desired quantity (e.g., "x" or "f" for atom properties, "f_id" for per-atom fix + * data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix + * property/atom vectors with count = 1, "d2_name" or "i2_name" for fix + * property/atom vectors with count > 1) + * \param type: 0 for ``int`` values, 1 for ``double`` values + * \param count: number of per-atom values (e.g., 1 for *type* or *charge*, 3 for *x* or *f*); + * use *count* = 3 with *image* if you want the image flags unpacked into + * (*x*,*y*,*z*) components. + * \param data: per-atom values packed into a one-dimensional array of length + * *natoms* \* *count*. + * + */ /* ---------------------------------------------------------------------- Contributing author: Thomas Swinburne (CNRS & CINaM, Marseille, France) gather the named atom-based entity for all atoms @@ -3633,6 +3718,41 @@ void lammps_gather_concat(void *handle, char *name, int type, int count, void *d END_CAPTURE } +/** Gather the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based entities + * from all processes for a subset of atoms. + * +\verbatim embed:rst + +This subroutine gathers data for the requested atom IDs and stores them in a one-dimensional array +allocated by the user. The data will be ordered by atom ID, but there is no requirement that the +IDs be consecutive. If you wish to return a similar array for *all* the atoms, use +:cpp:func:`lammps_gather` or :cpp:func:`lammps_gather_concat`. + +The *data* array will be in groups of *count* values, sorted by atom ID in the same order as the +array *ids* (e.g., if *name* is *x*, *count* = 3, and *ids* is {100, 57, 210}, then *data* might +look like {x[100][0], x[100][1], x[100][2], x[57][0], x[57][1], x[57][2], x[210][0], +:math:`\dots`}); *ids* must be provided by the user with length *ndata*, and *data* must be +pre-allocated by the caller to length (*count*\ :math:`{}\times{}`\ *ndata*). + +This function is not compatible with ``-DLAMMPS_BIGBIG``. + +\endverbatim + * + * \param handle: pointer to a previously created LAMMPS instance + * \param name desired quantity (e.g., "x" or "f" for atom properties, "f_id" for per-atom fix + * data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix + * property/atom vectors with *count* = 1, "d2_name" or "i2_name" for fix + * property/atom vectors with *count* > 1) + * \param type 0 for ``int`` values, 1 for ``double`` values + * \param count number of per-atom values (e.g., 1 for *type* or *charge*, 3 for *x* or *f*); + * use *count* = 3 with *image* if you want the image flags unpacked into + * (*x*,*y*,*z*) components. + * \param ndata: number of atoms for which to return data (can be all of them) + * \param ids: list of *ndata* atom IDs for which to return data + * \param data per-atom values packed into a one-dimensional array of length + * *ndata* \* *count*. + * + */ /* ---------------------------------------------------------------------- Contributing author: Thomas Swinburne (CNRS & CINaM, Marseille, France) gather the named atom-based entity for all atoms @@ -3884,6 +4004,37 @@ void lammps_gather_subset(void *handle, char *name, END_CAPTURE } +/** Scatter the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based + * entity in *data* to all processes. + * +\verbatim embed:rst + +This subroutine takes data stored in a one-dimensional array supplied by the user and scatters +them to all atoms on all processes. The data must be ordered by atom ID, with the requirement that +the IDs be consecutive. Use :cpp:func:`lammps_scatter_subset` to scatter data for some (or all) +atoms, unordered. + +The *data* array needs to be ordered in groups of *count* values, sorted by atom ID (e.g., if +*name* is *x* and *count* = 3, then *data* = {x[0][0], x[0][1], x[0][2], x[1][0], x[1][1], +x[1][2], x[2][0], :math:`\dots`}); *data* must be of length (*count* :math:`\times` *natoms*). + +This function is not compatible with ``-DLAMMPS_BIGBIG``. + +\endverbatim + * + * \param handle pointer to a previously created LAMMPS instance + * \param name desired quantity (e.g., "x" or "f" for atom properties, "f_id" for per-atom fix + * data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix + * property/atom vectors with *count* = 1, "d2_name" or "i2_name" for fix + * property/atom vectors with *count* > 1) + * \param type 0 for ``int`` values, 1 for ``double`` values + * \param count number of per-atom values (e.g., 1 for *type* or *charge*, + * 3 for *x* or *f*); use *count* = 3 with *image* if you have + * a single image flag packed into (*x*,*y*,*z*) components. + * \param data per-atom values packed in a one-dimensional array of length + * *natoms* \* *count*. + * + */ /* ---------------------------------------------------------------------- Contributing author: Thomas Swinburne (CNRS & CINaM, Marseille, France) scatter the named atom-based entity in data to all atoms @@ -4103,6 +4254,42 @@ void lammps_scatter(void *handle, char *name, int type, int count, void *data) END_CAPTURE } +/** Scatter the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based + * entities in *data* from a subset of atoms to all processes. + * +\verbatim embed:rst + +This subroutine takes data stored in a one-dimensional array supplied by the +user and scatters them to a subset of atoms on all processes. The array +*data* contains data associated with atom IDs, but there is no requirement that +the IDs be consecutive, as they are provided in a separate array. +Use :cpp:func:`lammps_scatter` to scatter data for all atoms, in order. + +The *data* array needs to be organized in groups of *count* values, with the +groups in the same order as the array *ids*. For example, if you want *data* +to be the array {x[1][0], x[1][1], x[1][2], x[100][0], x[100][1], x[100][2], +x[57][0], x[57][1], x[57][2]}, then *count* = 3, *ndata* = 3, and *ids* would +be {1, 100, 57}. + +This function is not compatible with ``-DLAMMPS_BIGBIG``. + +\endverbatim + * + * \param handle: pointer to a previously created LAMMPS instance + * \param name desired quantity (e.g., "x" or "f" for atom properties, "f_id" for per-atom fix + * data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix + * property/atom vectors with *count* = 1, "d2_name" or "i2_name" for fix + * property/atom vectors with *count* > 1) + * \param type: 0 for ``int`` values, 1 for ``double`` values + * \param count: number of per-atom values (e.g., 1 for *type* or *charge*, + * 3 for *x* or *f*); use *count* = 3 with "image" if you want + * single image flags unpacked into (*x*,*y*,*z*) + * \param ndata: number of atoms listed in *ids* and *data* arrays + * \param ids: list of *ndata* atom IDs to scatter data to + * \param data per-atom values packed in a 1-dimensional array of length + * *ndata* \* *count*. + * + */ /* ---------------------------------------------------------------------- Contributing author: Thomas Swinburne (CNRS & CINaM, Marseille, France) scatter the named atom-based entity in data to a subset of atoms From 2f321576c54c5beeffd782e1999222a94eb47f67 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Sat, 3 Dec 2022 16:03:29 -0600 Subject: [PATCH 12/19] Added documentation for scatter and gather; updated other docs --- doc/src/Fortran.rst | 134 ++++++++++++++++++++++++++------------------ fortran/lammps.f90 | 2 +- 2 files changed, 81 insertions(+), 55 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 80d3447437..a1f76a7267 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -39,21 +39,12 @@ can be found together with equivalent examples in C and C++ in the .. versionadded:: 9Oct2020 -.. admonition:: Work in Progress - :class: note - - This Fortran module is a work in progress, and only the documented - functionality is currently available. The final implementation should - cover the entire range of functionality available in the C and - Python library interfaces. - .. note:: - A contributed Fortran interface that more - closely resembles the C library interface is available in the - ``examples/COUPLE/fortran2`` folder. Please see the ``README`` file - in that folder for more information about it and how to contact its - author and maintainer. + A contributed Fortran interface that more closely resembles the C library + interface is available in the ``examples/COUPLE/fortran2`` folder. Please + see the ``README`` file in that folder for more information about it and how + to contact its author and maintainer. ---------- @@ -1280,11 +1271,11 @@ Procedures Bound to the :f:type:`lammps` Derived Type of atoms, see :f:func:`gather_atoms_subset`. The *data* array will be ordered in groups of *count* values, sorted by atom - ID (e.g., if *name* is *x* and *count* = 3, then *data* = *x*\ (1,1), + ID (e.g., if *name* is *x* and *count* = 3, then *data* = [*x*\ (1,1), *x*\ (2,1), *x*\ (3,1), *x*\ (1,2), *x*\ (2,2), *x*\ (3,2), *x*\ (1,3), - :math:`\dots`); *data* must be ``ALLOCATABLE`` and will be allocated to + :math:`\dots`]); *data* must be ``ALLOCATABLE`` and will be allocated to length (*count* :math:`\times` *natoms*), as queried by - :f:func:`extract_setting`. + :f:func:`get_natoms`. This function is not compatible with ``-DLAMMPS_BIGBIG``. @@ -1293,12 +1284,13 @@ Procedures Bound to the :f:type:`lammps` Derived Type (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use *count* = 3 with *image* if you want a single image flag unpacked into *x*/*y*/*z* components. - :p polymorphic data [dimension(:),allocatable]: array into which to store + :p data: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be reallocated to fit the length of the incoming data. It should have type ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if expecting floating-point data. + :ptype data: polymorphic,dimension(:),allocatable :to: :cpp:func:`lammps_gather_atoms` .. note:: @@ -1345,12 +1337,13 @@ Procedures Bound to the :f:type:`lammps` Derived Type (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use *count* = 3 with *image* if you want a single image flag unpacked into *x*/*y*/*z* components. - :p polymorphic data [dimension(:),allocatable]: array into which to store + :p data: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be reallocated to fit the length of the incoming data. It should have type ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if expecting floating-point data. + :ptype data: polymorphic,dimension(:),allocatable :to: :cpp:func:`lammps_gather_atoms_concat` -------- @@ -1386,12 +1379,13 @@ Procedures Bound to the :f:type:`lammps` Derived Type *x*/*y*/*z* components. :p integer(c_int) ids [dimension(:)]: atom IDs corresponding to the atoms to be gathered - :p polymorphic data [dimension(:),allocatable]: array into which to store + :p data: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be reallocated to fit the length of the incoming data. It should have type ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if expecting floating-point data. + :ptype data: polymorphic,dimension(:),allocatable :to: :cpp:func:`lammps_gather_atoms_subset` -------- @@ -1411,8 +1405,9 @@ Procedures Bound to the :f:type:`lammps` Derived Type The *data* array needs to be ordered in groups of *count* values, sorted by atom ID (e.g., if *name* is *x* and *count* = 3, then - *data* = [x(1,1), x(2,1), x(3,1), x(1,2), x(2,2), x(3,2), x(1,3), - :math:`\dots`]; *data* must be of length (*count* :math:`\times` *natoms*). + *data* = [*x*\ (1,1), *x*\ (2,1), *x*\ (3,1), *x*\ (1,2), *x*\ (2,2), + *x*\ (3,2), *x*\ (1,3), :math:`\dots`]); *data* must be of length *natoms* + or 3\*\ *natoms*. :p character(len=\*) name: quantity to be scattered (e.g., *x* or *charge*) :p data: per-atom values packed in a one-dimensional array @@ -1442,8 +1437,9 @@ Procedures Bound to the :f:type:`lammps` Derived Type The *data* array needs to be organized in groups of 1 or 3 values, depending on which quantity is being scattered, with the groups in the same order as the array *ids*. For example, if you want *data* to be the array - [x(1,1), x(2,1), x(3,1), x(1,100), x(2,100), x(3,100), x(1,57), x(2,57), - x(3,57)], then *ids* would be [1, 100, 57] and *name* would be *x*. + [*x*\ (1,1), *x*\ (2,1), *x*\ (3,1), *x*\ (1,100), *x*\ (2,100), + *x*\ (3,100), *x*\ (1,57), *x*\ (2,57), *x*\ (3,57)], then *ids* would be + [1, 100, 57] and *name* would be *x*. :p character(len=\*) name: quantity to be scattered (e.g., *x* or *charge*) :p integer(c_int) ids [dimension(:)]: atom IDs corresponding to the atoms @@ -1522,7 +1518,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type ID (e.g., if *name* is *x*, then *data* is [x(1,1), x(2,1), x(3,1), x(1,2), x(2,2), x(3,2), x(1,3), :math:`\dots`]); *data* must be ``ALLOCATABLE`` and will be allocated to length (*count*\ :math:`{}\times{}`\ *natoms*), as - queried by :f:func:`extract_setting`. + queried by :f:func:`get_natoms`. This function will return an error if fix or compute data are requested and the fix or compute ID given does not have per-atom data. See the note about @@ -1533,7 +1529,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p character(len=\*) name: desired quantity (e.g., "x" or "mask" for atom properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix property/atom vectors with *count* = 1, - "d2_name" or "i2_name" for fix propert/atom vectors with + "d2_name" or "i2_name" for fix property/atom vectors with *count*\ :math:`{}> 1`) :p integer(c_int) count: number of per-atom values (e.g., 1 for *type* or *charge*, 3 for *x* or *f*); use *count* = 3 with *image* if you want the @@ -1566,25 +1562,26 @@ Procedures Bound to the :f:type:`lammps` Derived Type then *data* might be something like [x(1,11), x(2,11), x(3,11), x(1,3), x(2,3), x(3,3), x(1,5), :math:`\dots`]); *data* must be ``ALLOCATABLE`` and will be allocated to length (*count* :math:`\times` *natoms*), as queried by - :f:func:`extract_setting`. + :f:func:`get_natoms`. This function is not compatible with ``-DLAMMPS_BIGBIG``. :p character(len=\*) name: desired quantity (e.g., "x" or "mask" for atom properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix property/atom vectors with *count* = 1, - "d2_name" or "i2_name" for fix propert/atom vectors with + "d2_name" or "i2_name" for fix property/atom vectors with *count*\ :math:`{}> 1`) :p integer(c_int) count: number of per-atom values you expect per atom (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use *count* = 3 with *image* if you want a single image flag unpacked into *x*/*y*/*z* components. - :p polymorphic data [dimension(:),allocatable]: array into which to store + :p data: array into which to store the data. Array *must* have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., ``DIMENSION(:)``). If this array is already allocated, it will be reallocated to fit the length of the incoming data. It should have type ``INTEGER(c_int)`` if expecting integer data and ``REAL(c_double)`` if expecting floating-point data. + :ptype data: polymorphic,dimension(:),allocatable :to: :cpp:func:`lammps_gather_concat` -------- @@ -1602,11 +1599,12 @@ Procedures Bound to the :f:type:`lammps` Derived Type The *data* array will be in groups of *count* values, sorted by atom ID in the same order as the array *ids* (e.g., if *name* is *x*, *count* = 3, and - *ids* is [100, 57, 210], then *data* might look like [x(1,100), x(2,100), - x(3,100), x(1,57), x(2,57), x(3,57), x(1,210), :math:`\dots`]); *ids* must - be provided by the user, and *data* must have the ``ALLOCATABLE`` attribute - and be of rank 1 (i.e., ``DIMENSION(:)``). If *data* is already allocated, - it will be reallocated to fit the length of the incoming data. + *ids* is [100, 57, 210], then *data* might look like [*x*\ (1,100), + *x*\ (2,100), *x*\ (3,100), *x*\ (1,57), *x*\ (2,57), *x*\ (3,57), + *x*\ (1,210), :math:`\dots`]); *ids* must be provided by the user, and + *data* must have the ``ALLOCATABLE`` attribute and be of rank 1 (i.e., + ``DIMENSION(:)``). If *data* is already allocated, it will be reallocated to + fit the length of the incoming data. This function is not compatible with ``-DLAMMPS_BIGBIG``. @@ -1617,24 +1615,52 @@ Procedures Bound to the :f:type:`lammps` Derived Type fix data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix property/atom vectors with *count* = 1, "d2_name" or "i2_name" for fix property/atom vectors with *count*\ :math:`{} > 1`) + :p integer(c_int) count: number of per-atom values you expect per atom + (e.g., 1 for *type*, *mask*, or *charge*; 3 for *x*, *v*, or *f*). Use + *count* = 3 with *image* if you want a single image flag unpacked into + *x*/*y*/*z* components. :p data: per-atom values packed into a one-dimensional array containing the - data to be scattered. This array must have either the same length as *ids* - (for *mask*, *type*, etc.) or three times its length (for *x*, *f*, etc.); - the array must be rank 1 and be of type ``INTEGER(c_int)`` (e.g., for - *mask* or *type*) or of type ``REAL(c_double)`` (e.g., for *charge*, *x*, - or *f*). - :ptype data: polymorphic,dimension(:) - :to: :cpp:func:`lammps_scatter_subset` + data to be scattered. This array must have the ``ALLOCATABLE`` attribute + and will be allocated either to the same length as *ids* + (for *mask*, *type*, etc.) or to three times its length (for *x*, *f*, + etc.); the array must be rank 1 and be of type ``INTEGER(c_int)`` (e.g., + for *mask* or *type*) or of type ``REAL(c_double)`` (e.g., for *charge*, + *x*, or *f*). + :ptype data: polymorphic,dimension(:),allocatable + :to: :cpp:func:`lammps_gather_subset` -------- .. f:subroutine:: scatter(name, data) This function calls :cpp:func:`lammps_scatter` to scatter the named - atom-based entities in *data* to all processes. + per-atom, per-atom fix, per-atom compute, or fix property/atom-based entity + in *data* to all processes. - *This function is not yet documented, as the underlying C routine has not - been thoroughly tested yet.* + .. versionadded:: TBD + + This subroutine takes data stored in a one-dimensional array supplied by the + user and scatters them to all atoms on all processes. The data must be + ordered by atom ID, with the requirement that the IDs be consecutive. Use + :f:subr:`scatter_subset` to scatter data for some (or all) atoms, unordered. + + The *data* array needs to be ordered in groups of *count* values, sorted by + atom ID (e.g., if *name* is *x* and *count* = 3, then *data* = [*x*\ (1,1), + *x*\ (2,1), *x*\ (3,1), *x*\ (1,2), *x*\ (2,2), *x*\ (3,2), *x*\ (1,3), + :math:`\dots`]); *data* must be of length (*count* :math:`\times` *natoms*). + + This function is not compatible with ``-DLAMMPS_BIGBIG``. + + :p character(len=*) name: desired quantity (e.g., "x" or "f" for atom + properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, + "d_name" or "i_name" for fix property/atom vectors with *count* = 1, + "d2_name" or "i2_name" for fix property/atom vectors with *count*\ + :math:`{} > 1`) + :p data: per-atom values packed in a one-dimensional array; *data* should be + of type ``INTEGER(c_int)`` or ``REAL(c_double)``, depending on the type of + data being scattered, and be of rank 1 (i.e., ``DIMENSION(:)``). + :ptype data: polymorphic,dimension(:) + :to: :cpp:func:`lammps_scatter` -------- @@ -1662,7 +1688,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p character(len=\*) name: desired quantity (e.g., "x" or "mask" for atom properties, "f_id" for per-atom fix data, "c_id" for per-atom compute data, "d_name" or "i_name" for fix property/atom vectors with *count* = 1, - "d2_name" or "i2_name" for fix propert/atom vectors with + "d2_name" or "i2_name" for fix property/atom vectors with *count*\ :math:`{}> 1`) :p integer(c_int) ids: list of atom IDs to scatter data for :p polymorphic data [dimension(:)]: per-atom values packed in a @@ -1681,7 +1707,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p integer(c_int) type [dimension(N)]: vector of :math:`N` atom types (required/see note below) - :p real(c_double) x [dimension(3N)]: vector of :math:`3N` x/y/z positions + :p real(c_double) x [dimension(3N)]: vector of :math:`3N\ x/y/z` positions of the new atoms, arranged as :math:`[x_1,y_1,z_1,x_2,y_2,\dotsc]` (required/see note below) :o integer(kind=\*) id [dimension(N),optional]: vector of :math:`N` atom @@ -2384,10 +2410,11 @@ Procedures Bound to the :f:type:`lammps` Derived Type MODULE stuff USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double, c_int64_t + USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : error_unit IMPLICIT NONE TYPE shield - REAL(c_double), DIMENSION(:), ALLOCATABLE :: k + REAL(c_double), DIMENSION(:,:), ALLOCATABLE :: k ! assume k gets allocated to dimension(3,nlocal) at some point ! and assigned values END TYPE shield @@ -2632,7 +2659,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type This function needs to be called **after** a call to :f:subr:`fix_external_set_vector_length` and **before** a run or minimize command. When running in parallel, it must be called from **all** MPI - processes with the **same**\ *idx* and *val* parameters. The variable + processes with the **same** *idx* and *val* parameters. The variable *val* is assumed to be extensive. .. note:: @@ -2686,7 +2713,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. versionadded:: 3Nov2022 This function can be used from signal handlers or multi-threaded - applications to cleanly terminate an ongoing run. + applications to terminate an ongoing run cleanly. :to: :cpp:func:`lammps_force_timeout` @@ -2725,12 +2752,11 @@ Procedures Bound to the :f:type:`lammps` Derived Type in the event of an error inside of LAMMPS that resulted in a :ref:`C++ exception `. A suitable buffer for a string has to be provided. If the internally-stored error message is longer than the - string and the string does not have ``ALLOCATABLE`` length, it will be - truncated accordingly. The optional argument *status* indicates the - kind of error: a "1" indicates an error that occurred on all MPI ranks and - is often recoverable, while a "2" indicates an abort that would happen only - in a single MPI rank and thus may not be recoverable, as other MPI ranks may - be waiting on the failing MPI rank(s) to send messages. + string, it will be truncated accordingly. The optional argument *status* + indicates the kind of error: a "1" indicates an error that occurred on all + MPI ranks and is often recoverable, while a "2" indicates an abort that + would happen only in a single MPI rank and thus may not be recoverable, as + other MPI ranks may be waiting on the failing MPI rank(s) to send messages. .. note:: diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index 4cf90e22a2..f8803e29e9 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -2098,7 +2098,7 @@ CONTAINS END SUBROUTINE lmp_scatter_int ! equivalent function to lammps_scatter (for doubles) - SUBROUTINE lmp_scatter_atoms_double(self, name, data) + SUBROUTINE lmp_scatter_double(self, name, data) CLASS(lammps), INTENT(IN) :: self CHARACTER(LEN=*), INTENT(IN) :: name REAL(c_double), DIMENSION(:), TARGET :: data From c0345845e86b73bc32d4174df37392d912e86daa Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Sat, 3 Dec 2022 20:38:42 -0600 Subject: [PATCH 13/19] unit test for gather and scatter; char* to const char* in library.* --- src/library.cpp | 32 +++--- src/library.h | 28 ++--- .../fortran/test_fortran_gather_scatter.f90 | 92 ++++++++++++++- unittest/fortran/wrap_gather_scatter.cpp | 108 ++++++++++++++++++ 4 files changed, 231 insertions(+), 29 deletions(-) diff --git a/src/library.cpp b/src/library.cpp index 9dc9a701f3..2c16cc6615 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -2303,7 +2303,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allreduce to sum vector into data across all procs ------------------------------------------------------------------------- */ -void lammps_gather_atoms(void *handle, char *name, int type, int count, void *data) +void lammps_gather_atoms(void *handle, const char *name, int type, int count, + void *data) { auto lmp = (LAMMPS *) handle; @@ -2460,7 +2461,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allgather Nlocal atoms from each proc into data ------------------------------------------------------------------------- */ -void lammps_gather_atoms_concat(void *handle, char *name, int type, int count, void *data) +void lammps_gather_atoms_concat(void *handle, const char *name, int type, + int count, void *data) { auto lmp = (LAMMPS *) handle; @@ -2627,8 +2629,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allreduce to sum vector into data across all procs ------------------------------------------------------------------------- */ -void lammps_gather_atoms_subset(void *handle, char *name, int type, int count, - int ndata, int *ids, void *data) +void lammps_gather_atoms_subset(void *handle, const char *name, int type, + int count, int ndata, int *ids, void *data) { auto lmp = (LAMMPS *) handle; @@ -2786,7 +2788,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. loop over Natoms, if I own atom ID, set its values from data ------------------------------------------------------------------------- */ -void lammps_scatter_atoms(void *handle, char *name, int type, int count, void *data) +void lammps_scatter_atoms(void *handle, const char *name, int type, int count, + void *data) { auto lmp = (LAMMPS *) handle; @@ -2938,8 +2941,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. loop over Ndata, if I own atom ID, set its values from data ------------------------------------------------------------------------- */ -void lammps_scatter_atoms_subset(void *handle, char *name, int type, int count, - int ndata, int *ids, void *data) +void lammps_scatter_atoms_subset(void *handle, const char *name, int type, + int count, int ndata, int *ids, void *data) { auto lmp = (LAMMPS *) handle; @@ -3219,7 +3222,7 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allreduce to sum vector into data across all procs ------------------------------------------------------------------------- */ -void lammps_gather(void *handle, char *name, int type, int count, void *data) +void lammps_gather(void *handle, const char *name, int type, int count, void *data) { auto lmp = (LAMMPS *) handle; @@ -3492,7 +3495,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allreduce to sum vector into data across all procs ------------------------------------------------------------------------- */ -void lammps_gather_concat(void *handle, char *name, int type, int count, void *data) +void lammps_gather_concat(void *handle, const char *name, int type, int count, + void *data) { auto lmp = (LAMMPS *) handle; @@ -3779,9 +3783,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allreduce to sum vector into data across all procs ------------------------------------------------------------------------- */ -void lammps_gather_subset(void *handle, char *name, - int type, int count, - int ndata, int *ids, void *data) +void lammps_gather_subset(void *handle, const char *name, int type, int count, + int ndata, int *ids, void *data) { auto lmp = (LAMMPS *) handle; @@ -4059,7 +4062,8 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. Allreduce to sum vector into data across all procs ------------------------------------------------------------------------- */ -void lammps_scatter(void *handle, char *name, int type, int count, void *data) +void lammps_scatter(void *handle, const char *name, int type, int count, + void *data) { auto lmp = (LAMMPS *) handle; @@ -4312,7 +4316,7 @@ This function is not compatible with ``-DLAMMPS_BIGBIG``. loop over Ndata, if I own atom ID, set its values from data ------------------------------------------------------------------------- */ -void lammps_scatter_subset(void *handle, char *name,int type, int count, +void lammps_scatter_subset(void *handle, const char *name,int type, int count, int ndata, int *ids, void *data) { auto lmp = (LAMMPS *) handle; diff --git a/src/library.h b/src/library.h index f0f8448b79..9cd2f34d73 100644 --- a/src/library.h +++ b/src/library.h @@ -181,23 +181,23 @@ int lammps_set_variable(void *, char *, char *); * Library functions for scatter/gather operations of data * ---------------------------------------------------------------------- */ -void lammps_gather_atoms(void *handle, char *name, int type, int count, void *data); -void lammps_gather_atoms_concat(void *handle, char *name, int type, int count, void *data); -void lammps_gather_atoms_subset(void *handle, char *name, int type, int count, int ndata, int *ids, - void *data); -void lammps_scatter_atoms(void *handle, char *name, int type, int count, void *data); -void lammps_scatter_atoms_subset(void *handle, char *name, int type, int count, int ndata, int *ids, - void *data); +void lammps_gather_atoms(void *handle, const char *name, int type, int count, void *data); +void lammps_gather_atoms_concat(void *handle, const char *name, int type, int count, void *data); +void lammps_gather_atoms_subset(void *handle, const char *name, int type, int count, int ndata, + int *ids, void *data); +void lammps_scatter_atoms(void *handle, const char *name, int type, int count, void *data); +void lammps_scatter_atoms_subset(void *handle, const char *name, int type, int count, int ndata, + int *ids, void *data); void lammps_gather_bonds(void *handle, void *data); -void lammps_gather(void *handle, char *name, int type, int count, void *data); -void lammps_gather_concat(void *handle, char *name, int type, int count, void *data); -void lammps_gather_subset(void *handle, char *name, int type, int count, int ndata, int *ids, - void *data); -void lammps_scatter(void *handle, char *name, int type, int count, void *data); -void lammps_scatter_subset(void *handle, char *name, int type, int count, int ndata, int *ids, - void *data); +void lammps_gather(void *handle, const char *name, int type, int count, void *data); +void lammps_gather_concat(void *handle, const char *name, int type, int count, void *data); +void lammps_gather_subset(void *handle, const char *name, int type, int count, int ndata, + int *ids, void *data); +void lammps_scatter(void *handle, const char *name, int type, int count, void *data); +void lammps_scatter_subset(void *handle, const char *name, int type, int count, int ndata, + int *ids, void *data); #if !defined(LAMMPS_BIGBIG) int lammps_create_atoms(void *handle, int n, const int *id, const int *type, const double *x, diff --git a/unittest/fortran/test_fortran_gather_scatter.f90 b/unittest/fortran/test_fortran_gather_scatter.f90 index 2115a10394..6df1421bac 100644 --- a/unittest/fortran/test_fortran_gather_scatter.f90 +++ b/unittest/fortran/test_fortran_gather_scatter.f90 @@ -24,13 +24,16 @@ END SUBROUTINE f_lammps_close SUBROUTINE f_lammps_setup_gather_scatter() BIND(C) USE LIBLAMMPS - USE keepstuff, ONLY : lmp, big_input, cont_input, more_input + USE keepstuff, ONLY : lmp, big_input, cont_input, more_input, pair_input IMPLICIT NONE CALL lmp%command('atom_modify map array') CALL lmp%commands_list(big_input) CALL lmp%commands_list(cont_input) CALL lmp%commands_list(more_input) + CALL lmp%commands_list(pair_input) + CALL lmp%command('mass 1 1.0') + CALL lmp%command("compute pe all pe/atom") END SUBROUTINE f_lammps_setup_gather_scatter FUNCTION f_lammps_gather_atoms_mask(i) BIND(C) @@ -262,3 +265,90 @@ FUNCTION f_lammps_test_gather_bonds_big() BIND(C) RESULT(success) success = 0_c_int END IF END FUNCTION f_lammps_test_gather_bonds_big + +FUNCTION f_lammps_gather_pe_atom(i) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i + REAL(c_double) :: f_lammps_gather_pe_atom + REAL(c_double), DIMENSION(:), ALLOCATABLE :: pe_atom + + CALL lmp%gather('c_pe', 1_c_int, pe_atom) + f_lammps_gather_pe_atom = pe_atom(i) +END FUNCTION f_lammps_gather_pe_atom + +FUNCTION f_lammps_gather_pe_atom_concat(i) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i + REAL(c_double) :: f_lammps_gather_pe_atom_concat + REAL(c_double), DIMENSION(:), ALLOCATABLE :: pe_atom + INTEGER(c_int), DIMENSION(:), ALLOCATABLE :: tag + INTEGER :: j + + CALL lmp%gather_concat('id', 1_c_int, tag) + CALL lmp%gather_concat('c_pe', 1_c_int, pe_atom) + DO j = 1, SIZE(tag) + IF (tag(j) == i) THEN + f_lammps_gather_pe_atom_concat = pe_atom(j) + EXIT + END IF + END DO + f_lammps_gather_pe_atom_concat = pe_atom(i) +END FUNCTION f_lammps_gather_pe_atom_concat + +SUBROUTINE f_lammps_gather_pe_atom_subset(ids, pe) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN) :: ids(2) + REAL(c_double), INTENT(OUT) :: pe(2) + REAL(c_double), DIMENSION(:), ALLOCATABLE :: pe_atom + INTEGER(c_int) :: natoms + + natoms = NINT(lmp%get_natoms(), c_int) + CALL lmp%gather_subset('c_pe', 1, ids, pe_atom) + pe(1:natoms) = pe_atom +END SUBROUTINE f_lammps_gather_pe_atom_subset + +SUBROUTINE f_lammps_scatter_compute() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + REAL(c_double), DIMENSION(:), ALLOCATABLE :: pe_atom + REAL(c_double) :: swap + + CALL lmp%gather('c_pe', 1_c_int, pe_atom) + + ! swap the computed energy of atoms 1 and 3 + swap = pe_atom(1) + pe_atom(1) = pe_atom(3) + pe_atom(3) = swap + + CALL lmp%scatter('c_pe', pe_atom) ! push the swap back to LAMMPS +END SUBROUTINE f_lammps_scatter_compute + +SUBROUTINE f_lammps_scatter_subset_compute() BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int, c_double + USE LIBLAMMPS + USE keepstuff, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), PARAMETER :: ids(2) = [3,1] + REAL(c_double), DIMENSION(:), ALLOCATABLE :: pe_atom + REAL(c_double) :: swap + + CALL lmp%gather_subset('c_pe', 1_c_int, ids, pe_atom) + + ! swap the computed energy of atoms 1 and 3 + swap = pe_atom(1) + pe_atom(1) = pe_atom(2) + pe_atom(2) = swap + + CALL lmp%scatter_subset('c_pe', ids, pe_atom) ! push the swap back to LAMMPS +END SUBROUTINE f_lammps_scatter_subset_compute diff --git a/unittest/fortran/wrap_gather_scatter.cpp b/unittest/fortran/wrap_gather_scatter.cpp index 9eb6082202..5e836c1759 100644 --- a/unittest/fortran/wrap_gather_scatter.cpp +++ b/unittest/fortran/wrap_gather_scatter.cpp @@ -3,6 +3,7 @@ #include "lammps.h" #include "library.h" +#include "atom.h" #include #include #include @@ -26,6 +27,11 @@ void f_lammps_scatter_atoms_positions(); void f_lammps_setup_gather_bonds(); int f_lammps_test_gather_bonds_small(); int f_lammps_test_gather_bonds_big(); +double f_lammps_gather_pe_atom(int); +double f_lammps_gather_pe_atom_concat(int); +void f_lammps_gather_pe_atom_subset(int*, double*); +void f_lammps_scatter_compute(); +void f_lammps_scatter_subset_compute(); } using namespace LAMMPS_NS; @@ -216,3 +222,105 @@ TEST_F(LAMMPS_gather_scatter, gather_bonds) #endif }; + +TEST_F(LAMMPS_gather_scatter, gather_compute) +{ +#ifdef LAMMPS_BIGBIG + GTEST_SKIP() +#endif + f_lammps_setup_gather_scatter(); + lammps_command(lmp, "run 0"); + int natoms = lmp->atom->natoms; + int *tag = lmp->atom->tag; + double *pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, + LMP_TYPE_VECTOR); + for (int i = 0; i < natoms; i++) + EXPECT_DOUBLE_EQ(f_lammps_gather_pe_atom(tag[i]), pe[i]); +}; + +TEST_F(LAMMPS_gather_scatter, gather_compute_concat) +{ +#ifdef LAMMPS_BIGBIG + GTEST_SKIP() +#endif + f_lammps_setup_gather_scatter(); + lammps_command(lmp, "run 0"); + int natoms = lmp->atom->natoms; + int *tag = lmp->atom->tag; + double *pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, + LMP_TYPE_VECTOR); + for (int i = 0; i < natoms; i++) + EXPECT_DOUBLE_EQ(f_lammps_gather_pe_atom(tag[i]), pe[i]); +}; + +TEST_F(LAMMPS_gather_scatter, gather_compute_subset) +{ + f_lammps_setup_gather_scatter(); + lammps_command(lmp, "run 0"); + int ids[2] = {3, 1}; + int *tag = lmp->atom->tag; + double pe[2] = {0.0, 0.0}; + int nlocal = lammps_extract_setting(lmp, "nlocal"); + double *pa_pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, + LMP_TYPE_VECTOR); + + for (int i = 0; i < nlocal; i++) { + if(tag[i] == ids[0]) pe[0] = pa_pe[i]; + if(tag[i] == ids[1]) pe[1] = pa_pe[i]; + } + + double ftn_pe[2]; + f_lammps_gather_pe_atom_subset(ids, ftn_pe); + EXPECT_DOUBLE_EQ(ftn_pe[0], pe[0]); + EXPECT_DOUBLE_EQ(ftn_pe[1], pe[1]); +}; + +TEST_F(LAMMPS_gather_scatter, scatter_compute) +{ +#ifdef LAMMPS_BIGBIG + GTEST_SKIP(); +#endif + f_lammps_setup_gather_scatter(); + int natoms = lmp->atom->natoms; + double *pe = new double[natoms]; + lammps_command(lmp, "run 0"); + lammps_gather(lmp, "c_pe", 1, 1, pe); + double *old_pe = new double[natoms]; + for (int i = 0; i < natoms; i++) + old_pe[i] = pe[i]; + EXPECT_DOUBLE_EQ(pe[0], old_pe[0]); + EXPECT_DOUBLE_EQ(pe[1], old_pe[1]); + EXPECT_DOUBLE_EQ(pe[2], old_pe[2]); + f_lammps_scatter_compute(); + lammps_gather(lmp, "c_pe", 1, 1, pe); + EXPECT_DOUBLE_EQ(pe[0], old_pe[2]); + EXPECT_DOUBLE_EQ(pe[1], old_pe[1]); + EXPECT_DOUBLE_EQ(pe[2], old_pe[0]); + delete[] old_pe; + delete[] pe; +}; + +TEST_F(LAMMPS_gather_scatter, scatter_subset_compute) +{ +#ifdef LAMMPS_BIGBIG + GTEST_SKIP(); +#endif + f_lammps_setup_gather_scatter(); + int natoms = lmp->atom->natoms; + double *pe = new double[natoms]; + lammps_command(lmp, "run 0"); + lammps_gather(lmp, "c_pe", 1, 1, pe); + double *old_pe = new double[natoms]; + for (int i = 0; i < natoms; i++) + old_pe[i] = pe[i]; + EXPECT_DOUBLE_EQ(pe[0], old_pe[0]); + EXPECT_DOUBLE_EQ(pe[1], old_pe[1]); + EXPECT_DOUBLE_EQ(pe[2], old_pe[2]); + f_lammps_scatter_subset_compute(); + lammps_gather(lmp, "c_pe", 1, 1, pe); + EXPECT_DOUBLE_EQ(pe[0], old_pe[2]); + EXPECT_DOUBLE_EQ(pe[1], old_pe[1]); + EXPECT_DOUBLE_EQ(pe[2], old_pe[0]); + delete[] old_pe; + delete[] pe; +}; From 411f9b450f4fd052de895d4d3c6cbef98ca3c381 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Sat, 3 Dec 2022 21:39:14 -0600 Subject: [PATCH 14/19] documented two overlooked functions; added NULL check to neighlist_element_neighbors --- doc/src/Fortran.rst | 39 +++++++++++++++++++++++++++++++++++++++ fortran/lammps.f90 | 7 ++++++- 2 files changed, 45 insertions(+), 1 deletion(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index a1f76a7267..27b2dc4a7f 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -277,6 +277,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype extract_fix: function :f extract_variable: :f:func:`extract_variable` :ftype extract_variable: function + :f set_variable: :f:subr:`set_variable` + :ftype set_variable: subroutine :f gather_atoms: :f:subr:`gather_atoms` :ftype gather_atoms: subroutine :f gather_atoms_concat: :f:subr:`gather_atoms_concat` @@ -309,6 +311,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype find_compute_neighlist: function :f neighlist_num_elements: :f:func:`neighlist_num_elements` :ftype neighlist_num_elements: function + :f neighlist_element_neighbors: :f:subr:`neighlist_element_neighbors` + :ftype neighlist_element_neighbors: subroutine :f version: :f:func:`version` :ftype version: function :f get_os_info: :f:subr:`get_os_info` @@ -1257,6 +1261,22 @@ Procedures Bound to the :f:type:`lammps` Derived Type -------- +.. f:subroutine:: set_variable(name, str) + + Set the value of a string-style variable. + + .. versionadded:: 3Nov2022 + + This function assigns a new value from the string *str* to the string-style + variable *name*\ . If *name* does not exist or is not a string-style + variable, an error is generated. + + :p character(len=*) name: name of the variable + :p character(len=*) str: new value to assign to the variable + :to: :cpp:func:`lammps_set_variable` + +-------- + .. f:subroutine:: gather_atoms(name, count, data) This function calls :cpp:func:`lammps_gather_atoms` to gather the named @@ -1837,6 +1857,25 @@ Procedures Bound to the :f:type:`lammps` Derived Type :r inum: number of entries in neighbor list, or :math:`-1` if *idx* is not a valid index. :rtype inum: integer(c_int) + :to: :cpp:func:`lammps_neighlist_num_elements` + +-------- + +.. f:subroutine:: neighlist_element_neighbors(idx, element, iatom, neighbors) + + Return atom local index, number of neighbors, and array of neighbor local + atom indices of a neighbor list entry. + + .. versionadded:: 3Nov2022 + + :p integer(c_int) idx: index of this neighbor list in the list of all + neighbor lists + :p integer(c_int) element: index of this neighbor list entry + :p integer(c_int) iatom: local atom index (i.e., in the range + [1,nlocal+nghost]; -1 if invalid or element value + :p integer(c_int) neighbors [dimension(:),pointer]: pointer to an array of + neighboring atom local indices + :to: :cpp:func:`lammps_neighlist_element_neighbors` -------- diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index f8803e29e9..e25b5ec206 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -2406,7 +2406,12 @@ CONTAINS CALL lammps_neighlist_element_neighbors(self%handle, idx, element, iatom, & numneigh, Cneighbors) - CALL C_F_POINTER(Cneighbors, neighbors, [numneigh]) + IF (C_ASSOCIATED(Cneighbors)) THEN + CALL C_F_POINTER(Cneighbors, neighbors, [numneigh]) + ELSE + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Pointer returned from lammps_neighlist_element_neighbors is NULL') + END IF END SUBROUTINE lmp_neighlist_element_neighbors ! equivalent function to lammps_version From f381a78c4676a146d63bd89929d06a56419e167b Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Sat, 3 Dec 2022 21:52:30 -0600 Subject: [PATCH 15/19] Added missing "call to" in Fortran docs --- doc/src/Fortran.rst | 3 +++ 1 file changed, 3 insertions(+) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 27b2dc4a7f..0f292af108 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -2432,6 +2432,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :ptype callback: external :p class(*) caller [optional]: object you wish to pass to the callback procedure (must be a scalar; see note) + :to: :cpp:func:`lammps_set_fix_external_callback` .. note:: @@ -2680,6 +2681,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p character(len=*) id: fix ID of fix external instance :p integer(c_int) length: length of the global vector to be stored with the fix + :to: :cpp:func:`lammps_fix_external_set_vector_length` -------- @@ -2715,6 +2717,7 @@ Procedures Bound to the :f:type:`lammps` Derived Type :p character(len=*) id: ID of fix external instance :p integer(c_int) idx: 1-based index in global vector :p integer(c_int) val: value to be stored in global vector at index *idx* + :to: :cpp:func:`lammps_fix_external_set_vector` -------- From 162f2f938446350f8a8322b23de1a51e0025220d Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Sun, 4 Dec 2022 00:27:48 -0600 Subject: [PATCH 16/19] whitespace; versionadded tags --- doc/src/Fortran.rst | 6 +- unittest/fortran/wrap_gather_scatter.cpp | 87 +++++++++++++----------- 2 files changed, 52 insertions(+), 41 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 0f292af108..4c450a0c08 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -1567,6 +1567,8 @@ Procedures Bound to the :f:type:`lammps` Derived Type Gather the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based entities from all processes, unordered. + .. versionadded:: TBD + This subroutine gathers data for all atoms and stores them in a one-dimensional allocatable array. The data will be a concatenation of chunks from each processor's owned atoms, in whatever order @@ -1610,7 +1612,9 @@ Procedures Bound to the :f:type:`lammps` Derived Type Gather the named per-atom, per-atom fix, per-atom compute, or fix property/atom-based entities from all processes for a subset of atoms. - + + .. versionadded:: TBD + This subroutine gathers data for the requested atom IDs and stores them in a one-dimensional allocatable array. The data will be ordered by atom ID, but there is no requirement that the IDs be consecutive. If you wish to return a diff --git a/unittest/fortran/wrap_gather_scatter.cpp b/unittest/fortran/wrap_gather_scatter.cpp index 5e836c1759..3aea3b8d44 100644 --- a/unittest/fortran/wrap_gather_scatter.cpp +++ b/unittest/fortran/wrap_gather_scatter.cpp @@ -214,72 +214,77 @@ TEST_F(LAMMPS_gather_scatter, scatter_atoms_subset_mask) TEST_F(LAMMPS_gather_scatter, gather_bonds) { - f_lammps_setup_gather_bonds(); + f_lammps_setup_gather_bonds(); #ifdef LAMMPS_BIGBIG - EXPECT_EQ(f_lammps_test_gather_bonds_big(), 1); + EXPECT_EQ(f_lammps_test_gather_bonds_big(), 1); #else - EXPECT_EQ(f_lammps_test_gather_bonds_small(), 1); + EXPECT_EQ(f_lammps_test_gather_bonds_small(), 1); #endif - }; TEST_F(LAMMPS_gather_scatter, gather_compute) { #ifdef LAMMPS_BIGBIG - GTEST_SKIP() + GTEST_SKIP(); +#else + f_lammps_setup_gather_scatter(); + lammps_command(lmp, "run 0"); + int natoms = lmp->atom->natoms; + int *tag = lmp->atom->tag; + double *pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, + LMP_TYPE_VECTOR); + for (int i = 0; i < natoms; i++) + EXPECT_DOUBLE_EQ(f_lammps_gather_pe_atom(tag[i]), pe[i]); #endif - f_lammps_setup_gather_scatter(); - lammps_command(lmp, "run 0"); - int natoms = lmp->atom->natoms; - int *tag = lmp->atom->tag; - double *pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, - LMP_TYPE_VECTOR); - for (int i = 0; i < natoms; i++) - EXPECT_DOUBLE_EQ(f_lammps_gather_pe_atom(tag[i]), pe[i]); }; TEST_F(LAMMPS_gather_scatter, gather_compute_concat) { #ifdef LAMMPS_BIGBIG - GTEST_SKIP() + GTEST_SKIP(); +#else + f_lammps_setup_gather_scatter(); + lammps_command(lmp, "run 0"); + int natoms = lmp->atom->natoms; + int *tag = lmp->atom->tag; + double *pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, + LMP_TYPE_VECTOR); + for (int i = 0; i < natoms; i++) + EXPECT_DOUBLE_EQ(f_lammps_gather_pe_atom(tag[i]), pe[i]); #endif - f_lammps_setup_gather_scatter(); - lammps_command(lmp, "run 0"); - int natoms = lmp->atom->natoms; - int *tag = lmp->atom->tag; - double *pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, - LMP_TYPE_VECTOR); - for (int i = 0; i < natoms; i++) - EXPECT_DOUBLE_EQ(f_lammps_gather_pe_atom(tag[i]), pe[i]); }; TEST_F(LAMMPS_gather_scatter, gather_compute_subset) { - f_lammps_setup_gather_scatter(); - lammps_command(lmp, "run 0"); - int ids[2] = {3, 1}; - int *tag = lmp->atom->tag; - double pe[2] = {0.0, 0.0}; - int nlocal = lammps_extract_setting(lmp, "nlocal"); - double *pa_pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, - LMP_TYPE_VECTOR); +#ifdef LAMMPS_BIGBIG + GTEST_SKIP(); +#else + f_lammps_setup_gather_scatter(); + lammps_command(lmp, "run 0"); + int ids[2] = {3, 1}; + int *tag = lmp->atom->tag; + double pe[2] = {0.0, 0.0}; + int nlocal = lammps_extract_setting(lmp, "nlocal"); + double *pa_pe = (double*) lammps_extract_compute(lmp, "pe", LMP_STYLE_ATOM, + LMP_TYPE_VECTOR); - for (int i = 0; i < nlocal; i++) { - if(tag[i] == ids[0]) pe[0] = pa_pe[i]; - if(tag[i] == ids[1]) pe[1] = pa_pe[i]; - } + for (int i = 0; i < nlocal; i++) { + if(tag[i] == ids[0]) pe[0] = pa_pe[i]; + if(tag[i] == ids[1]) pe[1] = pa_pe[i]; + } - double ftn_pe[2]; - f_lammps_gather_pe_atom_subset(ids, ftn_pe); - EXPECT_DOUBLE_EQ(ftn_pe[0], pe[0]); - EXPECT_DOUBLE_EQ(ftn_pe[1], pe[1]); + double ftn_pe[2]; + f_lammps_gather_pe_atom_subset(ids, ftn_pe); + EXPECT_DOUBLE_EQ(ftn_pe[0], pe[0]); + EXPECT_DOUBLE_EQ(ftn_pe[1], pe[1]); +#endif }; TEST_F(LAMMPS_gather_scatter, scatter_compute) { #ifdef LAMMPS_BIGBIG GTEST_SKIP(); -#endif +#else f_lammps_setup_gather_scatter(); int natoms = lmp->atom->natoms; double *pe = new double[natoms]; @@ -298,13 +303,14 @@ TEST_F(LAMMPS_gather_scatter, scatter_compute) EXPECT_DOUBLE_EQ(pe[2], old_pe[0]); delete[] old_pe; delete[] pe; +#endif }; TEST_F(LAMMPS_gather_scatter, scatter_subset_compute) { #ifdef LAMMPS_BIGBIG GTEST_SKIP(); -#endif +#else f_lammps_setup_gather_scatter(); int natoms = lmp->atom->natoms; double *pe = new double[natoms]; @@ -323,4 +329,5 @@ TEST_F(LAMMPS_gather_scatter, scatter_subset_compute) EXPECT_DOUBLE_EQ(pe[2], old_pe[0]); delete[] old_pe; delete[] pe; +#endif }; From 0984b11cb4f074ad1009c488a01ce60c13cfcad0 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 11 Dec 2022 17:50:47 -0500 Subject: [PATCH 17/19] skip gather_bonds test when atom style full is no available --- unittest/fortran/wrap_gather_scatter.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/unittest/fortran/wrap_gather_scatter.cpp b/unittest/fortran/wrap_gather_scatter.cpp index 3aea3b8d44..5eecf763a6 100644 --- a/unittest/fortran/wrap_gather_scatter.cpp +++ b/unittest/fortran/wrap_gather_scatter.cpp @@ -214,6 +214,7 @@ TEST_F(LAMMPS_gather_scatter, scatter_atoms_subset_mask) TEST_F(LAMMPS_gather_scatter, gather_bonds) { + if (!lammps_has_style(lmp, "atom", "full")) GTEST_SKIP(); f_lammps_setup_gather_bonds(); #ifdef LAMMPS_BIGBIG EXPECT_EQ(f_lammps_test_gather_bonds_big(), 1); From c0a39dc7b81adba44f6530908839aec12d5290fd Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 11 Dec 2022 17:50:03 -0500 Subject: [PATCH 18/19] add c wrapper to allow testing fortran interface w/o fortran MPI libs --- unittest/fortran/CMakeLists.txt | 52 ++++++++---------------- unittest/fortran/mpi_stubs.f90 | 20 --------- unittest/fortran/test_fortran_create.f90 | 29 ++++++++++--- unittest/fortran/wrap_create.cpp | 8 ++++ 4 files changed, 48 insertions(+), 61 deletions(-) delete mode 100644 unittest/fortran/mpi_stubs.f90 diff --git a/unittest/fortran/CMakeLists.txt b/unittest/fortran/CMakeLists.txt index d12b274765..08bd1b4356 100644 --- a/unittest/fortran/CMakeLists.txt +++ b/unittest/fortran/CMakeLists.txt @@ -21,20 +21,6 @@ if(CMAKE_Fortran_COMPILER) endif() get_filename_component(LAMMPS_FORTRAN_MODULE ${LAMMPS_SOURCE_DIR}/../fortran/lammps.f90 ABSOLUTE) - if(BUILD_MPI) - find_package(MPI REQUIRED) - if((NOT MPI_Fortran_FOUND) OR (NOT MPI_Fortran_HAVE_F77_HEADER)) - message(STATUS "Skipping Tests for the LAMMPS Fortran Module: no MPI support for Fortran") - return() - endif() - else() - add_library(fmpi_stubs STATIC mpi_stubs.f90) - get_filename_component(_tmp_fc ${CMAKE_Fortran_COMPILER} NAME) - if (_tmp_fc STREQUAL "flang") - target_link_libraries(fmpi_stubs PUBLIC gfortran) - endif() - add_library(MPI::MPI_Fortran ALIAS fmpi_stubs) - endif() add_library(flammps STATIC ${LAMMPS_FORTRAN_MODULE} keepstuff.f90) get_filename_component(_tmp_fc ${CMAKE_Fortran_COMPILER} NAME) @@ -42,66 +28,62 @@ if(CMAKE_Fortran_COMPILER) target_link_libraries(flammps PUBLIC gfortran) endif() - if(MPI_Fortran_HAVE_F90_MODULE) - add_executable(test_fortran_create wrap_create.cpp test_fortran_create.f90) - target_link_libraries(test_fortran_create PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) - target_include_directories(test_fortran_create PRIVATE "${LAMMPS_SOURCE_DIR}/../fortran") - add_test(NAME FortranOpen COMMAND test_fortran_create) - else() - message(STATUS "Skipping FortranOpen test since no working F90 MPI module was found") - endif() + add_executable(test_fortran_create wrap_create.cpp test_fortran_create.f90) + target_link_libraries(test_fortran_create PRIVATE flammps lammps GTest::GTestMain) + target_include_directories(test_fortran_create PRIVATE "${LAMMPS_SOURCE_DIR}/../fortran") + add_test(NAME FortranOpen COMMAND test_fortran_create) add_executable(test_fortran_commands wrap_commands.cpp test_fortran_commands.f90) - target_link_libraries(test_fortran_commands PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_commands PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranCommands COMMAND test_fortran_commands) add_executable(test_fortran_get_thermo wrap_get_thermo.cpp test_fortran_get_thermo.f90) - target_link_libraries(test_fortran_get_thermo PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_get_thermo PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranGetThermo COMMAND test_fortran_get_thermo) add_executable(test_fortran_box wrap_box.cpp test_fortran_box.f90) - target_link_libraries(test_fortran_box PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_box PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranBox COMMAND test_fortran_box) add_executable(test_fortran_properties wrap_properties.cpp test_fortran_properties.f90 test_fortran_commands.f90) - target_link_libraries(test_fortran_properties PRIVATE flammps lammps MPI::MPI_Fortran GTest::GMockMain) + target_link_libraries(test_fortran_properties PRIVATE flammps lammps GTest::GMockMain) add_test(NAME FortranProperties COMMAND test_fortran_properties) add_executable(test_fortran_extract_global wrap_extract_global.cpp test_fortran_extract_global.f90) - target_link_libraries(test_fortran_extract_global PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_extract_global PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranExtractGlobal COMMAND test_fortran_extract_global) add_executable(test_fortran_extract_atom wrap_extract_atom.cpp test_fortran_extract_atom.f90) - target_link_libraries(test_fortran_extract_atom PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_extract_atom PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranExtractAtom COMMAND test_fortran_extract_atom) add_executable(test_fortran_extract_compute wrap_extract_compute.cpp test_fortran_extract_compute.f90) - target_link_libraries(test_fortran_extract_compute PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_extract_compute PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranExtractCompute COMMAND test_fortran_extract_compute) add_executable(test_fortran_extract_fix wrap_extract_fix.cpp test_fortran_extract_fix.f90) - target_link_libraries(test_fortran_extract_fix PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_extract_fix PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranExtractFix COMMAND test_fortran_extract_fix) add_executable(test_fortran_extract_variable wrap_extract_variable.cpp test_fortran_extract_variable.f90) target_compile_definitions(test_fortran_extract_variable PRIVATE -DTEST_INPUT_FOLDER=${CMAKE_CURRENT_SOURCE_DIR}) - target_link_libraries(test_fortran_extract_variable PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_extract_variable PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranExtractVariable COMMAND test_fortran_extract_variable) add_executable(test_fortran_gather_scatter wrap_gather_scatter.cpp test_fortran_gather_scatter.f90) - target_link_libraries(test_fortran_gather_scatter PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_gather_scatter PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranGatherScatter COMMAND test_fortran_gather_scatter) add_executable(test_fortran_create_atoms wrap_create_atoms.cpp test_fortran_create_atoms.f90) - target_link_libraries(test_fortran_create_atoms PRIVATE flammps lammps MPI::MPI_Fortran GTest::GTestMain) + target_link_libraries(test_fortran_create_atoms PRIVATE flammps lammps GTest::GTestMain) add_test(NAME FortranCreateAtoms COMMAND test_fortran_create_atoms) add_executable(test_fortran_configuration wrap_configuration.cpp test_fortran_configuration.f90 test_fortran_commands.f90) - target_link_libraries(test_fortran_configuration PRIVATE flammps lammps MPI::MPI_Fortran GTest::GMockMain) + target_link_libraries(test_fortran_configuration PRIVATE flammps lammps GTest::GMockMain) add_test(NAME FortranConfiguration COMMAND test_fortran_configuration) add_executable(test_fortran_neighlist wrap_neighlist.cpp test_fortran_neighlist.f90) - target_link_libraries(test_fortran_neighlist PRIVATE flammps lammps MPI::MPI_Fortran GTest::GMockMain) + target_link_libraries(test_fortran_neighlist PRIVATE flammps lammps GTest::GMockMain) add_test(NAME FortranNeighlist COMMAND test_fortran_neighlist) add_executable(test_fortran_fixexternal wrap_fixexternal.cpp test_fortran_fixexternal.f90) diff --git a/unittest/fortran/mpi_stubs.f90 b/unittest/fortran/mpi_stubs.f90 deleted file mode 100644 index 8601f436d2..0000000000 --- a/unittest/fortran/mpi_stubs.f90 +++ /dev/null @@ -1,20 +0,0 @@ -MODULE MPI - IMPLICIT NONE - PRIVATE - - INTEGER, PARAMETER :: MPI_COMM_WORLD=0 - INTEGER, PARAMETER :: MPI_SUCCESS=0 - - PUBLIC :: MPI_COMM_WORLD, MPI_SUCCESS, & - mpi_comm_split - -CONTAINS - - SUBROUTINE mpi_comm_split(comm,color,key,newcomm,ierr) - INTEGER, INTENT(in) :: comm,color,key - INTEGER, INTENT(out) :: newcomm,ierr - - newcomm = comm + 1 - ierr = 0 - END SUBROUTINE mpi_comm_split -END MODULE MPI diff --git a/unittest/fortran/test_fortran_create.f90 b/unittest/fortran/test_fortran_create.f90 index 4ea2a33cfe..fed67064f0 100644 --- a/unittest/fortran/test_fortran_create.f90 +++ b/unittest/fortran/test_fortran_create.f90 @@ -1,3 +1,20 @@ +MODULE MYMPI + USE, INTRINSIC :: ISO_C_BINDING, ONLY: c_int + + IMPLICIT NONE + PRIVATE + PUBLIC :: lmp_comm_split + + INTERFACE + FUNCTION lmp_comm_split(color, key) BIND(C,name='create_mpi_comm_split') + IMPORT :: c_int + IMPLICIT NONE + INTEGER(c_int), VALUE, INTENT(IN) :: color, key + INTEGER(c_int) :: lmp_comm_split + END FUNCTION lmp_comm_split + END INTERFACE +END MODULE MYMPI + FUNCTION f_lammps_no_mpi_no_args() BIND(C, name="f_lammps_no_mpi_no_args") USE ISO_C_BINDING, ONLY: c_ptr USE liblammps @@ -25,35 +42,35 @@ END FUNCTION f_lammps_no_mpi_with_args FUNCTION f_lammps_open_no_args() BIND(C, name="f_lammps_open_no_args") USE ISO_C_BINDING, ONLY: c_ptr - USE MPI, ONLY: MPI_COMM_WORLD, mpi_comm_split + USE MYMPI, ONLY: lmp_comm_split USE liblammps USE keepstuff, ONLY: lmp,mycomm IMPLICIT NONE TYPE(c_ptr) :: f_lammps_open_no_args - INTEGER :: color, key, ierr + INTEGER :: color, key color = 1 key = 1 - CALL mpi_comm_split(MPI_COMM_WORLD, color, key, mycomm, ierr) + mycomm = lmp_comm_split(color, key) lmp = lammps(comm=mycomm) f_lammps_open_no_args = lmp%handle END FUNCTION f_lammps_open_no_args FUNCTION f_lammps_open_with_args() BIND(C, name="f_lammps_open_with_args") USE ISO_C_BINDING, ONLY: c_ptr - USE MPI, ONLY: MPI_COMM_WORLD, mpi_comm_split + USE MYMPI, ONLY: lmp_comm_split USE liblammps USE keepstuff, ONLY: lmp,mycomm IMPLICIT NONE TYPE(c_ptr) :: f_lammps_open_with_args - INTEGER :: color, key, ierr + INTEGER :: color, key CHARACTER(len=12), DIMENSION(4), PARAMETER :: args = & [ CHARACTER(len=12) :: 'liblammps', '-log', 'none', '-nocite' ] color = 2 key = 1 - CALL mpi_comm_split(MPI_COMM_WORLD, color, key, mycomm, ierr) + mycomm = lmp_comm_split(color, key) lmp = lammps(args,mycomm) f_lammps_open_with_args = lmp%handle END FUNCTION f_lammps_open_with_args diff --git a/unittest/fortran/wrap_create.cpp b/unittest/fortran/wrap_create.cpp index 04a62e4040..789e9fca6a 100644 --- a/unittest/fortran/wrap_create.cpp +++ b/unittest/fortran/wrap_create.cpp @@ -17,6 +17,14 @@ void f_lammps_close(); int f_lammps_get_comm(); } +// C wrapper to split MPI communicator w/o requiring a Fortran MPI lib +extern "C" int create_mpi_comm_split(int color, int key) +{ + MPI_Comm c_newcomm = MPI_COMM_NULL; + MPI_Comm_split(MPI_COMM_WORLD, color, key, &c_newcomm); + return MPI_Comm_c2f(c_newcomm); +} + TEST(open_no_mpi, no_args) { ::testing::internal::CaptureStdout(); From 126f597e71d90ee27ebfb70f62e3f41ffb37c160 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 11 Dec 2022 18:02:49 -0500 Subject: [PATCH 19/19] need no longer need to Fortran MPI library --- unittest/fortran/CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/unittest/fortran/CMakeLists.txt b/unittest/fortran/CMakeLists.txt index 08bd1b4356..93573eb297 100644 --- a/unittest/fortran/CMakeLists.txt +++ b/unittest/fortran/CMakeLists.txt @@ -87,7 +87,7 @@ if(CMAKE_Fortran_COMPILER) add_test(NAME FortranNeighlist COMMAND test_fortran_neighlist) add_executable(test_fortran_fixexternal wrap_fixexternal.cpp test_fortran_fixexternal.f90) - target_link_libraries(test_fortran_fixexternal PRIVATE flammps lammps MPI::MPI_Fortran GTest::GMockMain) + target_link_libraries(test_fortran_fixexternal PRIVATE flammps lammps GTest::GMockMain) add_test(NAME FortranFixExternal COMMAND test_fortran_fixexternal) else()