lammps/lammps

[BUG] Fortran scatter_atoms() does not update atom positions when using atom_style atomic

jwrm2 opened this issue · 3 comments

jwrm2 commented

Summary

I am attempting to use the Fortran interface. While the coordinates are read from the initial file correctly, any attempt to update the coordinates with scatter_atoms does nothing when using atom_style atomic.

LAMMPS Version and Platform

LAMMPS (7 Feb 2024 - Development - patch_7Feb2024_update1-384-g9ca47a36d2)

develop branch commit 9ca47a3 10th April

Debian Linux Bookworm, Intel i9-13900K

Compiler version: GNU Fortran (Debian 12.2.0-14) 12.2.0

LAMMPS built in path-to-lammps/build/ configured with

cmake ../cmake/ -DBUILD_MPI=no -DBUILD_OMP=no -DPKG_MOLECULE=yes

Fortran programme compiled with

gfortran path-to-lammps/fortran/lammps.f90 test.f90 path-to-lammps/build/liblammps.a -lstdc++

Expected Behavior

I would expect that calling scatter_atoms() to update the atom positions should update the atom positions. When atom_style is full, the positions are updated as expected.

Actual Behavior

When atom_style is atomic, the atom positions are not updated by scatter_atoms().

Steps to Reproduce

Compile LAMMPS with above settings and Fortran file test.f90 as above. Have attached input files in working directory. Run Fortran executable.

Further Information, Files, and Links

test.f90 - short Fortran programme to reproduce the bug

PROGRAM TEST

    USE LIBLAMMPS
    USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_DOUBLE

    IMPLICIT NONE

    TYPE(LAMMPS) :: LMP
    REAL(KIND=C_DOUBLE), ALLOCATABLE, DIMENSION(:) :: COORDS

    ! Read original coordinates from start_atomic
    LMP = LAMMPS([CHARACTER(LEN=10) :: 'liblammps', '-log', 'lammps.log', '-screen', 'none'])
    CALL LMP%FILE('in_atomic.txt')
    CALL LMP%GATHER_ATOMS('x', 3, COORDS)
    ! Prints 1.0 1.0 1.0 as expected
    PRINT *, COORDS(:)

    ! Try to update coordinates
    COORDS = (/2.0, 2.0, 2.0/)
    CALL LMP%SCATTER_ATOMS('x', COORDS)

    ! Read updated coordinates back
    COORDS(:) = 0.0
    CALL LMP%GATHER_ATOMS('x', 3, COORDS)
    ! Prints 1.0 1.0 1.0, but 2.0 2.0 2.0 expected
    PRINT *, COORDS(:)
    CALL LMP%CLOSE()

    ! Try again with atom_style full
    LMP = LAMMPS([CHARACTER(LEN=10) :: 'liblammps', '-log', 'lammps.log', '-screen', 'none'])
    CALL LMP%FILE('in_full.txt')
    CALL LMP%GATHER_ATOMS('x', 3, COORDS)
    ! Prints 1.0 1.0 1.0 as expected
    PRINT *, COORDS(:)

    ! Try to update coordinates
    COORDS = (/2.0, 2.0, 2.0/)
    CALL LMP%SCATTER_ATOMS('x', COORDS)

    ! Read updated coordinates back
    COORDS(:) = 0.0
    CALL LMP%GATHER_ATOMS('x', 3, COORDS)
    ! Prints 2.0 2.0 2.0 as expected
    PRINT *, COORDS(:)

END PROGRAM TEST

Output from running above programme

   1.0000000000000000        1.0000000000000000        1.0000000000000000     
   1.0000000000000000        1.0000000000000000        1.0000000000000000     
   1.0000000000000000        1.0000000000000000        1.0000000000000000     
   2.0000000000000000        2.0000000000000000        2.0000000000000000

in_atomic.txt
in_full.txt
start_atomic.txt
start_full.txt

akohlmey commented

@jwrm2 This is not a bug. In order to be able to scatter atoms, you need to know the atom's atom IDs, but they are not automatically enabled for atom style atomic (they are required for atom style full).
If you augment your in_atomic.txt file by adding atom_modify map yes at the top it will magically work.

The "scatter" and "gather" functions of the library interface are due for a refactor like the rest of the library interface has seen. In that process this requirement needs to be properly documented and tested for.

If you don't disable screen output you can see:

WARNING: Library error in lammps_scatter_atoms: ids must exist, be consecutive, and be mapped (src/library.cpp:3156)

This warning, for example, should be an error. But since we are now always using exceptions to trap errors, you need to check for error status after a command to see if the command was successful or triggered an error and then you can query for the error message.

jwrm2 commented

OK, thanks for explaining. atom_modify map yes does indeed make it work.

When I was trying to work out what was going on, I did check for an error condition after the scatter call, but LMP%HAS_ERROR() returns .FALSE..

akohlmey commented

When I was trying to work out what was going on, I did check for an error condition after the scatter call, but LMP%HAS_ERROR() returns .FALSE..

Yes. Because this part of the library interface has not yet been refactored to be consistent with the rest (as I pointed out). That would include turning the warning message into an error. At the moment, the screen output is the only indication of this condition.

  • 👍1