From bada1fb348b5a1bbd1d608b98d68f3db8c7053d6 Mon Sep 17 00:00:00 2001 From: Karl Hammond Date: Sun, 25 Sep 2022 23:54:18 -0500 Subject: [PATCH] Finished extract_compute and its unit tests and documentation --- doc/src/Fortran.rst | 143 +++++++++++++++--- fortran/lammps.f90 | 110 +++++++++----- .../fortran/test_fortran_extract_compute.f90 | 88 ++++++++++- 3 files changed, 276 insertions(+), 65 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 7726a43f6e..b24eb30d83 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -452,7 +452,7 @@ Procedures Bound to the lammps Derived Type .. note:: - The `MPI_F08` module, which defines Fortran 2008 bindings for MPI, + The ``MPI_F08`` module, which defines Fortran 2008 bindings for MPI, is not directly supported by this function. However, you should be able to convert between the two using the `MPI_VAL` member of the communicator. For example, @@ -461,12 +461,12 @@ Procedures Bound to the lammps Derived Type USE MPI_F08 USE LIBLAMMPS - TYPE (LAMMPS) :: lmp + TYPE (lammps) :: lmp TYPE (MPI_Comm) :: comm ! ... [commands to set up LAMMPS/etc.] comm%MPI_VAL = lmp%get_mpi_comm() - should assign an `MPI_F08` communicator properly. + should assign an ``MPI_F08`` communicator properly. -------- @@ -505,12 +505,12 @@ Procedures Bound to the lammps Derived Type .. code-block:: fortran PROGRAM demo - USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int64_t + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_int64_t, c_int, c_double USE LIBLAMMPS TYPE(lammps) :: lmp - INTEGER(C_int), POINTER :: nlocal => NULL() - INTEGER(C_int64_t), POINTER :: ntimestep => NULL() - REAL(C_double), POINTER :: dt => NULL() + INTEGER(c_int), POINTER :: nlocal => NULL() + INTEGER(c_int64_t), POINTER :: ntimestep => NULL() + REAL(c_double), POINTER :: dt => NULL() CHARACTER(LEN=10) :: units lmp = lammps() ! other commands @@ -540,7 +540,7 @@ Procedures Bound to the lammps Derived Type pointer (e.g., ``INTEGER (c_int), POINTER :: nlocal``) to the extracted property. If expecting vector data, the pointer should have dimension ":". -.. warning:: + .. warning:: Modifying the data in the location pointed to by the returned pointer may lead to inconsistent internal data and thus may cause failures, @@ -560,14 +560,14 @@ Procedures Bound to the 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 - 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. Pointers returned by this - function are generally persistent; therefore, it is not necessary to call - the function again unless the underlying LAMMPS data are destroyed, such as - through the :doc:`clear` command. + rank (e.g., ``INTEGER(c_int), DIMENSION(:)`` for "type", "mask", or "tag"; + ``INTEGER(c_int64_t), DIMENSION(:)`` for "tag" 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. Pointers returned by + this function are generally persistent; therefore, it is not necessary to + call the function again unless the underlying LAMMPS data are destroyed, + such as through the :doc:`clear` command. :p character(len=\*) name: string with the name of the property to extract :r polymorphic: pointer to LAMMPS data. The left-hand side of the assignment @@ -576,7 +576,7 @@ Procedures Bound to the lammps Derived Type property. If expecting vector data, the pointer should have dimension ":"; if expecting matrix data, the pointer should have dimension ":,:". - .. note:: + .. admonition:: Array index order Two-dimensional arrays returned from :f:func:`extract_atom` will be **transposed** from equivalent arrays in C, and they will be indexed @@ -596,12 +596,12 @@ Procedures Bound to the lammps Derived Type .. code-block:: Fortran TYPE(lammps) :: lmp - REAL(C_double), DIMENSION(:,:), POINTER :: x + REAL(c_double), DIMENSION(:,:), POINTER :: x => NULL() ! more code to setup, etc. x = lmp%extract_atom("x") print '(f0.6)', x(2,6) - will print the *y*-coordinate of the third atom on this processor + will print the *y*-coordinate of the sixth atom on this processor (note the transposition of the two indices). This is not a choice, but rather a consequence of the different conventions adopted by the Fortran and C standards decades ago. @@ -612,7 +612,7 @@ Procedures Bound to the lammps Derived Type .. code-block:: Fortran - REAL(C_double), DIMENSION(:,:), POINTER :: x, x0 + REAL(c_double), DIMENSION(:,:), POINTER :: x, x0 x = lmp%extract_atom("x") x0(0:,0:) => x @@ -621,6 +621,109 @@ Procedures Bound to the lammps Derived Type -------- +.. f:function:: extract_compute(id, style, type) + + This function calls :c:func:`lammps_extract_compute` and returns a pointer + to LAMMPS data tied to the :cpp:class:`Compute` class, specifically data + provided by the compute identified by *id*. Computes may provide global, + per-atom, or local data, and those data may be a scalar, a vector, or an + array. Since computes may provide multiple kinds of data, the user is + required to specify which set of data is to be returned through the + *style* and *type* variables. + + Note that this function actually does not return a value, but rather + 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 compatible rank. The pointer + being associated with LAMMPS data is type-, kind-, and rank-checked at + run-time via an overloaded assignment operator. + + For example, + + .. code-block:: Fortran + + TYPE(lammps) :: lmp + REAL(c_double), DIMENSION(:), POINTER :: COM + ! code to setup, create atoms, etc. + CALL lmp%compute('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 + your simulation. The vector in this case has length 3; the length (or, in + the case of array data, the number of rows and columns) is determined for + you based on data from the :cpp:class:`Compute` class. + + .. admonition:: Array index order + + Two-dimensional arrays returned from :f:func:`extract_compute` will be + **transposed** from equivalent arrays in C, and they will be indexed + from 1 instead of 0. See the similar note under + :f:func:`extract_atom` for further details. + + The following combinations are possible (assuming ``lmp`` is the name of + your LAMMPS instance): + + .. list-table:: + :header-rows: 1 + :widths: auto + + * - Style + - Type + - Pointer type to assign to + - Returned data + * - ``lmp%style%global`` + - ``lmp%type%scalar`` + - ``REAL(c_double), POINTER`` + - Global scalar + * - ``lmp%style%global`` + - ``lmp%type%vector`` + - ``REAL(c_double), DIMENSION(:), POINTER`` + - Global vector + * - ``lmp%style%global`` + - ``lmp%type%array`` + - ``REAL(c_double), DIMENSION(:,:), POINTER`` + - Global array + * - ``lmp%style%atom`` + - ``lmp%type%vector`` + - ``REAL(c_double), DIMENSION(:), POINTER`` + - Per-atom vector + * - ``lmp%style%atom`` + - ``lmp%type%array`` + - ``REAL(c_double), DIMENSION(:,:), POINTER`` + - Per-atom array + * - ``lmp%style%local`` + - ``lmp%type%vector`` + - ``REAL(c_double), DIMENSION(:), POINTER`` + - Local vector + * - ``lmp%style%local`` + - ``lmp%type%array`` + - ``REAL(c_double), DIMENSION(:,:), POINTER`` + - Local array + + :p character(len=\*) id: compute ID from which to extract data + :p integer(c_int) style: value indicating the style of data to extract + (global, per-atom, or local) + :p integer(c_int) type: value indicating the type of data to extract + (scalar, vector, or array) + + .. note:: + + If the compute's data are not already computed for the current step, the + compute will be invoked. LAMMPS cannot easily check at that time if it is + valid to invoke a compute, so it may fail with an error. The caller has + to check to avoid such an error. + + .. warning:: + + The pointers returned by this function are generally not persistent, + since the computed data may be re-distributed, re-allocated, and + re-ordered at every invocation. It is advisable to re-invoke this + function before the data are accessed or make a copy if the data are to + be used after other LAMMPS commands have been issued. Do **not** modify + the data returned by this function. + +-------- + .. f:function:: version() This method returns the numeric LAMMPS version like diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index d06de746d5..4486bd87c1 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -29,8 +29,9 @@ ! MODULE LIBLAMMPS - USE, INTRINSIC :: ISO_C_BINDING, ONLY: c_ptr, c_null_ptr, c_loc, & - c_int, c_int64_t, c_char, c_null_char, c_double, c_size_t, c_f_pointer + 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 USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : ERROR_UNIT IMPLICIT NONE @@ -44,22 +45,27 @@ MODULE LIBLAMMPS ! ! These are NOT part of the API (the part the user sees) INTEGER (c_int), PARAMETER :: & - LAMMPS_INT = 0, & ! 32-bit integer (array) - LAMMPS_INT_2D = 1, & ! two-dimensional 32-bit integer array - LAMMPS_DOUBLE = 2, & ! 64-bit double (array) - LAMMPS_DOUBLE_2D = 3, & ! two-dimensional 64-bit double array - LAMMPS_INT64 = 4, & ! 64-bit integer (array) - LAMMPS_INT64_2D = 5, & ! two-dimensional 64-bit integer array - LAMMPS_STRING = 6, & ! C-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 - LMP_TYPE_SCALAR = 0, & ! request scalar - LMP_TYPE_VECTOR = 1, & ! request vector - 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) + LAMMPS_INT = 0, & ! 32-bit integer (array) + LAMMPS_INT_2D = 1, & ! two-dimensional 32-bit integer array + LAMMPS_DOUBLE = 2, & ! 64-bit double (array) + LAMMPS_DOUBLE_2D = 3, & ! two-dimensional 64-bit double array + LAMMPS_INT64 = 4, & ! 64-bit integer (array) + LAMMPS_INT64_2D = 5, & ! two-dimensional 64-bit integer array + LAMMPS_STRING = 6, & ! C-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 + LMP_TYPE_SCALAR = 0, & ! request scalar + LMP_TYPE_VECTOR = 1, & ! request vector + 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_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) + LMP_ERROR_WORLD = 4, & ! error on comm->world + LMP_ERROR_UNIVERSE = 8 ! error on comm->universe ! "Constants" to use with extract_compute and friends TYPE lammps_style @@ -112,13 +118,13 @@ MODULE LIBLAMMPS ! pointers) TYPE lammps_data INTEGER(c_int) :: datatype - INTEGER(c_int), POINTER :: i32 - INTEGER(c_int), DIMENSION(:), POINTER :: i32_vec - INTEGER(c_int64_t), POINTER :: i64 - INTEGER(c_int64_t), DIMENSION(:), POINTER :: i64_vec - REAL(c_double), POINTER :: r64 - REAL(c_double), DIMENSION(:), POINTER :: r64_vec - REAL(c_double), DIMENSION(:,:), POINTER :: r64_mat + INTEGER(c_int), POINTER :: i32 => NULL() + INTEGER(c_int), DIMENSION(:), POINTER :: i32_vec => NULL() + INTEGER(c_int64_t), POINTER :: i64 => NULL() + INTEGER(c_int64_t), DIMENSION(:), POINTER :: i64_vec => NULL() + REAL(c_double), POINTER :: r64 => NULL() + REAL(c_double), DIMENSION(:), POINTER :: r64_vec => NULL() + REAL(c_double), DIMENSION(:,:), POINTER :: r64_mat => NULL() CHARACTER(LEN=:), ALLOCATABLE :: str END TYPE lammps_data @@ -683,9 +689,9 @@ CONTAINS FORALL ( I=1:length ) global_data%str(i:i) = Fptr(i) END FORALL - CASE DEFAULT - ! FIXME convert to use symbolic constants later - CALL lmp_error(self, 6, 'Unknown pointer type in extract_global') + CASE DEFAULT + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Unknown pointer type in extract_global') END SELECT END FUNCTION @@ -703,6 +709,7 @@ CONTAINS INTEGER :: nrows, ncols REAL(c_double), DIMENSION(:), POINTER :: dummy TYPE(c_ptr), DIMENSION(:), POINTER :: Catomptr + CHARACTER(LEN=:), ALLOCATABLE :: error_msg nmax = lmp_extract_setting(self, 'nmax') ntypes = lmp_extract_setting(self, 'ntypes') @@ -745,13 +752,18 @@ CONTAINS ! Catomptr(1) now points to the first element of the array CALL C_F_POINTER(Catomptr(1), peratom_data%r64_mat, [nrows,ncols]) CASE (-1) - WRITE(ERROR_UNIT,'(A)') 'ERROR: per-atom property "' // name // & - '" not found.' - STOP 2 + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'per-atom property ' // name // 'not found in extract_setting') +! WRITE(ERROR_UNIT,'(A)') 'ERROR: per-atom property "' // name // & +! '" not found.' +! STOP 2 CASE DEFAULT - WRITE(ERROR_UNIT,'(A,I0,A)') 'ERROR: return value ', datatype, & - ' from lammps_extract_atom_datatype not known' - STOP 1 + WRITE(error_msg,'(A,I0,A)') 'return value ', datatype, & + ' from lammps_extract_atom_datatype not known [Fortran/extract_atom]' + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, error_msg) +! WRITE(ERROR_UNIT,'(A,I0,A)') 'ERROR: return value ', datatype, & +! ' from lammps_extract_atom_datatype not known' +! STOP 1 END SELECT END FUNCTION lmp_extract_atom @@ -768,18 +780,26 @@ CONTAINS TYPE(c_ptr) :: Cid, Cptr, Ctemp INTEGER :: nrows, ncols, length INTEGER(c_int), POINTER :: temp + TYPE(c_ptr), DIMENSION(:), POINTER :: Ccomputeptr Cid = f2c_string(id) Cptr = lammps_extract_compute(self%handle, Cid, style, type) + IF ( .NOT. C_ASSOCIATED(Cptr) ) THEN + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Pointer from LAMMPS is NULL [Fortran/extract_compute]') + END IF + ! Remember that rows and columns in C are transposed in Fortran! SELECT CASE (type) CASE (LMP_TYPE_SCALAR) + compute_data%datatype = DATA_DOUBLE length = 1 nrows = 1 ncols = 1 CALL C_F_POINTER(Cptr, compute_data%r64) CASE (LMP_TYPE_VECTOR) + compute_data%datatype = DATA_DOUBLE_1D IF (style == LMP_STYLE_ATOM) THEN length = self%extract_setting('nmax') ELSE @@ -789,11 +809,12 @@ CONTAINS END IF CALL C_F_POINTER(Cptr, compute_data%r64_vec, [length]) CASE (LMP_TYPE_ARRAY) + compute_data%datatype = DATA_DOUBLE_2D IF (style == LMP_STYLE_ATOM) THEN - nrows = self%extract_setting('nmax') - Ctemp = lammps_extract_compute(self%handle,Cid,style,LMP_SIZE_ROWS) + ncols = self%extract_setting('nmax') + Ctemp = lammps_extract_compute(self%handle,Cid,style,LMP_SIZE_COLS) CALL C_F_POINTER(Ctemp, temp) - ncols = temp + nrows = temp ELSE Ctemp = lammps_extract_compute(self%handle,Cid,style,LMP_SIZE_ROWS) CALL C_F_POINTER(Ctemp, temp) @@ -802,11 +823,16 @@ CONTAINS CALL C_F_POINTER(Ctemp, temp) nrows = temp END IF - CALL C_F_POINTER(Cptr, compute_data%r64_mat, [nrows, ncols]) + ! First, we dereference the void** pointer to point to a void* pointer + CALL C_F_POINTER(Cptr, Ccomputeptr, [ncols]) + ! Ccomputeptr(1) now points to the first element of the array + CALL C_F_POINTER(Ccomputeptr(1), compute_data%r64_mat, [nrows, ncols]) CASE DEFAULT - WRITE(ERROR_UNIT,'(A,I0,A)') 'ERROR: unknown type value ', type, & - 'passed to extract_compute' - STOP 1 + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'unknown type value passed to extract_compute') + !WRITE(ERROR_UNIT,'(A,I0,A)') 'ERROR: unknown type value ', type, & + ! 'passed to extract_compute' + !STOP 1 END SELECT CALL lammps_free(Cid) END FUNCTION lmp_extract_compute @@ -944,6 +970,8 @@ CONTAINS CASE DEFAULT str1 = 'that type' END SELECT + !CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, 'cannot associate ' & + ! // str1 // ' with ' // type2 // ' [Fortran API]') WRITE (ERROR_UNIT,'(A)') 'ERROR (Fortran API): cannot associate ' & // str1 // ' with ' // type2 STOP ERROR_CODE diff --git a/unittest/fortran/test_fortran_extract_compute.f90 b/unittest/fortran/test_fortran_extract_compute.f90 index f358ec9d2c..64b5068171 100644 --- a/unittest/fortran/test_fortran_extract_compute.f90 +++ b/unittest/fortran/test_fortran_extract_compute.f90 @@ -6,10 +6,11 @@ MODULE keepcompute 'region box block 0 $x 0 3 0 4', & 'create_box 1 box', & 'create_atoms 1 single 1.0 1.0 ${zpos}' ] - CHARACTER(LEN=40), DIMENSION(2), PARAMETER :: cont_input = & + CHARACTER(LEN=40), DIMENSION(3), PARAMETER :: cont_input = & [ CHARACTER(len=40) :: & 'create_atoms 1 single &', & - ' 0.2 0.1 0.1' ] + ' 0.2 0.1 0.1', & + 'create_atoms 1 single 0.5 0.5 0.5' ] CHARACTER(LEN=40), DIMENSION(3), PARAMETER :: pair_input = & [ CHARACTER(LEN=40) :: & 'pair_style lj/cut 2.5', & @@ -51,11 +52,13 @@ SUBROUTINE f_lammps_setup_extract_compute () BIND(C) CALL lmp%commands_list(pair_input) CALL lmp%command("compute peratompe all pe/atom") ! per-atom vector call lmp%command("compute stress all stress/atom thermo_temp") ! per-atom array - CALL lmp%command("compute COM all com") ! global vector CALL lmp%command("compute totalpe all reduce sum c_peratompe") ! global scalar - CALL lmp%command("compute 1 all rdf 100") ! global array + CALL lmp%command("compute COM all com") ! global vector + CALL lmp%command("compute RDF all rdf 100") ! global array CALL lmp%command("compute pairdist all pair/local dist") ! local vector CALL lmp%command("compute pairlocal all pair/local dist dx dy dz") ! local array + CALL lmp%command("thermo_style custom step pe c_totalpe c_COM[1]") + CALL lmp%command("run 0") ! must be here, otherwise will SEGFAULT END SUBROUTINE f_lammps_setup_extract_compute FUNCTION f_lammps_extract_compute_peratom_vector (i) BIND(C) @@ -70,3 +73,80 @@ FUNCTION f_lammps_extract_compute_peratom_vector (i) BIND(C) vector = lmp%extract_compute('peratompe', lmp%style%atom, lmp%type%vector) f_lammps_extract_compute_peratom_vector = vector(i) END FUNCTION f_lammps_extract_compute_peratom_vector + +FUNCTION f_lammps_extract_compute_peratom_array (i,j) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_double, C_int + USE LIBLAMMPS + USE keepcompute, ONLY : lmp + IMPLICIT NONE + INTEGER(C_int), INTENT(IN), VALUE :: i, j + REAL(C_double) :: f_lammps_extract_compute_peratom_array + REAL(C_double), DIMENSION(:,:), POINTER :: array => NULL() + + array = lmp%extract_compute('stress', lmp%style%atom, lmp%type%array) + f_lammps_extract_compute_peratom_array = array(i,j) +END FUNCTION f_lammps_extract_compute_peratom_array + +FUNCTION f_lammps_extract_compute_global_scalar () BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_double, C_int + USE LIBLAMMPS + USE keepcompute, ONLY : lmp + IMPLICIT NONE + REAL(C_double) :: f_lammps_extract_compute_global_scalar + REAL(C_double), POINTER :: scalar + + scalar = lmp%extract_compute('totalpe', lmp%style%global, lmp%type%scalar) + f_lammps_extract_compute_global_scalar = scalar +END FUNCTION f_lammps_extract_compute_global_scalar + +FUNCTION f_lammps_extract_compute_global_vector (i) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_double, C_int + USE LIBLAMMPS + USE keepcompute, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i + REAL(C_double) :: f_lammps_extract_compute_global_vector + REAL(C_double), DIMENSION(:), POINTER :: vector + + vector = lmp%extract_compute('COM', lmp%style%global, lmp%type%vector) + f_lammps_extract_compute_global_vector = vector(i) +END FUNCTION f_lammps_extract_compute_global_vector + +FUNCTION f_lammps_extract_compute_global_array (i,j) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_double, C_int + USE LIBLAMMPS + USE keepcompute, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i, j + REAL(C_double) :: f_lammps_extract_compute_global_array + REAL(C_double), DIMENSION(:,:), POINTER :: array + + array = lmp%extract_compute('RDF', lmp%style%global, lmp%type%array) + f_lammps_extract_compute_global_array = array(i,j) +END FUNCTION f_lammps_extract_compute_global_array + +FUNCTION f_lammps_extract_compute_local_vector (i) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_double, C_int + USE LIBLAMMPS + USE keepcompute, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i + REAL(C_double) :: f_lammps_extract_compute_local_vector + REAL(C_double), DIMENSION(:), POINTER :: vector + + vector = lmp%extract_compute('pairdist', lmp%style%local, lmp%type%vector) + f_lammps_extract_compute_local_vector = vector(i) +END FUNCTION f_lammps_extract_compute_local_vector + +FUNCTION f_lammps_extract_compute_local_array (i, j) BIND(C) + USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_double, C_int + USE LIBLAMMPS + USE keepcompute, ONLY : lmp + IMPLICIT NONE + INTEGER(c_int), INTENT(IN), VALUE :: i, j + REAL(C_double) :: f_lammps_extract_compute_local_array + REAL(C_double), DIMENSION(:,:), POINTER :: array + + array = lmp%extract_compute('pairlocal', lmp%style%local, lmp%type%array) + f_lammps_extract_compute_local_array = array(i,j) +END FUNCTION f_lammps_extract_compute_local_array