diff --git a/.github/CODEOWNERS b/.github/CODEOWNERS index 59a937faba..d1617cbb25 100644 --- a/.github/CODEOWNERS +++ b/.github/CODEOWNERS @@ -50,6 +50,7 @@ src/PTM/* @pmla src/QMMM/* @akohlmey src/REACTION/* @jrgissing src/REAXFF/* @hasanmetin @stanmoore1 +src/RHEO/* @jtclemm src/SCAFACOS/* @rhalver src/SNAP/* @athomps src/SPIN/* @julient31 diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt index 224dc5cf9d..eac2cc87c2 100644 --- a/cmake/CMakeLists.txt +++ b/cmake/CMakeLists.txt @@ -887,7 +887,7 @@ endif() include(Testing) include(CodeCoverage) include(CodingStandard) -find_package(ClangFormat 8.0) +find_package(ClangFormat 11.0) if(ClangFormat_FOUND) add_custom_target(format-src diff --git a/cmake/Modules/FindClangFormat.cmake b/cmake/Modules/FindClangFormat.cmake index 3f0257f34f..24fad5b6a8 100644 --- a/cmake/Modules/FindClangFormat.cmake +++ b/cmake/Modules/FindClangFormat.cmake @@ -1,5 +1,8 @@ # Find clang-format find_program(ClangFormat_EXECUTABLE NAMES clang-format + clang-format-20.0 + clang-format-19.0 + clang-format-18.0 clang-format-17.0 clang-format-16.0 clang-format-15.0 diff --git a/cmake/Modules/LAMMPSUtils.cmake b/cmake/Modules/LAMMPSUtils.cmake index 2ec9d1b706..223577fe31 100644 --- a/cmake/Modules/LAMMPSUtils.cmake +++ b/cmake/Modules/LAMMPSUtils.cmake @@ -32,7 +32,13 @@ function(check_omp_h_include) set(CMAKE_REQUIRED_INCLUDES ${OpenMP_CXX_INCLUDE_DIRS}) set(CMAKE_REQUIRED_LINK_OPTIONS ${OpenMP_CXX_FLAGS}) set(CMAKE_REQUIRED_LIBRARIES ${OpenMP_CXX_LIBRARIES}) - check_include_file_cxx(omp.h _have_omp_h) + # there are all kinds of problems with finding omp.h + # for Clang and derived compilers so we pretend it is there. + if(CMAKE_CXX_COMPILER_ID STREQUAL "Clang") + set(_have_omp_h TRUE) + else() + check_include_file_cxx(omp.h _have_omp_h) + endif() else() set(_have_omp_h FALSE) endif() diff --git a/cmake/Modules/Packages/RHEO.cmake b/cmake/Modules/Packages/RHEO.cmake index be8c22877b..7639acd8bc 100644 --- a/cmake/Modules/Packages/RHEO.cmake +++ b/cmake/Modules/Packages/RHEO.cmake @@ -1,2 +1,2 @@ -find_package(GSL 2.7 REQUIRED) +find_package(GSL 2.6 REQUIRED) target_link_libraries(lammps PRIVATE GSL::gsl) diff --git a/doc/src/Build_extras.rst b/doc/src/Build_extras.rst index 4802c67420..eae247d66a 100644 --- a/doc/src/Build_extras.rst +++ b/doc/src/Build_extras.rst @@ -639,6 +639,9 @@ They must be specified in uppercase. * - AMD_GFX1100 - GPU - AMD GPU RX7900XTX + * - AMD_GFX1103 + - GPU + - AMD Phoenix APU with Radeon 740M/760M/780M/880M/890M * - INTEL_GEN - GPU - SPIR64-based devices, e.g. Intel GPUs, using JIT diff --git a/doc/src/Build_settings.rst b/doc/src/Build_settings.rst index 34100871ce..b6aa6e4114 100644 --- a/doc/src/Build_settings.rst +++ b/doc/src/Build_settings.rst @@ -414,8 +414,8 @@ Read or write compressed files If this option is enabled, large files can be read or written with compression by ``gzip`` or similar tools by several LAMMPS commands, including :doc:`read_data `, :doc:`rerun `, and -:doc:`dump `. Supported compression tools are currently -``gzip``, ``bzip2``, ``zstd``, and ``lzma``. +:doc:`dump `. Supported compression tools and algorithms are currently +``gzip``, ``bzip2``, ``zstd``, ``xz``, ``lz4``, and ``lzma`` (via xz). .. tabs:: diff --git a/doc/src/Howto_lammps_gui.rst b/doc/src/Howto_lammps_gui.rst index a91a46ad3c..1283bcfa5a 100644 --- a/doc/src/Howto_lammps_gui.rst +++ b/doc/src/Howto_lammps_gui.rst @@ -1,5 +1,5 @@ -Using the LAMMPS-GUI -==================== +Using LAMMPS-GUI +================ This document describes **LAMMPS-GUI version 1.6**. @@ -16,54 +16,101 @@ to the online LAMMPS documentation for known LAMMPS commands and styles. .. note:: - Pre-compiled, ready-to-use LAMMPS-GUI executables for Linux (Ubuntu - 20.04LTS or later and compatible), macOS (version 11 aka Big Sur or - later), and Windows (version 10 or later) :ref:`are available - ` for download. None-MPI LAMMPS executables of - the same LAMMPS version are included in these packages as well. The - source code for the LAMMPS-GUI is included in the LAMMPS source code + Pre-compiled, ready-to-use LAMMPS-GUI executables for Linux x86\_64 + (Ubuntu 20.04LTS or later and compatible), macOS (version 11 aka Big + Sur or later), and Windows (version 10 or later) :ref:`are available + ` for download. None-MPI LAMMPS executables for + running LAMMPS from the command line and :doc:`some LAMMPS tools ` + are also included. + + The source code for LAMMPS-GUI is included in the LAMMPS source code distribution and can be found in the ``tools/lammps-gui`` folder. It can be compiled alongside LAMMPS when :doc:`compiling with CMake `. LAMMPS-GUI tries to provide an experience similar to what people -traditionally would do to run LAMMPS using a command line window -but just rolled into a single executable: +traditionally would have running LAMMPS using a command line window +and the console LAMMPS executable but just rolled into a single executable: -- editing LAMMPS input files with a text editor +- writing & editing LAMMPS input files with a text editor - run LAMMPS on those input file with selected command line flags - use or extract data from the created files and visualize it with either a molecular visualization program or a plotting program - That procedure is quite effective for people proficient in using the command line, as that allows them to use tools for the individual steps that they are most comfortable with. It is often *required* to adopt this workflow when running LAMMPS simulations on high-performance computing facilities. -The main benefit of using the LAMMPS-GUI application instead is that -many basic tasks can be done directly from the GUI without switching to -a text console window or using external programs, let alone writing -scripts to extract data from the generated output. It also integrates -well with graphical desktop environments. +The main benefit of using LAMMPS-GUI is that many basic tasks can be +done directly from the GUI without switching to a text console window or +using external programs, let alone writing scripts to extract data from +the generated output. It also integrates well with graphical desktop +environments where the `.lmp` filename extension can be registered with +LAMMPS-GUI as the executable to launch when double clicking on such +files. Also, LAMMPS-GUI has support for drag-n-drop, i.e. an input +file can be selected and then moved and dropped on the LAMMPS-GUI +executable, and LAMMPS-GUI will launch and read the file into its +buffer. LAMMPS-GUI thus makes it easier for beginners to get started running simple LAMMPS simulations. It is very suitable for tutorials on LAMMPS since you only need to learn how to use a single program for most tasks -and thus time can be saved and people can focus on learning LAMMPS. It -is also designed to keep the barrier low when you decide to switch to a -full featured, standalone programming editor and more sophisticated -visualization and analysis tools, and run LAMMPS from the command line -or a batch script. +and thus time can be saved and people can focus on learning LAMMPS. +The tutorials at https://lammpstutorials.github.io/ were specifically +updated for use with LAMMPS-GUI. + +Another design goal is to keep the barrier low when replacing part of +the functionality of LAMMPS-GUI with external tools. The following text provides a detailed tour of the features and -functionality of the LAMMPS-GUI. Suggestions for new features and +functionality of LAMMPS-GUI. Suggestions for new features and reports of bugs are always welcome. You can use the :doc:`the same channels as for LAMMPS itself ` for that purpose. ----- +Installing Pre-compiled LAMMPS-GUI Packages +------------------------------------------- + +LAMMPS-GUI is available as pre-compiled binary packages for Linux +x86\_64, macOS 11 and later, and Windows 10 and later. Alternately, it +can be compiled from source. + +Windows 10 and later +^^^^^^^^^^^^^^^^^^^^ + +After downloading the ``LAMMPS-Win10-64bit-GUI-.exe`` installer +package, you need to execute it, and start the installation process. +Since those packages are currently unsigned, you have to enable "Developer Mode" +in the Windows System Settings to run the installer. + +MacOS 11 and later +^^^^^^^^^^^^^^^^^^ + +After downloading the ``LAMMPS-macOS-multiarch-GUI-.dmg`` +installer package, you need to double-click it and then, in the window +that opens, drag the app bundle as indicated into the "Applications" +folder. The follow the instructions in the "README.txt" file to +get access to the other included executables. + +Linux on x86\_64 +^^^^^^^^^^^^^^^^ + +After downloading and unpacking the +``LAMMPS-Linux-x86_64-GUI-.tar.gz`` package. You can switch +into the "LAMMPS_GUI" folder and execute "./lammps-gui" directly. + +Compiling from Source +^^^^^^^^^^^^^^^^^^^^^ + +There also are instructions for :ref:`compiling LAMMPS-GUI from source +code ` available elsewhere in the manual. +Compilation from source *requires* using CMake. + +----- + Starting LAMMPS-GUI ------------------- @@ -88,17 +135,24 @@ window is stored when exiting and restored when starting again. Opening Files ^^^^^^^^^^^^^ -The LAMMPS-GUI application tries to open the first command line argument -as a LAMMPS input script, further arguments are ignored. When no -argument is given, LAMMPS-GUI starts with an empty buffer. Files can -also be opened via the ``File`` menu or by drag-and-drop of a file from -a graphical file manager into the editor window. Only one file can be -edited at a time, so opening a new file with a filled buffer closes that -buffer. If the buffer has unsaved modifications, you are asked to -either cancel the operation, discard the changes, or save them. A -buffer with modifications can be saved any time from the "File" menu, by -the keyboard shortcut `Ctrl-S` (`Command-S` on macOS), or by clicking on -the "Save" button at the very left in the status bar. +The LAMMPS-GUI application can be launched without command line arguments +and then starts with an empty buffer in the *Editor* window. If arguments +are given LAMMPS will use first command line argument as the file name for +the *Editor* buffer and reads its contents into the buffer, if the file +exists. All further arguments are ignored. Files can also be opened via +the ``File`` menu, the `Ctrl-O` (`Command-O` on macOS) keyboard shortcut +or by drag-and-drop of a file from a graphical file manager into the editor +window. If a file extension (e.g. ``.lmp``) has been registered with the +graphical environment to launch LAMMPS-GUI, an existing input file can +be launched with LAMMPS-GUI through double clicking. + +Only one file can be edited at a time, so opening a new file with a +filled buffer closes that buffer. If the buffer has unsaved +modifications, you are asked to either cancel the operation, discard the +changes, or save them. A buffer with modifications can be saved any +time from the "File" menu, by the keyboard shortcut `Ctrl-S` +(`Command-S` on macOS), or by clicking on the "Save" button at the very +left in the status bar. Running LAMMPS ^^^^^^^^^^^^^^ @@ -235,20 +289,30 @@ run number that this chart window corresponds to. Same as for the *Output* window, the chart window is replaced on each new run, but the behavior can be changed in the preferences dialog. +.. versionadded:: 1.6 + + Support for YAML export added + From the ``File`` menu on the top left, it is possible to save an image of the currently displayed plot or export the data in either plain text columns (for use by plotting tools like `gnuplot `_ or `grace -`_), or as CSV data which can -be imported for further processing with Microsoft Excel or `pandas -`_ +`_), as CSV data which can be +imported for further processing with Microsoft Excel `LibreOffice Calc +`_ or with Python via `pandas +`_, or as YAML which can be imported into +Python with `PyYAML `_ or pandas. Thermo output data from successive run commands in the input script is combined into a single data set unless the format, number, or names of output columns are changed with a :doc:`thermo_style ` or a :doc:`thermo_modify ` command, or the current time step is reset with :doc:`reset_timestep `, or if a -:doc:`clear ` command is issued. +:doc:`clear ` command is issued. This is where the YAML export +from the *Charts* window differs from that of the *Output* window: +here you get the compounded data set starting with the last change of +output fields or timestep setting, while the export from the log will +contain *all* YAML output but *segmented* into individual runs. Image Slide Show ---------------- @@ -347,15 +411,16 @@ actual image size, high-quality (SSAO) rendering, anti-aliasing, view style, display of box or axes, zoom factor. The view of the system can be rotated horizontally and vertically. It is also possible to only display the atoms within a group defined in the input script (default is -"all"). After each change, the image is rendered again and the display -updated. The small palette icon on the top left is colored while LAMMPS -is running to render the new image; it is grayed out when LAMMPS is -finished. When there are many atoms to render and high quality images -with anti-aliasing are requested, re-rendering may take several seconds. -From the ``File`` menu of the image window, the current image can be -saved to a file (keyboard shortcut `Ctrl-S`) or copied to the clipboard -(keyboard shortcut `Ctrl-C`) for pasting the image into another -application. +"all"). The image can also be re-centered on the center of mass of the +selected group. After each change, the image is rendered again and the +display updated. The small palette icon on the top left is colored +while LAMMPS is running to render the new image; it is grayed out when +LAMMPS is finished. When there are many atoms to render and high +quality images with anti-aliasing are requested, re-rendering may take +several seconds. From the ``File`` menu of the image window, the +current image can be saved to a file (keyboard shortcut `Ctrl-S`) or +copied to the clipboard (keyboard shortcut `Ctrl-C`) for pasting the +image into another application. .. versionadded:: 1.6 @@ -427,7 +492,7 @@ Context Specific Help |gui-popup1| |gui-popup2| -A unique feature of the LAMMPS-GUI is the option to look up the +A unique feature of LAMMPS-GUI is the option to look up the LAMMPS documentation for the command in the current line. This can be done by either clicking the right mouse button or by using the `Ctrl-?` keyboard shortcut. When using the mouse, there are additional entries in the @@ -435,10 +500,16 @@ context menu that open the corresponding documentation page in the online LAMMPS documentation in a web browser window. When using the keyboard, the first of those entries is chosen. +.. versionadded:: 1.6 + If the word under the cursor is a file, then additionally the context menu has an entry to open the file in a read-only text viewer window. This is a convenient way to view the contents of files that are -referenced in the input. +referenced in the input. The file viewer also supports on-the-fly +decompression based on the file name suffix in a :ref:`similar fashion +as available with LAMMPS `. If the necessary decompression +program is missing or the file cannot be decompressed, the viewer window +will contain a corresponding message. Menu ---- @@ -458,7 +529,7 @@ The ``File`` menu offers the usual options: - ``New`` clears the current buffer and resets the file name to ``*unknown*`` - ``Open`` opens a dialog to select a new file for editing in the *Editor* -- ``View`` opens a dialog to select a file for viewing in a *separate* window (read-only) +- ``View`` opens a dialog to select a file for viewing in a *separate* window (read-only) with support for on-the-fly decompression as explained above. - ``Save`` saves the current file; if the file name is ``*unknown*`` a dialog will open to select a new file name - ``Save As`` opens a dialog to select and new file name (and folder, if @@ -531,12 +602,12 @@ in an ``Image Viewer`` window. The ``View in OVITO`` entry will launch `OVITO `_ with a :doc:`data file ` containing the current state of -the system. This option is only available if the LAMMPS-GUI can find +the system. This option is only available if LAMMPS-GUI can find the OVITO executable in the system path. The ``View in VMD`` entry will launch VMD with a :doc:`data file ` containing the current state of the system. This option -is only available if the LAMMPS-GUI can find the VMD executable in the +is only available if LAMMPS-GUI can find the VMD executable in the system path. View @@ -559,6 +630,9 @@ a minimal description of LAMMPS-GUI. The ``LAMMPS-GUI Howto`` entry will open this documentation page from the online documentation in a web browser window. The ``LAMMPS Manual`` entry will open the main page of the LAMMPS online documentation in a web browser window. +The ``LAMMPS Tutorial`` entry will open the main page of the set of +LAMMPS tutorials authored and maintained by Simon Gravelle at +https://lammpstutorials.github.io/ in a web browser window. ----- @@ -566,8 +640,8 @@ Preferences ----------- The ``Preferences`` dialog allows customization of the behavior and -look of the LAMMPS-GUI application. The settings are grouped and each -group is displayed within a tab. +look of LAMMPS-GUI. The settings are grouped and each group is +displayed within a tab. .. |guiprefs1| image:: JPG/lammps-gui-prefs-general.png :width: 24% @@ -744,12 +818,12 @@ available (On macOS use the Command key instead of Ctrl/Control). - Reformat line - Shift+TAB - Show Completions - * - Ctrl+Shift+Enter + * - Ctrl+Shift+T + - LAMMPS Tutorial + - Ctrl+Shift+Enter - Run File - - - - - - Further editing keybindings `are documented with the Qt documentation `_. In diff --git a/doc/src/Howto_pylammps.rst b/doc/src/Howto_pylammps.rst index bce37d5ac7..5ef3248e1d 100644 --- a/doc/src/Howto_pylammps.rst +++ b/doc/src/Howto_pylammps.rst @@ -6,19 +6,22 @@ PyLammps Tutorial Overview -------- -``PyLammps`` is a Python wrapper class for LAMMPS which can be created -on its own or use an existing lammps Python object. It creates a simpler, +:py:class:`PyLammps ` is a Python wrapper class for +LAMMPS which can be created on its own or use an existing +:py:class:`lammps Python ` object. It creates a simpler, more "pythonic" interface to common LAMMPS functionality, in contrast to -the ``lammps`` wrapper for the C-style LAMMPS library interface which -is written using `Python ctypes `_. The ``lammps`` wrapper -is discussed on the :doc:`Python_head` doc page. +the :py:class:`lammps ` wrapper for the LAMMPS :ref:`C +language library interface API ` which is written using +`Python ctypes `_. The :py:class:`lammps ` +wrapper is discussed on the :doc:`Python_head` doc page. -Unlike the flat ``ctypes`` interface, PyLammps exposes a discoverable -API. It no longer requires knowledge of the underlying C++ code -implementation. Finally, the ``IPyLammps`` wrapper builds on top of -``PyLammps`` and adds some additional features for -`IPython integration `_ into `Jupyter notebooks `_, -e.g. for embedded visualization output from :doc:`dump style image `. +Unlike the flat `ctypes `_ interface, PyLammps exposes a +discoverable API. It no longer requires knowledge of the underlying C++ +code implementation. Finally, the :py:class:`IPyLammps +` wrapper builds on top of :py:class:`PyLammps +` and adds some additional features for `IPython +integration `_ into `Jupyter notebooks `_, e.g. for +embedded visualization output from :doc:`dump style image `. .. _ctypes: https://docs.python.org/3/library/ctypes.html .. _ipython: https://ipython.org/ @@ -30,19 +33,22 @@ Comparison of lammps and PyLammps interfaces lammps.lammps """"""""""""" -* uses ``ctypes`` -* direct memory access to native C++ data +* uses `ctypes `_ +* direct memory access to native C++ data with optional support for NumPy arrays * provides functions to send and receive data to LAMMPS +* interface modeled after the LAMMPS :ref:`C language library interface API ` * requires knowledge of how LAMMPS internally works (C pointers, etc) +* full support for running Python with MPI using `mpi4py `_ lammps.PyLammps """"""""""""""" -* higher-level abstraction built on top of original ctypes interface +* higher-level abstraction built on *top* of original :py:class:`ctypes based interface ` * manipulation of Python objects * communication with LAMMPS is hidden from API user * shorter, more concise Python * better IPython integration, designed for quick prototyping +* designed for serial execution Quick Start ----------- @@ -506,14 +512,26 @@ inside of the IPython notebook. Using PyLammps and mpi4py (Experimental) ---------------------------------------- -PyLammps can be run in parallel using mpi4py. This python package can be installed using +PyLammps can be run in parallel using `mpi4py +`_. This python package can be installed +using .. code-block:: bash pip install mpi4py -The following is a short example which reads in an existing LAMMPS input file and -executes it in parallel. You can find in.melt in the examples/melt folder. +.. warning:: + + Usually, any :py:class:`PyLammps ` command must be + executed by *all* MPI processes. However, evaluations and querying + the system state is only available on MPI rank 0. Using these + functions from other MPI ranks will raise an exception. + +The following is a short example which reads in an existing LAMMPS input +file and executes it in parallel. You can find in.melt in the +examples/melt folder. Please take note that the +:py:meth:`PyLammps.eval() ` is called only from +MPI rank 0. .. code-block:: python @@ -535,10 +553,6 @@ following mpirun command: mpirun -np 4 python melt.py -.. warning:: - - Any command must be executed by all MPI processes. However, evaluations and querying the system state is only available on rank 0. - Feedback and Contributing ------------------------- diff --git a/doc/src/JPG/lammps-gui-image.png b/doc/src/JPG/lammps-gui-image.png index fb3fe609e0..48b9fd1ce3 100644 Binary files a/doc/src/JPG/lammps-gui-image.png and b/doc/src/JPG/lammps-gui-image.png differ diff --git a/doc/src/Tools.rst b/doc/src/Tools.rst index 4505142459..060fd9f89d 100644 --- a/doc/src/Tools.rst +++ b/doc/src/Tools.rst @@ -484,23 +484,22 @@ are in the :doc:`Howto_lammps_gui` tutorial Howto page. Here are a few highlights of LAMMPS-GUI -- Text editor with syntax highlighting customized for LAMMPS -- Text editor features command completion for known commands and styles -- Text editor will switch working directory to folder of file in buffer -- Text editor will remember up to 5 recent files +- Text editor with line numbers and syntax highlighting customized for LAMMPS +- Text editor features command completion and auto-indentation for known commands and styles +- Text editor will switch its working directory to folder of file in buffer +- Many adjustable settings and preferences that are persistent including the 5 most recent files - Context specific LAMMPS command help via online documentation - LAMMPS is running in a concurrent thread, so the GUI remains responsive -- Support for most accelerator packages - Progress bar indicates how far a run command is completed -- LAMMPS can be started and stopped with a hotkey -- Screen output is captured in a Log Window -- Thermodynamic output is captured and displayed as line graph in a Chart Window +- LAMMPS can be started and stopped with a mouse click or a hotkey +- Screen output is captured in an *Output* Window +- Thermodynamic output is captured and displayed as line graph in a *Chart* Window - Indicator for currently executed command - Indicator for line that caused an error - Visualization of current state in Image Viewer (via calling :doc:`write_dump image `) - Capture of images created via :doc:`dump image ` in Slide show window -- Many adjustable settings and preferences that are persistent - Dialog to set variables, similar to the LAMMPS command line flag '-v' / '-var' +- Support for GPU, INTEL, KOKKOS/OpenMP, OPENMAP, and OPT and accelerator packages Parallelization ^^^^^^^^^^^^^^^ @@ -542,6 +541,8 @@ variable so the executables will be found automatically. The LAMMPS-GUI executable is called ``lammps-gui`` and either takes no arguments or attempts to load the first argument as LAMMPS input file. +.. _lammps_gui_compilation: + Compilation ^^^^^^^^^^^ diff --git a/doc/src/angle_hybrid.rst b/doc/src/angle_hybrid.rst index 5eba42cd25..fc213c13e9 100644 --- a/doc/src/angle_hybrid.rst +++ b/doc/src/angle_hybrid.rst @@ -82,6 +82,10 @@ for specific angle types. ---------- +.. include:: accel_styles.rst + +---------- + Restrictions """""""""""" @@ -90,8 +94,9 @@ MOLECULE package. See the :doc:`Build package ` doc page for more info. Unlike other angle styles, the hybrid angle style does not store angle -coefficient info for individual sub-styles in a :doc:`binary restart files `. Thus when restarting a simulation from a restart -file, you need to re-specify :doc:`angle_coeff ` commands. +coefficient info for individual sub-styles in :doc:`binary restart files +` or :doc:`data files `. Thus when restarting a +simulation, you need to re-specify the angle_coeff commands. Related commands """""""""""""""" diff --git a/doc/src/bond_hybrid.rst b/doc/src/bond_hybrid.rst index ea90ed2b57..fce393008e 100644 --- a/doc/src/bond_hybrid.rst +++ b/doc/src/bond_hybrid.rst @@ -75,8 +75,9 @@ package. See the :doc:`Build package ` page for more info. Unlike other bond styles, the hybrid bond style does not store bond -coefficient info for individual sub-styles in a :doc:`binary restart files `. Thus when restarting a simulation from a restart -file, you need to re-specify bond_coeff commands. +coefficient info for individual sub-styles in :doc:`binary restart files +` or :doc:`data files `. Thus when restarting a +simulation, you need to re-specify the bond_coeff commands. Related commands """""""""""""""" diff --git a/doc/src/dihedral_hybrid.rst b/doc/src/dihedral_hybrid.rst index b22fecef44..460c329de9 100644 --- a/doc/src/dihedral_hybrid.rst +++ b/doc/src/dihedral_hybrid.rst @@ -83,6 +83,10 @@ for specific dihedral types. ---------- +.. include:: accel_styles.rst + +---------- + Restrictions """""""""""" @@ -91,8 +95,10 @@ MOLECULE package. See the :doc:`Build package ` doc page for more info. Unlike other dihedral styles, the hybrid dihedral style does not store -dihedral coefficient info for individual sub-styles in a :doc:`binary restart files `. Thus when restarting a simulation from a -restart file, you need to re-specify dihedral_coeff commands. +dihedral coefficient info for individual sub-styles in :doc:`binary +restart files ` or :doc:`data files `. Thus when +restarting a simulation, you need to re-specify the dihedral_coeff +commands. Related commands """""""""""""""" diff --git a/doc/src/fix_meso_move.rst b/doc/src/fix_meso_move.rst index 281a405ab0..64b451b7f1 100644 --- a/doc/src/fix_meso_move.rst +++ b/doc/src/fix_meso_move.rst @@ -247,7 +247,7 @@ defined by the :doc:`atom_style sph ` command. All particles in the group must be mesoscopic SPH/SDPD particles. -versionchanged:: TBD +.. versionchanged:: TBD This fix is incompatible with deformation controls that remap velocity, for instance the *remap v* option of :doc:`fix deform `. diff --git a/doc/src/fix_mvv_dpd.rst b/doc/src/fix_mvv_dpd.rst index 740785d562..e64a162bf4 100644 --- a/doc/src/fix_mvv_dpd.rst +++ b/doc/src/fix_mvv_dpd.rst @@ -97,7 +97,7 @@ These fixes are part of the DPD-MESO package. They are only enabled if LAMMPS was built with that package. See the :doc:`Build package ` page for more info. -versionchanged:: TBD +.. versionchanged:: TBD This fix is incompatible with deformation controls that remap velocity, for instance the *remap v* option of :doc:`fix deform `. diff --git a/doc/src/fix_rigid_meso.rst b/doc/src/fix_rigid_meso.rst index 933369788e..3f734e3fef 100644 --- a/doc/src/fix_rigid_meso.rst +++ b/doc/src/fix_rigid_meso.rst @@ -353,7 +353,7 @@ defined by the :doc:`atom_style sph ` command. All particles in the group must be mesoscopic SPH/SDPD particles. -versionchanged:: TBD +.. versionchanged:: TBD This fix is incompatible with deformation controls that remap velocity, for instance the *remap v* option of :doc:`fix deform `. diff --git a/doc/src/fix_smd_integrate_tlsph.rst b/doc/src/fix_smd_integrate_tlsph.rst index e714c91a17..44d4bab3a5 100644 --- a/doc/src/fix_smd_integrate_tlsph.rst +++ b/doc/src/fix_smd_integrate_tlsph.rst @@ -53,7 +53,7 @@ Restrictions This fix is part of the MACHDYN package. It is only enabled if LAMMPS was built with that package. See the :doc:`Build package ` page for more info. -versionchanged:: TBD +.. versionchanged:: TBD This fix is incompatible with deformation controls that remap velocity, for instance the *remap v* option of :doc:`fix deform `. diff --git a/doc/src/fix_smd_integrate_ulsph.rst b/doc/src/fix_smd_integrate_ulsph.rst index 60c185db5f..6b1e070763 100644 --- a/doc/src/fix_smd_integrate_ulsph.rst +++ b/doc/src/fix_smd_integrate_ulsph.rst @@ -61,7 +61,7 @@ Restrictions This fix is part of the MACHDYN package. It is only enabled if LAMMPS was built with that package. See the :doc:`Build package ` page for more info. -versionchanged:: TBD +.. versionchanged:: TBD This fix is incompatible with deformation controls that remap velocity, for instance the *remap v* option of :doc:`fix deform `. diff --git a/doc/src/improper_hybrid.rst b/doc/src/improper_hybrid.rst index c2f80bdd25..d7b9a425f3 100644 --- a/doc/src/improper_hybrid.rst +++ b/doc/src/improper_hybrid.rst @@ -82,6 +82,10 @@ types. ---------- +.. include:: accel_styles.rst + +---------- + Restrictions """""""""""" @@ -90,9 +94,10 @@ MOLECULE package. See the :doc:`Build package ` doc page for more info. Unlike other improper styles, the hybrid improper style does not store -improper coefficient info for individual sub-styles in a :doc:`binary restart files `. -Thus when restarting a simulation from a -restart file, you need to re-specify improper_coeff commands. +improper coefficient info for individual sub-styles in :doc:`binary +restart files ` or :doc:`data files `. Thus when +restarting a simulation, you need to re-specify the improper_coeff +commands. Related commands """""""""""""""" diff --git a/doc/src/pair_born_gauss.rst b/doc/src/pair_born_gauss.rst index 97f0042ee7..58796438ad 100644 --- a/doc/src/pair_born_gauss.rst +++ b/doc/src/pair_born_gauss.rst @@ -19,7 +19,7 @@ Examples .. code-block:: LAMMPS pair_style born/gauss 10.0 - pair_coeff 1 1 1 1 8.2464e13 12.48 0.042644277 0.44 3.56 + pair_coeff 1 1 8.2464e13 12.48 0.042644277 0.44 3.56 Description """"""""""" diff --git a/doc/src/pair_hybrid.rst b/doc/src/pair_hybrid.rst index 634717f8a8..617b0c4372 100644 --- a/doc/src/pair_hybrid.rst +++ b/doc/src/pair_hybrid.rst @@ -479,11 +479,12 @@ For the hybrid pair styles, the list of sub-styles and their respective settings are written to :doc:`binary restart files `, so a :doc:`pair_style ` command does not need to specified in an input script that reads a restart file. However, the coefficient -information is not stored in the restart file. Thus, pair_coeff -commands need to be re-specified in the restart input script. For pair -style *hybrid/scaled* also the names of any variables used as scale -factors are restored, but not the variables themselves, so those may -need to be redefined when continuing from a restart. +information is not stored in the restart file. The same is true for +:doc:`data files `. Thus, pair_coeff commands need to be +re-specified in the restart input script. For pair style +*hybrid/scaled* also the names of any variables used as scale factors +are restored, but not the variables themselves, so those may need to be +redefined when continuing from a restart. These pair styles support the use of the *inner*, *middle*, and *outer* keywords of the :doc:`run_style respa ` command, if diff --git a/doc/src/write_data.rst b/doc/src/write_data.rst index e1e5b04953..937ab5ba02 100644 --- a/doc/src/write_data.rst +++ b/doc/src/write_data.rst @@ -51,10 +51,12 @@ value. The write_data command may not always write all coefficient settings to the corresponding Coeff sections of the data file. This can have - one of multiple reasons. 1) A few styles may be missing the code that - would write those sections (if you come across one, please notify - the LAMMPS developers). 2) Some pair styles require a single pair_coeff - statement and those are not compatible with data files. 3) The + one of multiple reasons. 1) The style may be a hybrid style. In that + case *no* coeff information is written. 2) A few styles may be + missing the code that would write those sections (This is rare these + days, but if you come across one, please notify the LAMMPS + developers). 3) Some pair styles require a single pair_coeff + statement and those are not compatible with data files. 4) The default for write_data is to write a PairCoeff section, which has only entries for atom types i == j. The remaining coefficients would be inferred through the currently selected mixing rule. If there has diff --git a/doc/utils/requirements.txt b/doc/utils/requirements.txt index 83d4999dd7..6787b9d682 100644 --- a/doc/utils/requirements.txt +++ b/doc/utils/requirements.txt @@ -1,4 +1,4 @@ -Sphinx >= 5.3.0, <8.0 +Sphinx >= 5.3.0, <7.5 sphinxcontrib-spelling sphinxcontrib-jquery sphinx-design diff --git a/doc/utils/sphinx-config/false_positives.txt b/doc/utils/sphinx-config/false_positives.txt index d18ab66181..002e861dc7 100644 --- a/doc/utils/sphinx-config/false_positives.txt +++ b/doc/utils/sphinx-config/false_positives.txt @@ -1355,6 +1355,7 @@ Grama grana granregion graphene +Gravelle Greathouse greenyellow Greffet @@ -3080,6 +3081,7 @@ qx qy qz Rackers +Radeon radi radialscreened radialscreenedspin diff --git a/examples/charmmfsw/data.charmmfsw.gz b/examples/charmmfsw/data.charmmfsw.gz index 0fe7d2df0b..8d5741d0f9 100644 Binary files a/examples/charmmfsw/data.charmmfsw.gz and b/examples/charmmfsw/data.charmmfsw.gz differ diff --git a/examples/charmmfsw/in.charmmfsw b/examples/charmmfsw/in.charmmfsw index fe7f486006..204c2627d3 100644 --- a/examples/charmmfsw/in.charmmfsw +++ b/examples/charmmfsw/in.charmmfsw @@ -28,7 +28,16 @@ neighbor 2 bin neigh_modify delay 2 every 1 fix 1 all shake 1e-6 100 100 m 1.008 a 142 -fix 2 all nvt temp 303.15 303.15 100.0 +fix 2 all nvt temp 303.15 303.15 100.0 + +# for visualization with LAMMPS-GUI +group water type 18 60 +group nowater subtract all water +group ions type 63 64 +group other subtract all water ions + +# dump viz other image 10 myimage-*.ppm element type size 800 800 zoom 2.82954 shiny 0.5 fsaa yes bond none none view 20 10 box no 0.0 axes no 0.0 0.0 center s 0.521318 0.489856 0.489856 +# dump_modify viz pad 9 boxcolor darkblue backcolor darkgray element H H H H H H H H H H H H H H H H H H C C C C C C C C C C C C C C C C C C C N N N N N N N N N N N N N N O O O O O O O O O P S Cl K adiam 1 1.92 adiam 2 1.92 adiam 3 1.92 adiam 4 1.92 adiam 5 1.92 adiam 6 1.92 adiam 7 1.92 adiam 8 1.92 adiam 9 1.92 adiam 10 1.92 adiam 11 1.92 adiam 12 1.92 adiam 13 1.92 adiam 14 1.92 adiam 15 1.92 adiam 16 1.92 adiam 17 1.92 adiam 18 1.92 adiam 19 2.72 adiam 20 2.72 adiam 21 2.72 adiam 22 2.72 adiam 23 2.72 adiam 24 2.72 adiam 25 2.72 adiam 26 2.72 adiam 27 2.72 adiam 28 2.72 adiam 29 2.72 adiam 30 2.72 adiam 31 2.72 adiam 32 2.72 adiam 33 2.72 adiam 34 2.72 adiam 35 2.72 adiam 36 2.72 adiam 37 2.72 adiam 38 2.48 adiam 39 2.48 adiam 40 2.48 adiam 41 2.48 adiam 42 2.48 adiam 43 2.48 adiam 44 2.48 adiam 45 2.48 adiam 46 2.48 adiam 47 2.48 adiam 48 2.48 adiam 49 2.48 adiam 50 2.48 adiam 51 2.48 adiam 52 2.432 adiam 53 2.432 adiam 54 2.432 adiam 55 2.432 adiam 56 2.432 adiam 57 2.432 adiam 58 2.432 adiam 59 2.432 adiam 60 2.432 adiam 61 2.88 adiam 62 2.88 adiam 63 3.632 adiam 64 2.816 thermo 10 thermo_style custom step etotal evdwl ecoul elong edihed pe ke temp press diff --git a/lib/gpu/lal_dpd_coul_slater_long.cu b/lib/gpu/lal_dpd_coul_slater_long.cu index d9679ccdfd..4835b8777b 100644 --- a/lib/gpu/lal_dpd_coul_slater_long.cu +++ b/lib/gpu/lal_dpd_coul_slater_long.cu @@ -186,15 +186,13 @@ __kernel void k_dpd_coul_slater_long(const __global numtyp4 *restrict x_, atom_info(t_per_atom,ii,tid,offset); __local numtyp sp_cl[4]; - ///local_allocate_store_charge(); - sp_cl[0]=sp_cl_in[0]; sp_cl[1]=sp_cl_in[1]; sp_cl[2]=sp_cl_in[2]; sp_cl[3]=sp_cl_in[3]; int n_stride; - local_allocate_store_pair(); + local_allocate_store_charge(); acctyp3 f; f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0; @@ -332,8 +330,7 @@ __kernel void k_dpd_coul_slater_long(const __global numtyp4 *restrict x_, } // for nbor } // if ii - store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,vflag, - ans,engv); + store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,vflag,ans,engv); } __kernel void k_dpd_coul_slater_long_fast(const __global numtyp4 *restrict x_, @@ -378,7 +375,7 @@ __kernel void k_dpd_coul_slater_long_fast(const __global numtyp4 *restrict x_, int n_stride; - local_allocate_store_pair(); + local_allocate_store_charge(); acctyp3 f; f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0; @@ -517,7 +514,6 @@ __kernel void k_dpd_coul_slater_long_fast(const __global numtyp4 *restrict x_, } // for nbor } // if ii - store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,vflag, - ans,engv); + store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,vflag,ans,engv); } diff --git a/lib/gpu/lal_sw.cpp b/lib/gpu/lal_sw.cpp index 9687a0352d..60d388f100 100644 --- a/lib/gpu/lal_sw.cpp +++ b/lib/gpu/lal_sw.cpp @@ -52,12 +52,12 @@ int SWT::init(const int ntypes, const int nlocal, const int nall, double ***costheta, const int *map, int ***e2param) { _lj_types=ntypes; - int oldparam=-1; int onetype=-1; int onetype3=0; int spq=1; - int mtypes=0; #ifdef USE_OPENCL + int oldparam=-1; + int mtypes=0; for (int ii=1; ii, local_offset_value = element_values(team_id, i - 1); // FIXME_OPENMPTARGET We seem to access memory illegaly on AMD GPUs #if defined(KOKKOS_ARCH_AMD_GPU) && !defined(KOKKOS_ARCH_AMD_GFX1030) && \ - !defined(KOKKOS_ARCH_AMD_GFX1100) + !defined(KOKKOS_ARCH_AMD_GFX1100) && !defined(KOKKOS_ARCH_AMD_GFX1103) if constexpr (Analysis::Reducer::has_join_member_function()) { if constexpr (std::is_void_v) a_functor_reducer.get_functor().join(local_offset_value, diff --git a/lib/kokkos/core/src/impl/Kokkos_Core.cpp b/lib/kokkos/core/src/impl/Kokkos_Core.cpp index 4a69652616..c7addbe337 100644 --- a/lib/kokkos/core/src/impl/Kokkos_Core.cpp +++ b/lib/kokkos/core/src/impl/Kokkos_Core.cpp @@ -750,6 +750,9 @@ void pre_initialize_internal(const Kokkos::InitializationSettings& settings) { #elif defined(KOKKOS_ARCH_AMD_GFX1100) declare_configuration_metadata("architecture", "GPU architecture", "AMD_GFX1100"); +#elif defined(KOKKOS_ARCH_AMD_GFX1103) + declare_configuration_metadata("architecture", "GPU architecture", + "AMD_GFX1103"); #else declare_configuration_metadata("architecture", "GPU architecture", "none"); diff --git a/lib/kokkos/generate_makefile.bash b/lib/kokkos/generate_makefile.bash index 25370daa3f..70dd61f9af 100755 --- a/lib/kokkos/generate_makefile.bash +++ b/lib/kokkos/generate_makefile.bash @@ -164,6 +164,7 @@ display_help_text() { echo " AMD_GFX942 = AMD GPU MI300 GFX942" echo " AMD_GFX1030 = AMD GPU V620/W6800 GFX1030" echo " AMD_GFX1100 = AMD GPU RX 7900 XT(X) GFX1100" + echo " AMD_GFX1103 = AMD APU Radeon 740M/760M/780M/880M/890M GFX1103" echo " [ARM]" echo " ARMV80 = ARMv8.0 Compatible CPU" echo " ARMV81 = ARMv8.1 Compatible CPU" diff --git a/python/lammps/pylammps.py b/python/lammps/pylammps.py index 96384255c2..49a2c6cb09 100644 --- a/python/lammps/pylammps.py +++ b/python/lammps/pylammps.py @@ -463,13 +463,19 @@ class PyLammps(object): self.lmp = lammps(name=name,cmdargs=cmdargs,ptr=ptr,comm=comm) else: self.lmp = lammps(name=name,cmdargs=cmdargs,ptr=None,comm=comm) - print("LAMMPS output is captured by PyLammps wrapper") + self.comm_nprocs = self.lmp.extract_setting("world_size") + self.comm_me = self.lmp.extract_setting("world_rank") + if self.comm_me == 0: + print("LAMMPS output is captured by PyLammps wrapper") + if self.comm_nprocs > 1: + print("WARNING: Using PyLammps with multiple MPI ranks is experimental. Not all functionality is supported.") self._cmd_history = [] self._enable_cmd_history = False self.runs = [] if not self.lmp.has_package("PYTHON"): - print("WARNING: run thermo data not captured since PYTHON LAMMPS package is not enabled") + if self.comm_me == 0: + print("WARNING: run thermo data not captured since PYTHON LAMMPS package is not enabled") def __enter__(self): return self @@ -727,7 +733,15 @@ class PyLammps(object): def eval(self, expr): """ - Evaluate expression + Evaluate LAMMPS input file expression. + + This is equivalent to using immediate variable expressions in the format "$(...)" + in the LAMMPS input and will return the result of that expression. + + .. warning:: + + This function is only supported on MPI rank 0. Calling it from a different + MPI rank will raise an exception. :param expr: the expression string that should be evaluated inside of LAMMPS :type expr: string @@ -735,6 +749,9 @@ class PyLammps(object): :return: the value of the evaluated expression :rtype: float if numeric, string otherwise """ + if self.comm_me > 0: + raise Exception("PyLammps.eval() may only be used on MPI rank 0") + value = self.lmp_print('"$(%s)"' % expr).strip() try: return float(value) diff --git a/src/AMOEBA/amoeba_convolution.h b/src/AMOEBA/amoeba_convolution.h index 60825bb8b6..c89d214cc5 100644 --- a/src/AMOEBA/amoeba_convolution.h +++ b/src/AMOEBA/amoeba_convolution.h @@ -14,8 +14,8 @@ #ifndef LMP_AMOEBA_CONVOLUTION_H #define LMP_AMOEBA_CONVOLUTION_H -#include "pointers.h" #include "lmpfftsettings.h" +#include "pointers.h" namespace LAMMPS_NS { @@ -27,7 +27,7 @@ class AmoebaConvolution : protected Pointers { int nxlo_in, nxhi_in, nylo_in, nyhi_in, nzlo_in, nzhi_in; int nxlo_out, nxhi_out, nylo_out, nyhi_out, nzlo_out, nzhi_out; int nxlo_fft, nxhi_fft, nylo_fft, nyhi_fft, nzlo_fft, nzhi_fft; - bigint nfft_global; // nx * ny * nz + bigint nfft_global; // nx * ny * nz FFT_SCALAR *grid_brick_start; // lower left corner of (c)grid_brick data AmoebaConvolution(class LAMMPS *, class Pair *, int, int, int, int, int); @@ -40,11 +40,11 @@ class AmoebaConvolution : protected Pointers { double time_fft; protected: - int which; // caller name for convolution being performed - int flag3d; // 1 if using 3d grid_brick, 0 for 4d cgrid_brick - int nbrick_owned; // owned grid points in brick decomp - int nbrick_ghosts; // owned + ghost brick grid points - int ngrid_either; // max of nbrick_owned or nfft_owned + int which; // caller name for convolution being performed + int flag3d; // 1 if using 3d grid_brick, 0 for 4d cgrid_brick + int nbrick_owned; // owned grid points in brick decomp + int nbrick_ghosts; // owned + ghost brick grid points + int ngrid_either; // max of nbrick_owned or nfft_owned class Pair *amoeba; class FFT3d *fft1, *fft2; diff --git a/src/DPD-REACT/fix_eos_table.cpp b/src/DPD-REACT/fix_eos_table.cpp index 42567119df..b287a6b0dd 100644 --- a/src/DPD-REACT/fix_eos_table.cpp +++ b/src/DPD-REACT/fix_eos_table.cpp @@ -307,8 +307,8 @@ void FixEOStable::param_extract(Table *tb, Table *tb2, char *line) while (word) { if (strcmp(word,"N") == 0) { word = strtok(nullptr," \t\n\r\f"); - tb->ninput = atoi(word); - tb2->ninput = atoi(word); + tb->ninput = std::stoi(word); + tb2->ninput = std::stoi(word); } else { error->one(FLERR,"Invalid keyword in fix eos/table parameters"); } diff --git a/src/DPD-REACT/fix_eos_table_rx.cpp b/src/DPD-REACT/fix_eos_table_rx.cpp index bf71b502f0..bd7f4bd8cc 100644 --- a/src/DPD-REACT/fix_eos_table_rx.cpp +++ b/src/DPD-REACT/fix_eos_table_rx.cpp @@ -134,11 +134,11 @@ FixEOStableRX::FixEOStableRX(LAMMPS *lmp, int narg, char **arg) : } if (rx_flag) read_file(arg[7]); - else dHf[0] = atof(arg[7]); + else dHf[0] = std::stod(arg[7]); if (narg==10) { - energyCorr[0] = atof(arg[8]); - tempCorrCoeff[0] = atof(arg[9]); + energyCorr[0] = std::stod(arg[8]); + tempCorrCoeff[0] = std::stod(arg[9]); } comm_forward = 3; @@ -373,11 +373,11 @@ void FixEOStableRX::read_file(char *file) if (strcmp(words[0],&atom->dvname[ispecies][0]) == 0) break; if (ispecies < nspecies) { - dHf[ispecies] = atof(words[1]); + dHf[ispecies] = std::stod(words[1]); if (nwords > min_params_per_line+1) { - energyCorr[ispecies] = atof(words[2]); - tempCorrCoeff[ispecies] = atof(words[3]); - moleculeCorrCoeff[ispecies] = atof(words[4]); + energyCorr[ispecies] = std::stod(words[2]); + tempCorrCoeff[ispecies] = std::stod(words[3]); + moleculeCorrCoeff[ispecies] = std::stod(words[4]); } } } @@ -482,7 +482,7 @@ void FixEOStableRX::read_table(Table *tb, Table *tb2, char *file, char *keyword) word = strtok(line," \t\n\r\f"); word = strtok(nullptr," \t\n\r\f"); - rtmp = atof(word); + rtmp = std::stod(word); for (int icolumn=0; icolumn < ncolumn; icolumn++) { ispecies = eosSpecies[icolumn]; @@ -491,7 +491,7 @@ void FixEOStableRX::read_table(Table *tb, Table *tb2, char *file, char *keyword) Table *tbl2 = &tables2[ispecies]; word = strtok(nullptr," \t\n\r\f"); - tmpE = atof(word); + tmpE = std::stod(word); tbl->rfile[i] = rtmp; tbl->efile[i] = tmpE; @@ -574,7 +574,7 @@ void FixEOStableRX::param_extract(Table *tb, char *line) char *word = strtok(line," \t\n\r\f"); if (strcmp(word,"N") == 0) { word = strtok(nullptr," \t\n\r\f"); - tb->ninput = atoi(word); + tb->ninput = std::stoi(word); } else error->one(FLERR,"Invalid keyword in fix eos/table/rx parameters"); word = strtok(nullptr," \t\n\r\f"); diff --git a/src/DPD-REACT/fix_rx.cpp b/src/DPD-REACT/fix_rx.cpp index a7e9e4ea77..fe7538bd10 100644 --- a/src/DPD-REACT/fix_rx.cpp +++ b/src/DPD-REACT/fix_rx.cpp @@ -863,7 +863,7 @@ void FixRX::read_file(char *file) word = strtok(line," \t\n\r\f"); while (word != nullptr) { - tmpStoich = atof(word); + tmpStoich = std::stod(word); word = strtok(nullptr, " \t\n\r\f"); for (ispecies = 0; ispecies < nspecies; ispecies++) { if (strcmp(word,&atom->dvname[ispecies][0]) == 0) { @@ -886,13 +886,13 @@ void FixRX::read_file(char *file) if (strcmp(word,"=") == 0) sign = 1.0; if (strcmp(word,"+") != 0 && strcmp(word,"=") != 0) { if (word==nullptr) error->all(FLERR,"Missing parameters in reaction kinetic equation"); - Arr[nreactions] = atof(word); + Arr[nreactions] = std::stod(word); word = strtok(nullptr, " \t\n\r\f"); if (word==nullptr) error->all(FLERR,"Missing parameters in reaction kinetic equation"); - nArr[nreactions] = atof(word); + nArr[nreactions] = std::stod(word); word = strtok(nullptr, " \t\n\r\f"); if (word==nullptr) error->all(FLERR,"Missing parameters in reaction kinetic equation"); - Ea[nreactions] = atof(word); + Ea[nreactions] = std::stod(word); sign = -1.0; break; } diff --git a/src/DPD-REACT/pair_exp6_rx.cpp b/src/DPD-REACT/pair_exp6_rx.cpp index e0ac9c0b27..0cd9dd3a05 100644 --- a/src/DPD-REACT/pair_exp6_rx.cpp +++ b/src/DPD-REACT/pair_exp6_rx.cpp @@ -889,15 +889,15 @@ void PairExp6rx::read_file2(char *file) if (strcmp(words[0],"alpha") == 0) { for (int ii=1; iininput = atoi(word); + tb->ninput = std::stoi(word); } else if (strcmp(word,"R") == 0 || strcmp(word,"RSQ") == 0) { if (strcmp(word,"R") == 0) tb->rflag = RLINEAR; else if (strcmp(word,"RSQ") == 0) tb->rflag = RSQ; word = strtok(nullptr," \t\n\r\f"); - tb->rlo = atof(word); + tb->rlo = std::stod(word); word = strtok(nullptr," \t\n\r\f"); - tb->rhi = atof(word); + tb->rhi = std::stod(word); } else if (strcmp(word,"FP") == 0) { tb->fpflag = 1; word = strtok(nullptr," \t\n\r\f"); - tb->fplo = atof(word); + tb->fplo = std::stod(word); word = strtok(nullptr," \t\n\r\f"); - tb->fphi = atof(word); + tb->fphi = std::stod(word); } else { printf("WORD: %s\n",word); error->one(FLERR,"Invalid keyword in pair table parameters"); diff --git a/src/DPD-REACT/pair_multi_lucy_rx.cpp b/src/DPD-REACT/pair_multi_lucy_rx.cpp index c248d92694..178b5c3b85 100644 --- a/src/DPD-REACT/pair_multi_lucy_rx.cpp +++ b/src/DPD-REACT/pair_multi_lucy_rx.cpp @@ -612,20 +612,20 @@ void PairMultiLucyRX::param_extract(Table *tb, char *line) while (word) { if (strcmp(word,"N") == 0) { word = strtok(nullptr," \t\n\r\f"); - tb->ninput = atoi(word); + tb->ninput = std::stoi(word); } else if (strcmp(word,"R") == 0 || strcmp(word,"RSQ") == 0) { if (strcmp(word,"R") == 0) tb->rflag = RLINEAR; else if (strcmp(word,"RSQ") == 0) tb->rflag = RSQ; word = strtok(nullptr," \t\n\r\f"); - tb->rlo = atof(word); + tb->rlo = std::stod(word); word = strtok(nullptr," \t\n\r\f"); - tb->rhi = atof(word); + tb->rhi = std::stod(word); } else if (strcmp(word,"FP") == 0) { tb->fpflag = 1; word = strtok(nullptr," \t\n\r\f"); - tb->fplo = atof(word); + tb->fplo = std::stod(word); word = strtok(nullptr," \t\n\r\f"); - tb->fphi = atof(word); + tb->fphi = std::stod(word); } else { printf("WORD: %s\n",word); error->one(FLERR,"Invalid keyword in pair table parameters"); diff --git a/src/Depend.sh b/src/Depend.sh index e55b100975..85542b21c0 100755 --- a/src/Depend.sh +++ b/src/Depend.sh @@ -195,6 +195,10 @@ if (test $1 = "ML-SNAP") then depend INTEL fi +if (test $1 = "ML-UF3") then + depend KOKKOS +fi + if (test $1 = "CG-SPICA") then depend GPU depend KOKKOS diff --git a/src/ELECTRODE/electrode_vector.cpp b/src/ELECTRODE/electrode_vector.cpp index fc2cca5e46..1ffd8c1317 100644 --- a/src/ELECTRODE/electrode_vector.cpp +++ b/src/ELECTRODE/electrode_vector.cpp @@ -37,8 +37,7 @@ using namespace LAMMPS_NS; using namespace MathConst; ElectrodeVector::ElectrodeVector(LAMMPS *lmp, int sensor_group, int source_group, double eta, - bool invert_source) : - Pointers(lmp) + bool invert_source) : Pointers(lmp) { igroup = sensor_group; // group of all atoms at which we calculate potential this->source_group = source_group; // group of all atoms influencing potential diff --git a/src/ELECTRODE/fix_electrode_conp.cpp b/src/ELECTRODE/fix_electrode_conp.cpp index 5e92af7a82..565a750af5 100644 --- a/src/ELECTRODE/fix_electrode_conp.cpp +++ b/src/ELECTRODE/fix_electrode_conp.cpp @@ -931,7 +931,8 @@ void FixElectrodeConp::update_charges() dot_old = dot_new; } recompute_potential(std::move(b), q_local); - if ((delta > cg_threshold) && (comm->me == 0)) error->warning(FLERR, "CG threshold not reached"); + if ((delta > cg_threshold) && (comm->me == 0)) + error->warning(FLERR, "CG threshold not reached"); } else { error->all(FLERR, "This algorithm is not implemented, yet"); } diff --git a/src/H5MD/dump_h5md.cpp b/src/H5MD/dump_h5md.cpp index 21c02bfc23..f8b7ecfbd4 100644 --- a/src/H5MD/dump_h5md.cpp +++ b/src/H5MD/dump_h5md.cpp @@ -47,7 +47,7 @@ static int element_args(int narg, char **arg, int *every) while (iarg narg) error->all(FLERR,"Illegal package intel command"); - _offload_tpc = atoi(arg[iarg+1]); + _offload_tpc = std::stoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"tptask") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal package intel command"); - _offload_threads = atoi(arg[iarg+1]); + _offload_threads = std::stoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"no_affinity") == 0) { no_affinity = 1; @@ -150,7 +150,7 @@ FixIntel::FixIntel(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg) iarg++; } else if (strcmp(arg[iarg],"buffers") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal package intel command"); - _allow_separate_buffers = atoi(arg[iarg+1]); + _allow_separate_buffers = std::stoi(arg[iarg+1]); iarg += 2; } else error->all(FLERR,"Illegal package intel command"); } @@ -1198,7 +1198,7 @@ int FixIntel::set_host_affinity(const int nomp) if (p == nullptr) return -1; ncores = 0; while (fgets(readbuf, 512, p)) { - proc_list[ncores] = atoi(readbuf); + proc_list[ncores] = std::stoi(readbuf); ncores++; } pclose(p); @@ -1218,7 +1218,7 @@ int FixIntel::set_host_affinity(const int nomp) if (nthreads == 0) { estring = getenv("OMP_NUM_THREADS"); if (estring != nullptr) { - nthreads = atoi(estring); + nthreads = std::stoi(estring); if (nthreads < 2) nthreads = 1; } else nthreads = 1; @@ -1251,7 +1251,7 @@ int FixIntel::set_host_affinity(const int nomp) if (p == nullptr) return -1; while (fgets(readbuf, 512, p)) { - lwp = atoi(readbuf); + lwp = std::stoi(readbuf); int first = coi_cores + node_rank * mpi_cores; CPU_ZERO(&cpuset); for (int i = first; i < first + mpi_cores; i++) @@ -1287,7 +1287,7 @@ int FixIntel::set_host_affinity(const int nomp) if (p == nullptr) return -1; while (fgets(readbuf, 512, p)) { - lwp = atoi(readbuf); + lwp = std::stoi(readbuf); nlwp++; if (nlwp <= plwp) continue; diff --git a/src/KIM/kim_param.cpp b/src/KIM/kim_param.cpp index c50474fe67..7f2c891507 100644 --- a/src/KIM/kim_param.cpp +++ b/src/KIM/kim_param.cpp @@ -204,14 +204,14 @@ void KimParam::command(int narg, char **arg) if (npos != std::string::npos) { argtostr[npos] = ' '; auto words = utils::split_words(argtostr); - nlbound = atoi(words[0].c_str()); - nubound = atoi(words[1].c_str()); + nlbound = std::stoi(words[0]); + nubound = std::stoi(words[1]); if ((nubound < 1) || (nubound > extent) || (nlbound < 1) || (nlbound > nubound)) error->all(FLERR, "Illegal index_range '{}-{}' for '{}' parameter with the " "extent of '{}'",nlbound, nubound, paramname, extent); } else { - nlbound = atoi(argtostr.c_str()); + nlbound = std::stoi(argtostr); if (nlbound < 1 || nlbound > extent) error->all(FLERR, "Illegal index '{}' for '{}' parameter with the extent of '{}'", diff --git a/src/KIM/pair_kim.cpp b/src/KIM/pair_kim.cpp index e0fc90ea46..ded7221d76 100644 --- a/src/KIM/pair_kim.cpp +++ b/src/KIM/pair_kim.cpp @@ -491,14 +491,14 @@ void PairKIM::coeff(int narg, char **arg) if (npos != std::string::npos) { argtostr[npos] = ' '; auto words = utils::split_words(argtostr); - nlbound = atoi(words[0].c_str()); - nubound = atoi(words[1].c_str()); + nlbound = std::stoi(words[0]); + nubound = std::stoi(words[1]); if ((nubound < 1) || (nubound > extent) || (nlbound < 1) || (nlbound > nubound)) error->all(FLERR,"Illegal index_range '{}-{}' for '{}' parameter with the extent " "of '{}'", nlbound, nubound, paramname, extent); } else { - nlbound = atoi(argtostr.c_str()); + nlbound = std::stoi(argtostr); if ((nlbound < 1) || (nlbound > extent)) error->all(FLERR,"Illegal index '{}' for '{}' parameter with the extent of '{}'", diff --git a/src/KOKKOS/Install.sh b/src/KOKKOS/Install.sh index 94effc5d68..31359d4e4a 100755 --- a/src/KOKKOS/Install.sh +++ b/src/KOKKOS/Install.sh @@ -391,6 +391,8 @@ action pair_tersoff_mod_kokkos.cpp pair_tersoff_mod.cpp action pair_tersoff_mod_kokkos.h pair_tersoff_mod.h action pair_tersoff_zbl_kokkos.cpp pair_tersoff_zbl.cpp action pair_tersoff_zbl_kokkos.h pair_tersoff_zbl.h +action pair_uf3_kokkos.cpp pair_uf3.cpp +action pair_uf3_kokkos.h pair_uf3.h action pair_vashishta_kokkos.cpp pair_vashishta.cpp action pair_vashishta_kokkos.h pair_vashishta.h action pair_yukawa_kokkos.cpp diff --git a/src/KOKKOS/dihedral_charmmfsw_kokkos.cpp b/src/KOKKOS/dihedral_charmmfsw_kokkos.cpp index 3414a02ec4..a3c000b3cf 100644 --- a/src/KOKKOS/dihedral_charmmfsw_kokkos.cpp +++ b/src/KOKKOS/dihedral_charmmfsw_kokkos.cpp @@ -380,16 +380,17 @@ void DihedralCharmmfswKokkos::operator()(TagDihedralCharmmfswCompute const F_FLOAT dely = x(i1,1) - x(i4,1); const F_FLOAT delz = x(i1,2) - x(i4,2); const F_FLOAT rsq = delx*delx + dely*dely + delz*delz; + const F_FLOAT r = sqrt(rsq); const F_FLOAT r2inv = 1.0/rsq; const F_FLOAT r6inv = r2inv*r2inv*r2inv; F_FLOAT forcecoul; if (implicit) forcecoul = qqrd2e * q[i1]*q[i4]*r2inv; - else forcecoul = qqrd2e * q[i1]*q[i4]*sqrt(r2inv); + else if (dihedflag) forcecoul = qqrd2e * q[i1]*q[i4]*sqrt(r2inv); + else forcecoul = qqrd2e * q[i1]*q[i4]*(sqrt(r2inv) - r*cut_coulinv14*cut_coulinv14); const F_FLOAT forcelj = r6inv * (d_lj14_1(itype,jtype)*r6inv - d_lj14_2(itype,jtype)); const F_FLOAT fpair = d_weight[type] * (forcelj+forcecoul)*r2inv; - const F_FLOAT r = sqrt(rsq); F_FLOAT ecoul = 0.0; F_FLOAT evdwl = 0.0; F_FLOAT evdwl14_12, evdwl14_6; diff --git a/src/KOKKOS/dihedral_hybrid_kokkos.h b/src/KOKKOS/dihedral_hybrid_kokkos.h index 29a3d29689..63a59505af 100644 --- a/src/KOKKOS/dihedral_hybrid_kokkos.h +++ b/src/KOKKOS/dihedral_hybrid_kokkos.h @@ -11,7 +11,7 @@ See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ -#ifdef BOND_CLASS +#ifdef DIHEDRAL_CLASS // clang-format off DihedralStyle(hybrid/kk,DihedralHybridKokkos); DihedralStyle(hybrid/kk/device,DihedralHybridKokkos); @@ -20,8 +20,8 @@ DihedralStyle(hybrid/kk/host,DihedralHybridKokkos); #else // clang-format off -#ifndef LMP_BOND_HYBRID_KOKKOS_H -#define LMP_BOND_HYBRID_KOKKOS_H +#ifndef LMP_DIHEDRAL_HYBRID_KOKKOS_H +#define LMP_DIHEDRAL_HYBRID_KOKKOS_H #include "dihedral_hybrid.h" #include "kokkos_type.h" @@ -40,13 +40,13 @@ class DihedralHybridKokkos : public DihedralHybrid { double memory_usage() override; private: - int maxbond_all; + int maxdihedral_all; class NeighborKokkos *neighborKK; - DAT::tdual_int_1d k_map; // which style each bond type points to - DAT::tdual_int_1d k_nbondlist; // # of bonds in sub-style bondlists - DAT::tdual_int_3d k_bondlist; // bondlist for each sub-style + DAT::tdual_int_1d k_map; // which style each dihedral type points to + DAT::tdual_int_1d k_ndihedrallist; // # of dihedrals in sub-style dihedrallists + DAT::tdual_int_3d k_dihedrallist; // dihedrallist for each sub-style void allocate() override; void deallocate() override; diff --git a/src/KOKKOS/improper_hybrid_kokkos.h b/src/KOKKOS/improper_hybrid_kokkos.h index f2a80f6a0c..78bafe5df6 100644 --- a/src/KOKKOS/improper_hybrid_kokkos.h +++ b/src/KOKKOS/improper_hybrid_kokkos.h @@ -11,7 +11,7 @@ See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ -#ifdef BOND_CLASS +#ifdef IMPROPER_CLASS // clang-format off ImproperStyle(hybrid/kk,ImproperHybridKokkos); ImproperStyle(hybrid/kk/device,ImproperHybridKokkos); @@ -20,8 +20,8 @@ ImproperStyle(hybrid/kk/host,ImproperHybridKokkos); #else // clang-format off -#ifndef LMP_BOND_HYBRID_KOKKOS_H -#define LMP_BOND_HYBRID_KOKKOS_H +#ifndef LMP_IMPROPER_HYBRID_KOKKOS_H +#define LMP_IMPROPER_HYBRID_KOKKOS_H #include "improper_hybrid.h" #include "kokkos_type.h" diff --git a/src/KOKKOS/kokkos.cpp b/src/KOKKOS/kokkos.cpp index 58b9436af6..f4ba967c9a 100644 --- a/src/KOKKOS/kokkos.cpp +++ b/src/KOKKOS/kokkos.cpp @@ -104,14 +104,14 @@ KokkosLMP::KokkosLMP(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) int set_flag = 0; char *str; if ((str = getenv("SLURM_LOCALID"))) { - int local_rank = atoi(str); + int local_rank = std::stoi(str); device = local_rank % ngpus; if (device >= skip_gpu) device++; set_flag = 1; } if ((str = getenv("FLUX_TASK_LOCAL_ID"))) { if (ngpus > 0) { - int local_rank = atoi(str); + int local_rank = std::stoi(str); device = local_rank % ngpus; if (device >= skip_gpu) device++; set_flag = 1; @@ -119,7 +119,7 @@ KokkosLMP::KokkosLMP(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) } if ((str = getenv("MPT_LRANK"))) { if (ngpus > 0) { - int local_rank = atoi(str); + int local_rank = std::stoi(str); device = local_rank % ngpus; if (device >= skip_gpu) device++; set_flag = 1; @@ -127,7 +127,7 @@ KokkosLMP::KokkosLMP(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) } if ((str = getenv("MV2_COMM_WORLD_LOCAL_RANK"))) { if (ngpus > 0) { - int local_rank = atoi(str); + int local_rank = std::stoi(str); device = local_rank % ngpus; if (device >= skip_gpu) device++; set_flag = 1; @@ -135,7 +135,7 @@ KokkosLMP::KokkosLMP(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) } if ((str = getenv("OMPI_COMM_WORLD_LOCAL_RANK"))) { if (ngpus > 0) { - int local_rank = atoi(str); + int local_rank = std::stoi(str); device = local_rank % ngpus; if (device >= skip_gpu) device++; set_flag = 1; @@ -143,7 +143,7 @@ KokkosLMP::KokkosLMP(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) } if ((str = getenv("PMI_LOCAL_RANK"))) { if (ngpus > 0) { - int local_rank = atoi(str); + int local_rank = std::stoi(str); device = local_rank % ngpus; if (device >= skip_gpu) device++; set_flag = 1; diff --git a/src/KOKKOS/pair_exp6_rx_kokkos.cpp b/src/KOKKOS/pair_exp6_rx_kokkos.cpp index dad7413669..7eddd7527e 100644 --- a/src/KOKKOS/pair_exp6_rx_kokkos.cpp +++ b/src/KOKKOS/pair_exp6_rx_kokkos.cpp @@ -1771,9 +1771,9 @@ void PairExp6rxKokkos::read_file(char *file) params[nparams].potential = utils::strdup(words[1]); if (strcmp(params[nparams].potential,"exp6") == 0) { - params[nparams].alpha = atof(words[2]); - params[nparams].epsilon = atof(words[3]); - params[nparams].rm = atof(words[4]); + params[nparams].alpha = std::stod(words[2]); + params[nparams].epsilon = std::stod(words[3]); + params[nparams].rm = std::stod(words[4]); if (params[nparams].epsilon <= 0.0 || params[nparams].rm <= 0.0 || params[nparams].alpha < 0.0) error->all(FLERR,"Illegal exp6/rx parameters. Rm and Epsilon must be greater than zero. Alpha cannot be negative."); diff --git a/src/KOKKOS/pair_uf3_kokkos.cpp b/src/KOKKOS/pair_uf3_kokkos.cpp index 4a3bcc164b..59112ddab0 100644 --- a/src/KOKKOS/pair_uf3_kokkos.cpp +++ b/src/KOKKOS/pair_uf3_kokkos.cpp @@ -1655,7 +1655,7 @@ double PairUF3Kokkos::single(int /*i*/, int /*j*/, int itype, int jt namespace LAMMPS_NS { template class PairUF3Kokkos; -#ifdef KOKKOS_ENABLE_CUDA +#ifdef KOKKOS_ENABLE_GPU template class PairUF3Kokkos; #endif } // namespace LAMMPS_NS diff --git a/src/KOKKOS/pair_uf3_kokkos.h b/src/KOKKOS/pair_uf3_kokkos.h index 15c2832da1..952e2aba25 100644 --- a/src/KOKKOS/pair_uf3_kokkos.h +++ b/src/KOKKOS/pair_uf3_kokkos.h @@ -22,6 +22,7 @@ // clang-format off PairStyle(uf3/kk,PairUF3Kokkos) PairStyle(uf3/kk/device,PairUF3Kokkos) +PairStyle(uf3/kk/host,PairUF3Kokkos) // clang-format on #else @@ -44,7 +45,7 @@ template class PairUF3Kokkos : public PairUF3 { void compute(int, int) override; void settings(int, char **) override; void coeff(int, char **) override; - void allocate(); + void allocate() override; void init_style() override; void init_list(int, class NeighList *) override; // needed for ptr to full neigh list double init_one(int, int) override; // needed for cutoff radius for neighbour list @@ -117,9 +118,11 @@ template class PairUF3Kokkos : public PairUF3 { std::vector get_dncoefficients(const double *knots, const double coefficient) const; template + KOKKOS_INLINE_FUNCTION void twobody(const int itype, const int jtype, const F_FLOAT r, F_FLOAT &evdwl, F_FLOAT &fpair) const; template + KOKKOS_INLINE_FUNCTION void threebody(const int itype, const int jtype, const int ktype, const F_FLOAT value_rij, const F_FLOAT value_rik, const F_FLOAT value_rjk, F_FLOAT &evdwl3, F_FLOAT (&fforce)[3]) const; diff --git a/src/MANYBODY/pair_extep.cpp b/src/MANYBODY/pair_extep.cpp index 7ed65f0f71..607106b4ce 100644 --- a/src/MANYBODY/pair_extep.cpp +++ b/src/MANYBODY/pair_extep.cpp @@ -640,7 +640,7 @@ void PairExTeP::read_file(char *file) if (!utils::is_integer(kname)) continue; - int Ni = atoi(kname.c_str()); + int Ni = std::stoi(kname); int Nj = values.next_int(); double spline_val = values.next_double(); double spline_derx = values.next_double(); diff --git a/src/MGPT/mgpt_readpot.cpp b/src/MGPT/mgpt_readpot.cpp index 3cc5d5d693..3bb39ab353 100644 --- a/src/MGPT/mgpt_readpot.cpp +++ b/src/MGPT/mgpt_readpot.cpp @@ -531,7 +531,7 @@ void potdata::readpot(const char *parmin_file,const char *potin_file,const doubl /* int main(int argc,char *argv[]) { - double vol = atof(argv[3]); + double vol = std::stod(argv[3]); int n = 25,i; printf("%% parmin = %s\n%% potin = %s\n%% vol = %15.5e\n", diff --git a/src/ML-IAP/mliap_unified.cpp b/src/ML-IAP/mliap_unified.cpp index 3f0d0b2a8f..9998ff2af4 100644 --- a/src/ML-IAP/mliap_unified.cpp +++ b/src/ML-IAP/mliap_unified.cpp @@ -253,8 +253,10 @@ void LAMMPS_NS::update_pair_energy(MLIAPData *data, double *eij) double e = 0.5 * eij[ii]; // must not count any contribution where i is not a local atom - data->eatoms[i] += e; - e_total += e; + if (i < data->nlocal) { + data->eatoms[i] += e; + e_total += e; + } } data->energy = e_total; } diff --git a/src/ML-UF3/pair_uf3.cpp b/src/ML-UF3/pair_uf3.cpp index a952403287..1f91f25f5b 100644 --- a/src/ML-UF3/pair_uf3.cpp +++ b/src/ML-UF3/pair_uf3.cpp @@ -120,7 +120,6 @@ void PairUF3::settings(int narg, char **arg) "Invalid number of arguments for pair_style uf3" " Are you using a 2-body or 2 & 3-body UF potential?"); nbody_flag = utils::numeric(FLERR, arg[0], true, lmp); - const int num_of_elements = atom->ntypes; if (nbody_flag == 2) { pot_3b = false; manybody_flag = 0; diff --git a/src/MOLFILE/reader_molfile.cpp b/src/MOLFILE/reader_molfile.cpp index 43154d658b..81f2fa3d10 100644 --- a/src/MOLFILE/reader_molfile.cpp +++ b/src/MOLFILE/reader_molfile.cpp @@ -322,7 +322,7 @@ void ReaderMolfile::read_atoms(int n, int nfield, double **fields) ++nid; if (mf->property(MFI::P_TYPE,nid-1,buf) != MFI::P_NONE) { - mytype = atoi(buf); + mytype = std::stoi(buf); } else mytype = 0; for (m = 0; m < nfield; m++) { diff --git a/src/OPENMP/reaxff_forces_omp.cpp b/src/OPENMP/reaxff_forces_omp.cpp index 26922add1e..ce8ad0716f 100644 --- a/src/OPENMP/reaxff_forces_omp.cpp +++ b/src/OPENMP/reaxff_forces_omp.cpp @@ -155,7 +155,7 @@ namespace ReaxFF { /* ---------------------------------------------------------------------- */ static void Validate_ListsOMP(reax_system *system, reax_list **lists, - int step, int n, int N, int numH) + int step, int N, int numH) { int comp, Hindex; reax_list *bonds, *hbonds; @@ -195,7 +195,7 @@ namespace ReaxFF { #if defined(_OPENMP) #pragma omp for schedule(guided) #endif - for (int i = 0; i < n; ++i) { + for (int i = 0; i < N; ++i) { Hindex = system->my_atoms[i].Hindex; if (Hindex > -1) { system->my_atoms[i].num_hbonds = @@ -457,8 +457,7 @@ namespace ReaxFF { workspace->realloc.num_bonds = num_bonds; workspace->realloc.num_hbonds = num_hbonds; - Validate_ListsOMP(system, lists, data->step, - system->n, system->N, system->numH); + Validate_ListsOMP(system, lists, data->step, system->N, system->numH); } /* ---------------------------------------------------------------------- */ diff --git a/src/PYTHON/python_compat.h b/src/PYTHON/python_compat.h index 775435ffd8..d1194056d8 100644 --- a/src/PYTHON/python_compat.h +++ b/src/PYTHON/python_compat.h @@ -20,13 +20,12 @@ #if PY_MAJOR_VERSION == 2 #if defined(_MSC_VER) || defined(__MINGW32__) #define PY_INT_FROM_LONG(X) PyLong_FromLongLong(X) +#define PY_INT_AS_LONG(X) PyLong_AsLongLong(X) +#define PY_LONG_FROM_STRING(X) std::stoll(X) #else #define PY_INT_FROM_LONG(X) PyInt_FromLong(X) -#endif -#if defined(_MSC_VER) || defined(__MINGW32__) -#define PY_INT_AS_LONG(X) PyLong_AsLongLong(X) -#else #define PY_INT_AS_LONG(X) PyInt_AsLong(X) +#define PY_LONG_FROM_STRING(X) std::stol(X) #endif #define PY_STRING_FROM_STRING(X) PyString_FromString(X) #define PY_VOID_POINTER(X) PyCObject_FromVoidPtr((void *) X, nullptr) @@ -35,13 +34,12 @@ #elif PY_MAJOR_VERSION == 3 #if defined(_MSC_VER) || defined(__MINGW32__) #define PY_INT_FROM_LONG(X) PyLong_FromLongLong(X) +#define PY_INT_AS_LONG(X) PyLong_AsLongLong(X) +#define PY_LONG_FROM_STRING(X) std::stoll(X) #else #define PY_INT_FROM_LONG(X) PyLong_FromLong(X) -#endif -#if defined(_MSC_VER) || defined(__MINGW32__) -#define PY_INT_AS_LONG(X) PyLong_AsLongLong(X) -#else #define PY_INT_AS_LONG(X) PyLong_AsLong(X) +#define PY_LONG_FROM_STRING(X) std::stol(X) #endif #define PY_STRING_FROM_STRING(X) PyUnicode_FromString(X) #define PY_VOID_POINTER(X) PyCapsule_New((void *) X, nullptr, nullptr) diff --git a/src/PYTHON/python_impl.cpp b/src/PYTHON/python_impl.cpp index 87a57187bf..91e0a7abe5 100644 --- a/src/PYTHON/python_impl.cpp +++ b/src/PYTHON/python_impl.cpp @@ -341,7 +341,7 @@ void PythonImpl::invoke_function(int ifunc, char *result) if (!str) error->all(FLERR, "Could not evaluate Python function {} input variable: {}", pfuncs[ifunc].name, pfuncs[ifunc].svalue[i]); - pValue = PY_INT_FROM_LONG(atoi(str)); + pValue = PY_INT_FROM_LONG(PY_LONG_FROM_STRING(str)); } else { pValue = PY_INT_FROM_LONG(pfuncs[ifunc].ivalue[i]); } @@ -351,7 +351,7 @@ void PythonImpl::invoke_function(int ifunc, char *result) if (!str) error->all(FLERR, "Could not evaluate Python function {} input variable: {}", pfuncs[ifunc].name, pfuncs[ifunc].svalue[i]); - pValue = PyFloat_FromDouble(atof(str)); + pValue = PyFloat_FromDouble(std::stod(str)); } else { pValue = PyFloat_FromDouble(pfuncs[ifunc].dvalue[i]); } diff --git a/src/RHEO/bond_rheo_shell.cpp b/src/RHEO/bond_rheo_shell.cpp index 258d047086..dbe65e2b2f 100644 --- a/src/RHEO/bond_rheo_shell.cpp +++ b/src/RHEO/bond_rheo_shell.cpp @@ -45,7 +45,8 @@ using namespace RHEO_NS; /* ---------------------------------------------------------------------- */ BondRHEOShell::BondRHEOShell(LAMMPS *_lmp) : - BondBPM(_lmp), compute_surface(nullptr), k(nullptr), ecrit(nullptr), gamma(nullptr) + BondBPM(_lmp), k(nullptr), ecrit(nullptr), gamma(nullptr), dbond(nullptr), nbond(nullptr), + id_fix(nullptr), compute_surface(nullptr) { partial_flag = 1; comm_reverse = 1; diff --git a/src/RHEO/compute_rheo_grad.cpp b/src/RHEO/compute_rheo_grad.cpp index cdf90e1dc5..fdddea18b5 100644 --- a/src/RHEO/compute_rheo_grad.cpp +++ b/src/RHEO/compute_rheo_grad.cpp @@ -42,8 +42,9 @@ enum{COMMGRAD, COMMFIELD}; /* ---------------------------------------------------------------------- */ ComputeRHEOGrad::ComputeRHEOGrad(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg), fix_rheo(nullptr), list(nullptr), rho0(nullptr), compute_interface(nullptr), compute_kernel(nullptr), - gradv(nullptr), gradr(nullptr), grade(nullptr), gradn(nullptr) + Compute(lmp, narg, arg), gradv(nullptr), gradr(nullptr), grade(nullptr), gradn(nullptr), + fix_rheo(nullptr), rho0(nullptr), compute_kernel(nullptr), compute_interface(nullptr), + list(nullptr) { if (narg < 4) error->all(FLERR,"Illegal compute rheo/grad command"); @@ -130,10 +131,9 @@ void ComputeRHEOGrad::compute_peratom() { int i, j, k, ii, jj, jnum, itype, jtype, a, b, fluidi, fluidj; double xtmp, ytmp, ztmp, delx, dely, delz; - double rsq, imass, jmass; - double rhoi, rhoj, Voli, Volj, drho, de, deta; + double rsq, rhoi, rhoj, Voli, Volj, drho, de, deta; double vi[3], vj[3], vij[3]; - double wp, *dWij, *dWji; + double *dWij, *dWji; int inum, *ilist, *numneigh, **firstneigh; int *jlist; @@ -147,6 +147,7 @@ void ComputeRHEOGrad::compute_peratom() int *status = atom->rheo_status; int *type = atom->type; double *mass = atom->mass; + double *rmass = atom->rmass; int newton = force->newton; int dim = domain->dimension; @@ -225,8 +226,13 @@ void ComputeRHEOGrad::compute_peratom() } } - Voli = mass[itype] / rhoi; - Volj = mass[jtype] / rhoj; + if (rmass) { + Voli = rmass[i] / rhoi; + Volj = rmass[j] / rhoj; + } else { + Voli = mass[itype] / rhoi; + Volj = mass[jtype] / rhoj; + } vij[0] = vi[0] - vj[0]; vij[1] = vi[1] - vj[1]; @@ -236,7 +242,7 @@ void ComputeRHEOGrad::compute_peratom() if (energy_flag) de = energy[i] - energy[j]; if (eta_flag) deta = viscosity[i] - viscosity[j]; - wp = compute_kernel->calc_dw(i, j, delx, dely, delz, sqrt(rsq)); + compute_kernel->calc_dw(i, j, delx, dely, delz, sqrt(rsq)); dWij = compute_kernel->dWij; dWji = compute_kernel->dWji; diff --git a/src/RHEO/compute_rheo_interface.cpp b/src/RHEO/compute_rheo_interface.cpp index 8ccd4e6a3b..2d2a368926 100644 --- a/src/RHEO/compute_rheo_interface.cpp +++ b/src/RHEO/compute_rheo_interface.cpp @@ -45,8 +45,10 @@ static constexpr double EPSILON = 1e-1; /* ---------------------------------------------------------------------- */ ComputeRHEOInterface::ComputeRHEOInterface(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg), fix_rheo(nullptr), compute_kernel(nullptr), fp_store(nullptr), - rho0(nullptr), norm(nullptr), normwf(nullptr), chi(nullptr), id_fix_pa(nullptr) + Compute(lmp, narg, arg), chi(nullptr), fp_store(nullptr), fix_rheo(nullptr), + rho0(nullptr), norm(nullptr), normwf(nullptr), id_fix_pa(nullptr), list(nullptr), + compute_kernel(nullptr), fix_pressure(nullptr) + { if (narg != 3) error->all(FLERR,"Illegal compute rheo/interface command"); @@ -228,12 +230,11 @@ void ComputeRHEOInterface::compute_peratom() int ComputeRHEOInterface::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i, j, k, m; - m = 0; + int m = 0; double *rho = atom->rho; - for (i = 0; i < n; i++) { - j = list[i]; + for (int i = 0; i < n; i++) { + int j = list[i]; if (comm_stage == 0) { buf[m++] = fp_store[j][0]; buf[m++] = fp_store[j][1]; @@ -250,11 +251,10 @@ int ComputeRHEOInterface::pack_forward_comm(int n, int *list, double *buf, int / void ComputeRHEOInterface::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last; double *rho = atom->rho; - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { if (comm_stage == 0) { fp_store[i][0] = buf[m++]; fp_store[i][1] = buf[m++]; @@ -270,12 +270,10 @@ void ComputeRHEOInterface::unpack_forward_comm(int n, int first, double *buf) int ComputeRHEOInterface::pack_reverse_comm(int n, int first, double *buf) { - int i, k, m, last; double *rho = atom->rho; - - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { buf[m++] = norm[i]; buf[m++] = chi[i]; buf[m++] = normwf[i]; @@ -288,12 +286,11 @@ int ComputeRHEOInterface::pack_reverse_comm(int n, int first, double *buf) void ComputeRHEOInterface::unpack_reverse_comm(int n, int *list, double *buf) { - int i, k, j, m; double *rho = atom->rho; int *status = atom->rheo_status; - m = 0; - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; norm[j] += buf[m++]; chi[j] += buf[m++]; if (status[j] & PHASECHECK){ @@ -338,6 +335,7 @@ void ComputeRHEOInterface::store_forces() int *type = atom->type; int *mask = atom->mask; double *mass = atom->mass; + double *rmass = atom->rmass; double **f = atom->f; // When this is called, fp_store stores the pressure force @@ -349,7 +347,8 @@ void ComputeRHEOInterface::store_forces() if (fixlist.size() != 0) { for (const auto &fix : fixlist) { for (int i = 0; i < atom->nlocal; i++) { - minv = 1.0 / mass[type[i]]; + if (rmass) minv = 1.0 / rmass[i]; + else minv = 1.0 / mass[type[i]]; if (mask[i] & fix->groupbit) for (int a = 0; a < 3; a++) fp_store[i][a] = f[i][a] * minv; @@ -359,10 +358,18 @@ void ComputeRHEOInterface::store_forces() } } } else { - for (int i = 0; i < atom->nlocal; i++) { - minv = 1.0 / mass[type[i]]; - for (int a = 0; a < 3; a++) - fp_store[i][a] = (f[i][a] - fp_store[i][a]) * minv; + if (rmass) { + for (int i = 0; i < atom->nlocal; i++) { + minv = 1.0 / rmass[i]; + for (int a = 0; a < 3; a++) + fp_store[i][a] = (f[i][a] - fp_store[i][a]) * minv; + } + } else { + for (int i = 0; i < atom->nlocal; i++) { + minv = 1.0 / mass[type[i]]; + for (int a = 0; a < 3; a++) + fp_store[i][a] = (f[i][a] - fp_store[i][a]) * minv; + } } } diff --git a/src/RHEO/compute_rheo_kernel.cpp b/src/RHEO/compute_rheo_kernel.cpp index 4558ddccc8..108ec58b09 100644 --- a/src/RHEO/compute_rheo_kernel.cpp +++ b/src/RHEO/compute_rheo_kernel.cpp @@ -48,13 +48,11 @@ using namespace RHEO_NS; using namespace MathConst; using namespace MathExtra; -static constexpr int DELTA = 2000; - /* ---------------------------------------------------------------------- */ ComputeRHEOKernel::ComputeRHEOKernel(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg), - list(nullptr), C(nullptr), C0(nullptr), coordination(nullptr), compute_interface(nullptr) + Compute(lmp, narg, arg), coordination(nullptr), fix_rheo(nullptr), C(nullptr), C0(nullptr), + list(nullptr), compute_interface(nullptr) { if (narg != 4) error->all(FLERR,"Illegal compute rheo/kernel command"); @@ -190,7 +188,7 @@ double ComputeRHEOKernel::calc_w_self(int i, int j) double ComputeRHEOKernel::calc_w(int i, int j, double delx, double dely, double delz, double r) { - double w; + double w = 0.0; int corrections_i, corrections_j, corrections; if (kernel_style == WENDLANDC4) @@ -282,8 +280,6 @@ double ComputeRHEOKernel::calc_w_quintic(int i, int j, double delx, double dely, double ComputeRHEOKernel::calc_dw_quintic(int i, int j, double delx, double dely, double delz, double r, double *dW1, double *dW2) { double wp, tmp1, tmp2, tmp3, tmp1sq, tmp2sq, tmp3sq, s, wprinv; - double *mass = atom->mass; - int *type = atom->type; s = r * 3.0 * cutinv; @@ -326,9 +322,9 @@ double ComputeRHEOKernel::calc_w_wendlandc4(int i, int j, double delx, double de double w, tmp6, s; s = r * cutinv; - if (s > 1.0) { - w = 0.0; - } else { + if (s > 1.0) { + w = 0.0; + } else { tmp6 = (1.0 - s) * (1.0 - s); tmp6 *= tmp6 * tmp6; w = tmp6 * (1.0 + 6.0 * s + 35.0 * THIRD * s * s); @@ -347,14 +343,12 @@ double ComputeRHEOKernel::calc_w_wendlandc4(int i, int j, double delx, double de double ComputeRHEOKernel::calc_dw_wendlandc4(int i, int j, double delx, double dely, double delz, double r, double *dW1, double *dW2) { double wp, tmp1, tmp5, tmp6, s, wprinv; - double *mass = atom->mass; - int *type = atom->type; s = r * cutinv; - if (s > 1.0) { - wp = 0.0; - } else { + if (s > 1.0) { + wp = 0.0; + } else { tmp1 = 1.0 - s; tmp5 = tmp1 * tmp1; tmp5 = tmp5 * tmp5 * tmp1; @@ -395,7 +389,7 @@ double ComputeRHEOKernel::calc_w_rk0(int i, int j, double delx, double dely, dou double ComputeRHEOKernel::calc_w_rk1(int i, int j, double delx, double dely, double delz, double r) { int b; - double w, wR, dx[3], H[Mdim]; + double w, dx[3], H[Mdim]; dx[0] = delx; dx[1] = dely; @@ -437,7 +431,7 @@ double ComputeRHEOKernel::calc_w_rk1(int i, int j, double delx, double dely, dou double ComputeRHEOKernel::calc_w_rk2(int i, int j, double delx, double dely, double delz, double r) { int b; - double w, wR, dx[3], H[Mdim]; + double w, dx[3], H[Mdim]; dx[0] = delx; dx[1] = dely; dx[2] = delz; @@ -574,7 +568,7 @@ void ComputeRHEOKernel::compute_peratom() if (kernel_style == QUINTIC) return; corrections_calculated = 1; - int i, j, ii, jj, inum, jnum, itype, g, a, b, gsl_error; + int i, j, ii, jj, inum, jnum, a, b, gsl_error; double xtmp, ytmp, ztmp, r, rsq, w, vj, rhoj; double dx[3]; gsl_matrix_view gM; @@ -585,6 +579,7 @@ void ComputeRHEOKernel::compute_peratom() double **x = atom->x; int *type = atom->type; double *mass = atom->mass; + double *rmass = atom->rmass; double *rho = atom->rho; int *status = atom->rheo_status; tagint *tag = atom->tag; @@ -630,7 +625,8 @@ void ComputeRHEOKernel::compute_peratom() if (status[j] & PHASECHECK) rhoj = compute_interface->correct_rho(j, i); - vj = mass[type[j]] / rhoj; + if (rmass) vj = rmass[j] / rhoj; + else vj = mass[type[j]] / rhoj; M += w * vj; } } @@ -651,7 +647,6 @@ void ComputeRHEOKernel::compute_peratom() jlist = firstneigh[i]; jnum = numneigh[i]; - itype = type[i]; // Zero upper-triangle M and cut (will be symmetric): for (a = 0; a < Mdim; a++) { @@ -679,7 +674,8 @@ void ComputeRHEOKernel::compute_peratom() if (status[j] & PHASECHECK) rhoj = compute_interface->correct_rho(j, i); - vj = mass[type[j]] / rhoj; + if (rmass) vj = rmass[j] / rhoj; + else vj = mass[type[j]] / rhoj; //Populate the H-vector of polynomials (2D) if (dim == 2) { @@ -854,19 +850,17 @@ void ComputeRHEOKernel::grow_arrays(int nmax) int ComputeRHEOKernel::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i,j,k,m,a,b; - m = 0; - - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; if (comm_stage == 0) { buf[m++] = coordination[j]; } else { if (kernel_style == RK0) { buf[m++] = C0[j]; } else { - for (a = 0; a < ncor; a++) - for (b = 0; b < Mdim; b++) + for (int a = 0; a < ncor; a++) + for (int b = 0; b < Mdim; b++) buf[m++] = C[j][a][b]; } } @@ -878,19 +872,17 @@ int ComputeRHEOKernel::pack_forward_comm(int n, int *list, double *buf, void ComputeRHEOKernel::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last,a,b; - m = 0; - last = first + n; - - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { if (comm_stage == 0) { coordination[i] = buf[m++]; } else { if (kernel_style == RK0) { C0[i] = buf[m++]; } else { - for (a = 0; a < ncor; a++) - for (b = 0; b < Mdim; b++) + for (int a = 0; a < ncor; a++) + for (int b = 0; b < Mdim; b++) C[i][a][b] = buf[m++]; } } diff --git a/src/RHEO/compute_rheo_property_atom.cpp b/src/RHEO/compute_rheo_property_atom.cpp index 7a450e7708..cd0ff6692d 100644 --- a/src/RHEO/compute_rheo_property_atom.cpp +++ b/src/RHEO/compute_rheo_property_atom.cpp @@ -46,9 +46,10 @@ using namespace RHEO_NS; /* ---------------------------------------------------------------------- */ ComputeRHEOPropertyAtom::ComputeRHEOPropertyAtom(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg), fix_rheo(nullptr), fix_pressure(nullptr), fix_thermal(nullptr), compute_interface(nullptr), - compute_kernel(nullptr), compute_surface(nullptr), compute_vshift(nullptr), compute_grad(nullptr), - avec_index(nullptr), pack_choice(nullptr), col_index(nullptr) + Compute(lmp, narg, arg), avec_index(nullptr), col_index(nullptr), col_t_index(nullptr), + buf(nullptr), pack_choice(nullptr), fix_rheo(nullptr), fix_pressure(nullptr), + fix_thermal(nullptr), compute_interface(nullptr), compute_kernel(nullptr), + compute_surface(nullptr), compute_vshift(nullptr), compute_grad(nullptr) { if (narg < 4) utils::missing_cmd_args(FLERR, "compute property/atom", error); @@ -88,7 +89,6 @@ ComputeRHEOPropertyAtom::ComputeRHEOPropertyAtom(LAMMPS *lmp, int narg, char **a col_t_index = new int[nvalues]; int i = 0; - int index, a, b; for (int iarg = 3; iarg < narg; iarg++) { if (strcmp(arg[iarg], "phase") == 0) { pack_choice[i] = &ComputeRHEOPropertyAtom::pack_phase; @@ -536,7 +536,7 @@ int ComputeRHEOPropertyAtom::add_tensor_component(char* option, int i, FnPtrPack } shift = dim * dim; } else { - int index; + int index = -1; int dim_error = 0; if (utils::strmatch(option, "xx$")) { @@ -592,7 +592,7 @@ int ComputeRHEOPropertyAtom::add_vector_component(char* option, int i, FnPtrPack } shift = dim; } else { - int index; + int index = -1; if (utils::strmatch(option, "x$")) { index = 0; } else if (utils::strmatch(option, "y$")) { diff --git a/src/RHEO/compute_rheo_rho_sum.cpp b/src/RHEO/compute_rheo_rho_sum.cpp index 6e25b35374..8067bae21b 100644 --- a/src/RHEO/compute_rheo_rho_sum.cpp +++ b/src/RHEO/compute_rheo_rho_sum.cpp @@ -72,7 +72,7 @@ void ComputeRHEORhoSum::init_list(int /*id*/, NeighList *ptr) void ComputeRHEORhoSum::compute_peratom() { - int i, j, ii, jj, inum, jnum, itype, jtype; + int i, j, ii, jj, inum, jnum; double xtmp, ytmp, ztmp, delx, dely, delz; int *ilist, *jlist, *numneigh, **firstneigh; double rsq, w; @@ -81,12 +81,11 @@ void ComputeRHEORhoSum::compute_peratom() double **x = atom->x; double *rho = atom->rho; - int *type = atom->type; double *mass = atom->mass; + double *rmass = atom->rmass; + int *type = atom->type; int newton = force->newton; - double jmass; - inum = list->inum; ilist = list->ilist; numneigh = list->numneigh; @@ -96,7 +95,8 @@ void ComputeRHEORhoSum::compute_peratom() // initialize arrays, local with quintic self-contribution, ghosts are zeroed for (i = 0; i < nlocal; i++) { w = compute_kernel->calc_w_self(i, i); - rho[i] = w * mass[type[i]]; + if (rmass) rho[i] = w * rmass[i]; + else rho[i] = w * mass[type[i]]; } for (i = nlocal; i < nall; i++) rho[i] = 0.0; @@ -106,7 +106,6 @@ void ComputeRHEORhoSum::compute_peratom() xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; - itype = type[i]; jlist = firstneigh[i]; jnum = numneigh[i]; @@ -121,14 +120,26 @@ void ComputeRHEORhoSum::compute_peratom() if (rsq < cutsq) { w = compute_kernel->calc_w(i, j, delx, dely, delz, sqrt(rsq)); - if (self_mass_flag) { - rho[i] += w * mass[type[i]]; - if (newton || j < nlocal) - rho[j] += w * mass[type[j]]; + if (rmass) { + if (self_mass_flag) { + rho[i] += w * rmass[i]; + if (newton || j < nlocal) + rho[j] += w * rmass[j]; + } else { + rho[i] += w * rmass[j]; + if (newton || j < nlocal) + rho[j] += w * rmass[i]; + } } else { - rho[i] += w * mass[type[j]]; - if (newton || j < nlocal) - rho[j] += w * mass[type[i]]; + if (self_mass_flag) { + rho[i] += w * mass[type[i]]; + if (newton || j < nlocal) + rho[j] += w * mass[type[j]]; + } else { + rho[i] += w * mass[type[j]]; + if (newton || j < nlocal) + rho[j] += w * mass[type[i]]; + } } } } @@ -143,7 +154,7 @@ void ComputeRHEORhoSum::compute_peratom() int ComputeRHEORhoSum::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i, j, k, m; + int i, j, m; double *rho = atom->rho; m = 0; @@ -157,7 +168,7 @@ int ComputeRHEORhoSum::pack_forward_comm(int n, int *list, double *buf, /* ---------------------------------------------------------------------- */ void ComputeRHEORhoSum::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last; + int i, m, last; double *rho = atom->rho; m = 0; @@ -171,7 +182,7 @@ void ComputeRHEORhoSum::unpack_forward_comm(int n, int first, double *buf) int ComputeRHEORhoSum::pack_reverse_comm(int n, int first, double *buf) { - int i, k, m, last; + int i, m, last; double *rho = atom->rho; m = 0; @@ -186,7 +197,7 @@ int ComputeRHEORhoSum::pack_reverse_comm(int n, int first, double *buf) void ComputeRHEORhoSum::unpack_reverse_comm(int n, int *list, double *buf) { - int i, k, j, m; + int i, j, m; double *rho = atom->rho; m = 0; diff --git a/src/RHEO/compute_rheo_surface.cpp b/src/RHEO/compute_rheo_surface.cpp index c3a3774cdc..314d2d2ef5 100644 --- a/src/RHEO/compute_rheo_surface.cpp +++ b/src/RHEO/compute_rheo_surface.cpp @@ -43,8 +43,9 @@ static constexpr double EPSILON = 1e-10; /* ---------------------------------------------------------------------- */ ComputeRHEOSurface::ComputeRHEOSurface(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg), fix_rheo(nullptr), list(nullptr), rho0(nullptr), compute_kernel(nullptr), compute_interface(nullptr), - B(nullptr), gradC(nullptr), nsurface(nullptr), divr(nullptr), rsurface(nullptr) + Compute(lmp, narg, arg), nsurface(nullptr), rsurface(nullptr), divr(nullptr), fix_rheo(nullptr), + rho0(nullptr), B(nullptr), gradC(nullptr), list(nullptr), compute_kernel(nullptr), + compute_interface(nullptr) { if (narg != 3) error->all(FLERR,"Illegal compute RHEO/SURFACE command"); @@ -98,8 +99,8 @@ void ComputeRHEOSurface::init_list(int /*id*/, NeighList *ptr) void ComputeRHEOSurface::compute_peratom() { - int i, j, ii, jj, inum, jnum, a, b, itype, jtype, fluidi, fluidj; - double xtmp, ytmp, ztmp, rsq, Voli, Volj, rhoi, rhoj, wp; + int i, j, ii, jj, inum, jnum, a, itype, jtype, fluidi, fluidj; + double xtmp, ytmp, ztmp, rsq, Voli, Volj, rhoi, rhoj; double dWij[3], dWji[3], dx[3]; int *ilist, *jlist, *numneigh, **firstneigh; @@ -112,6 +113,7 @@ void ComputeRHEOSurface::compute_peratom() int *mask = atom->mask; int *type = atom->type; double *mass = atom->mass; + double *rmass = atom->rmass; double *rho = atom->rho; int *coordination = compute_kernel->coordination; @@ -172,10 +174,14 @@ void ComputeRHEOSurface::compute_peratom() } } - Voli = mass[itype] / rhoi; - Volj = mass[jtype] / rhoj; - - wp = compute_kernel->calc_dw_quintic(i, j, dx[0], dx[1], dx[2], sqrt(rsq), dWij, dWji); + if (rmass) { + Voli = rmass[i] / rhoi; + Volj = rmass[j] / rhoj; + } else { + Voli = mass[itype] / rhoi; + Volj = mass[jtype] / rhoj; + } + compute_kernel->calc_dw_quintic(i, j, dx[0], dx[1], dx[2], sqrt(rsq), dWij, dWji); for (a = 0; a < dim; a++) { divr[i] -= dWij[a] * dx[a] * Volj; @@ -310,17 +316,15 @@ void ComputeRHEOSurface::compute_peratom() int ComputeRHEOSurface::pack_reverse_comm(int n, int first, double *buf) { - int i,a,b,k,m,last; int dim = domain->dimension; int *status = atom->rheo_status; - - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { if (comm_stage == 0) { buf[m++] = divr[i]; - for (a = 0; a < dim; a ++ ) - for (b = 0; b < dim; b ++) + for (int a = 0; a < dim; a ++ ) + for (int b = 0; b < dim; b ++) buf[m++] = gradC[i][a * dim + b]; } else if (comm_stage == 1) { buf[m++] = (double) status[i]; @@ -334,27 +338,23 @@ int ComputeRHEOSurface::pack_reverse_comm(int n, int first, double *buf) void ComputeRHEOSurface::unpack_reverse_comm(int n, int *list, double *buf) { - int i,a,b,k,j,m; int dim = domain->dimension; int *status = atom->rheo_status; - int tmp1; - double tmp2; - - m = 0; - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; if (comm_stage == 0) { divr[j] += buf[m++]; - for (a = 0; a < dim; a ++ ) - for (b = 0; b < dim; b ++) + for (int a = 0; a < dim; a ++ ) + for (int b = 0; b < dim; b ++) gradC[j][a * dim + b] += buf[m++]; } else if (comm_stage == 1) { - tmp1 = (int) buf[m++]; + auto tmp1 = (int) buf[m++]; if ((status[j] & STATUS_BULK) && (tmp1 & STATUS_LAYER)) { status[j] &= SURFACEMASK; status[j] |= STATUS_LAYER; } - tmp2 = buf[m++]; + auto tmp2 = buf[m++]; rsurface[j] = MIN(rsurface[j], tmp2); } } @@ -366,12 +366,10 @@ void ComputeRHEOSurface::unpack_reverse_comm(int n, int *list, double *buf) int ComputeRHEOSurface::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i,j,a,b,k,m; int *status = atom->rheo_status; - m = 0; - - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; if (comm_stage == 0) { buf[m++] = divr[j]; } else if (comm_stage == 1) { @@ -386,12 +384,10 @@ int ComputeRHEOSurface::pack_forward_comm(int n, int *list, double *buf, void ComputeRHEOSurface::unpack_forward_comm(int n, int first, double *buf) { - int i, k, a, b, m, last; int *status = atom->rheo_status; - - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { if (comm_stage == 0) { divr[i] = buf[m++]; } else if (comm_stage == 1) { diff --git a/src/RHEO/compute_rheo_vshift.cpp b/src/RHEO/compute_rheo_vshift.cpp index c06ef533ac..8d2f04d7f6 100644 --- a/src/RHEO/compute_rheo_vshift.cpp +++ b/src/RHEO/compute_rheo_vshift.cpp @@ -41,8 +41,8 @@ using namespace RHEO_NS; /* ---------------------------------------------------------------------- */ ComputeRHEOVShift::ComputeRHEOVShift(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg), list(nullptr), vshift(nullptr), fix_rheo(nullptr), - compute_kernel(nullptr), compute_interface(nullptr), compute_surface(nullptr) + Compute(lmp, narg, arg), vshift(nullptr), fix_rheo(nullptr), rho0(nullptr), list(nullptr), + compute_interface(nullptr), compute_kernel(nullptr), compute_surface(nullptr) { if (narg != 3) error->all(FLERR,"Illegal compute RHEO/VShift command"); @@ -90,7 +90,7 @@ void ComputeRHEOVShift::init_list(int /*id*/, NeighList *ptr) void ComputeRHEOVShift::compute_peratom() { - int i, j, a, b, ii, jj, jnum, itype, jtype; + int i, j, a, ii, jj, jnum, itype, jtype; int fluidi, fluidj; double xtmp, ytmp, ztmp, rsq, r, rinv; double w, wp, dr, w0, w4, vmag, prefactor; @@ -108,6 +108,7 @@ void ComputeRHEOVShift::compute_peratom() double **v = atom->v; double *rho = atom->rho; double *mass = atom->mass; + double *rmass = atom->rmass; int nlocal = atom->nlocal; int nall = nlocal + atom->nghost; @@ -140,7 +141,8 @@ void ComputeRHEOVShift::compute_peratom() itype = type[i]; jlist = firstneigh[i]; jnum = numneigh[i]; - imass = mass[itype]; + if (rmass) imass = rmass[i]; + else imass = mass[itype]; fluidi = !(status[i] & PHASECHECK); for (jj = 0; jj < jnum; jj++) { @@ -160,7 +162,8 @@ void ComputeRHEOVShift::compute_peratom() if (rsq < cutsq) { jtype = type[j]; - jmass = mass[jtype]; + if (rmass) jmass = rmass[j]; + else jmass = mass[jtype]; r = sqrt(rsq); rinv = 1 / r; @@ -229,8 +232,6 @@ void ComputeRHEOVShift::correct_surfaces() { if (!surface_flag) return; - int i, a, b; - int *status = atom->rheo_status; int *mask = atom->mask; double **nsurface = compute_surface->nsurface; @@ -239,7 +240,7 @@ void ComputeRHEOVShift::correct_surfaces() int dim = domain->dimension; double nx, ny, nz, vx, vy, vz, dot; - for (i = 0; i < nlocal; i++) { + for (int i = 0; i < nlocal; i++) { if (mask[i] & groupbit) { if (status[i] & PHASECHECK) continue; @@ -283,11 +284,11 @@ void ComputeRHEOVShift::correct_surfaces() int ComputeRHEOVShift::pack_reverse_comm(int n, int first, double *buf) { - int i,m,last; + int m, last; m = 0; last = first + n; - for (i = first; i < last; i++) { + for (int i = first; i < last; i++) { buf[m++] = vshift[i][0]; buf[m++] = vshift[i][1]; buf[m++] = vshift[i][2]; diff --git a/src/RHEO/fix_rheo.cpp b/src/RHEO/fix_rheo.cpp index f70b9e121f..b093eff681 100644 --- a/src/RHEO/fix_rheo.cpp +++ b/src/RHEO/fix_rheo.cpp @@ -50,8 +50,9 @@ static const char cite_rheo[] = /* ---------------------------------------------------------------------- */ FixRHEO::FixRHEO(LAMMPS *lmp, int narg, char **arg) : - Fix(lmp, narg, arg), compute_grad(nullptr), compute_kernel(nullptr), compute_surface(nullptr), - compute_interface(nullptr), compute_rhosum(nullptr), compute_vshift(nullptr), rho0(nullptr), csq(nullptr) + Fix(lmp, narg, arg), rho0(nullptr), csq(nullptr), compute_grad(nullptr), compute_kernel(nullptr), + compute_interface(nullptr), compute_surface(nullptr), compute_rhosum(nullptr), + compute_vshift(nullptr) { time_integrate = 1; @@ -457,7 +458,6 @@ void FixRHEO::final_integrate() double dtfm, divu; int i, a; - double **x = atom->x; double **v = atom->v; double **f = atom->f; double **gradv = compute_grad->gradv; @@ -469,7 +469,6 @@ void FixRHEO::final_integrate() int *mask = atom->mask; int *status = atom->rheo_status; - int rmass_flag = atom->rmass_flag; int dim = domain->dimension; // Update velocity @@ -477,7 +476,7 @@ void FixRHEO::final_integrate() if (mask[i] & groupbit) { if (status[i] & STATUS_NO_INTEGRATION) continue; - if (rmass_flag) { + if (rmass) { dtfm = dtf / rmass[i]; } else { dtfm = dtf / mass[type[i]]; diff --git a/src/RHEO/fix_rheo_oxidation.cpp b/src/RHEO/fix_rheo_oxidation.cpp index a51f2feb95..699f3428c8 100644 --- a/src/RHEO/fix_rheo_oxidation.cpp +++ b/src/RHEO/fix_rheo_oxidation.cpp @@ -155,7 +155,6 @@ void FixRHEOOxidation::post_integrate() double delx, dely, delz, rsq; tagint tagi, tagj; - int nlocal = atom->nlocal; int newton_bond = force->newton_bond; tagint *tag = atom->tag; @@ -267,7 +266,7 @@ void FixRHEOOxidation::post_force(int /*vflag*/) int FixRHEOOxidation::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i, j, k, m; + int i, j, m; double **x = atom->x; m = 0; @@ -284,7 +283,7 @@ int FixRHEOOxidation::pack_forward_comm(int n, int *list, double *buf, void FixRHEOOxidation::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last; + int i, m, last; double **x = atom->x; m = 0; last = first + n; diff --git a/src/RHEO/fix_rheo_pressure.cpp b/src/RHEO/fix_rheo_pressure.cpp index 82adf52012..84b0825570 100644 --- a/src/RHEO/fix_rheo_pressure.cpp +++ b/src/RHEO/fix_rheo_pressure.cpp @@ -39,7 +39,9 @@ static constexpr double SEVENTH = 1.0 / 7.0; /* ---------------------------------------------------------------------- */ FixRHEOPressure::FixRHEOPressure(LAMMPS *lmp, int narg, char **arg) : - Fix(lmp, narg, arg), fix_rheo(nullptr), rho0(nullptr), csq(nullptr), rho0inv(nullptr), csqinv(nullptr), c_cubic(nullptr), tpower(nullptr), pbackground(nullptr), pressure_style(nullptr) + Fix(lmp, narg, arg), c_cubic(nullptr), csq(nullptr), csqinv(nullptr), rho0(nullptr), + rho0inv(nullptr), tpower(nullptr), pbackground(nullptr), pressure_style(nullptr), + fix_rheo(nullptr) { if (narg < 4) error->all(FLERR,"Illegal fix command"); @@ -160,9 +162,6 @@ void FixRHEOPressure::setup_pre_force(int /*vflag*/) void FixRHEOPressure::pre_force(int /*vflag*/) { - int i; - double dr, rr3, rho_ratio; - int *mask = atom->mask; int *type = atom->type; double *rho = atom->rho; @@ -170,7 +169,7 @@ void FixRHEOPressure::pre_force(int /*vflag*/) int nlocal = atom->nlocal; - for (i = 0; i < nlocal; i++) + for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) pressure[i] = calc_pressure(rho[i], type[i]); @@ -182,12 +181,10 @@ void FixRHEOPressure::pre_force(int /*vflag*/) int FixRHEOPressure::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i,j,k,m; double *pressure = atom->pressure; - m = 0; - - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; buf[m++] = pressure[j]; } return m; @@ -197,12 +194,10 @@ int FixRHEOPressure::pack_forward_comm(int n, int *list, double *buf, void FixRHEOPressure::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last; double *pressure = atom->pressure; - - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { pressure[i] = buf[m++]; } } @@ -211,7 +206,8 @@ void FixRHEOPressure::unpack_forward_comm(int n, int first, double *buf) double FixRHEOPressure::calc_pressure(double rho, int type) { - double p, dr, rr3, rho_ratio; + double p = 0.0; + double dr, rr3, rho_ratio; if (pressure_style[type] == LINEAR) { p = csq[type] * (rho - rho0[type]); @@ -234,7 +230,7 @@ double FixRHEOPressure::calc_pressure(double rho, int type) double FixRHEOPressure::calc_rho(double p, int type) { - double rho, dr, rr3, rho_ratio; + double rho = 0.0; if (pressure_style[type] == LINEAR) { rho = csqinv[type] * p + rho0[type]; diff --git a/src/RHEO/fix_rheo_thermal.cpp b/src/RHEO/fix_rheo_thermal.cpp index c8c1877e32..a7d9fc8df9 100644 --- a/src/RHEO/fix_rheo_thermal.cpp +++ b/src/RHEO/fix_rheo_thermal.cpp @@ -60,10 +60,9 @@ static const char cite_rheo_oxide[] = /* ---------------------------------------------------------------------- */ FixRHEOThermal::FixRHEOThermal(LAMMPS *lmp, int narg, char **arg) : - Fix(lmp, narg, arg), fix_rheo(nullptr), compute_grad(nullptr), compute_vshift(nullptr), - Tc(nullptr), kappa(nullptr), cv(nullptr), L(nullptr), - Tc_style(nullptr), kappa_style(nullptr), cv_style(nullptr), L_style(nullptr), - fix_update_special_bonds(nullptr) + Fix(lmp, narg, arg), cv(nullptr), Tc(nullptr), kappa(nullptr), L(nullptr), cv_style(nullptr), + Tc_style(nullptr), kappa_style(nullptr), L_style(nullptr), list(nullptr), fix_rheo(nullptr), + compute_grad(nullptr), compute_vshift(nullptr), fix_update_special_bonds(nullptr) { if (narg < 4) error->all(FLERR,"Illegal fix command"); @@ -189,13 +188,14 @@ FixRHEOThermal::FixRHEOThermal(LAMMPS *lmp, int narg, char **arg) : cut_bond = utils::numeric(FLERR, arg[iarg + 1], false, lmp); btype = utils::numeric(FLERR, arg[iarg + 2], false, lmp); comm_forward = 4; - if (cut_bond <= 0.0) error->all(FLERR, "Illegal max bond length must be greater than zero");\ - if (btype < 1 || btype > atom->nbondtypes) error->all(FLERR, "Illegal value for bond type"); + if (cut_bond <= 0.0) error->all(FLERR, "Illegal max bond length must be greater than zero"); + if ((btype < 1) || (btype > atom->nbondtypes)) + error->all(FLERR, "Illegal value {} for bond type", btype); cutsq_bond = cut_bond * cut_bond; iarg += 3; } else { - error->all(FLERR,"Illegal fix command, {}", arg[iarg]); + error->all(FLERR,"Unknown fix rheo/thermal keyword: {}", arg[iarg]); } } @@ -460,20 +460,16 @@ void FixRHEOThermal::post_neighbor() void FixRHEOThermal::pre_force(int /*vflag*/) { - int i, itype; double cvi, Tci, Ti, Li; double *energy = atom->esph; double *temperature = atom->temperature; int *type = atom->type; - int *status = atom->rheo_status; - - int nlocal = atom->nlocal; - int nall = nlocal + atom->nghost; + int nall = atom->nlocal + atom->nghost; // Calculate temperature - for (i = 0; i < nall; i++) { - itype = type[i]; + for (int i = 0; i < nall; i++) { + int itype = type[i]; cvi = calc_cv(i, itype); temperature[i] = energy[i] / cvi; @@ -529,7 +525,6 @@ void FixRHEOThermal::break_bonds() int nbondlist = neighbor->nbondlist; int nlocal = atom->nlocal; - int nall = nlocal + atom->nghost; // Delete all bonds for local atoms that melt of a given type for (int i = 0; i < nlocal; i++) { @@ -753,13 +748,11 @@ double FixRHEOThermal::calc_L(int i, int itype) int FixRHEOThermal::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i, j, k, m; int *status = atom->rheo_status; double **x = atom->x; - m = 0; - - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; buf[m++] = ubuf(status[j]).d; buf[m++] = x[j][0]; buf[m++] = x[j][1]; @@ -772,12 +765,11 @@ int FixRHEOThermal::pack_forward_comm(int n, int *list, double *buf, void FixRHEOThermal::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last; int *status = atom->rheo_status; double **x = atom->x; - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { status[i] = (int) ubuf(buf[m++]).i; x[i][0] = buf[m++]; x[i][1] = buf[m++]; diff --git a/src/RHEO/fix_rheo_viscosity.cpp b/src/RHEO/fix_rheo_viscosity.cpp index 3ca7fd8d13..5f7a1208c5 100644 --- a/src/RHEO/fix_rheo_viscosity.cpp +++ b/src/RHEO/fix_rheo_viscosity.cpp @@ -38,8 +38,8 @@ enum {NONE, CONSTANT, POWER}; /* ---------------------------------------------------------------------- */ FixRHEOViscosity::FixRHEOViscosity(LAMMPS *lmp, int narg, char **arg) : - Fix(lmp, narg, arg), fix_rheo(nullptr), compute_grad(nullptr), eta(nullptr), - npow(nullptr), K(nullptr), gd0(nullptr), tau0(nullptr), viscosity_style(nullptr) + Fix(lmp, narg, arg), eta(nullptr), npow(nullptr), K(nullptr), gd0(nullptr), tau0(nullptr), + viscosity_style(nullptr), fix_rheo(nullptr), compute_grad(nullptr) { if (narg < 4) error->all(FLERR, "Illegal fix command"); @@ -225,12 +225,10 @@ void FixRHEOViscosity::pre_force(int /*vflag*/) int FixRHEOViscosity::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) { - int i, j, k, m; double *viscosity = atom->viscosity; - m = 0; - - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; buf[m++] = viscosity[j]; } return m; @@ -240,12 +238,10 @@ int FixRHEOViscosity::pack_forward_comm(int n, int *list, double *buf, void FixRHEOViscosity::unpack_forward_comm(int n, int first, double *buf) { - int i, k, m, last; double *viscosity = atom->viscosity; - - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { viscosity[i] = buf[m++]; } } diff --git a/src/RHEO/pair_rheo.cpp b/src/RHEO/pair_rheo.cpp index 9ebf884b6e..520a1c4470 100644 --- a/src/RHEO/pair_rheo.cpp +++ b/src/RHEO/pair_rheo.cpp @@ -47,8 +47,9 @@ static constexpr double EPSILON = 1e-2; /* ---------------------------------------------------------------------- */ PairRHEO::PairRHEO(LAMMPS *lmp) : - Pair(lmp), compute_kernel(nullptr), compute_grad(nullptr), compute_interface(nullptr), fix_rheo(nullptr), - fix_pressure(nullptr), rho0(nullptr), csq(nullptr), cs(nullptr) + Pair(lmp), csq(nullptr), rho0(nullptr), cs(nullptr), compute_kernel(nullptr), + compute_grad(nullptr), compute_interface(nullptr), fix_rheo(nullptr), + fix_pressure(nullptr) { restartinfo = 0; single_enable = 0; @@ -80,10 +81,10 @@ void PairRHEO::compute(int eflag, int vflag) int i, j, a, b, ii, jj, inum, jnum, itype, jtype; int pair_force_flag, pair_rho_flag, pair_avisc_flag; int fluidi, fluidj; - double xtmp, ytmp, ztmp, w, wp, Ti, Tj, dT, csq_ave, cs_ave; + double xtmp, ytmp, ztmp, wp, Ti, Tj, dT, csq_ave, cs_ave; double rhoi, rhoj, rho0i, rho0j, voli, volj, Pi, Pj, etai, etaj, kappai, kappaj, eta_ave, kappa_ave, dT_prefactor; double mu, q, fp_prefactor, drho_damp, fmag, psi_ij, Fij; - double *dWij, *dWji, *dW1ij, *dW1ji; + double *dWij, *dWji; double dx[3], du[3], dv[3], fv[3], dfp[3], fsolid[3], ft[3], vi[3], vj[3]; int *ilist, *jlist, *numneigh, **firstneigh; @@ -102,16 +103,15 @@ void PairRHEO::compute(int eflag, int vflag) double **f = atom->f; double *rho = atom->rho; double *mass = atom->mass; + double *rmass = atom->rmass; double *drho = atom->drho; double *pressure = atom->pressure; double *viscosity = atom->viscosity; double *conductivity = atom->conductivity; double *temperature = atom->temperature; double *heatflow = atom->heatflow; - double *special_lj = force->special_lj; int *type = atom->type; int *status = atom->rheo_status; - tagint *tag = atom->tag; double **fp_store, *chi; if (compute_interface) { @@ -150,7 +150,8 @@ void PairRHEO::compute(int eflag, int vflag) itype = type[i]; jlist = firstneigh[i]; jnum = numneigh[i]; - imass = mass[itype]; + if (rmass) imass = rmass[i]; + else imass = mass[itype]; etai = viscosity[i]; fluidi = !(status[i] & PHASECHECK); if (thermal_flag) { @@ -172,7 +173,8 @@ void PairRHEO::compute(int eflag, int vflag) r = sqrt(rsq); rinv = 1 / r; - jmass = mass[jtype]; + if (rmass) jmass = rmass[i]; + else jmass = mass[jtype]; etaj = viscosity[j]; fluidj = !(status[j] & PHASECHECK); if (thermal_flag) { @@ -461,11 +463,13 @@ void PairRHEO::setup() { auto fixes = modify->get_fix_by_style("rheo"); if (fixes.size() == 0) error->all(FLERR, "Need to define fix rheo to use pair rheo"); + if (fixes.size() > 1) error->all(FLERR, "Must have only one fix rheo defined"); fix_rheo = dynamic_cast(fixes[0]); // Currently only allow one instance of fix rheo/pressure fixes = modify->get_fix_by_style("rheo/pressure"); if (fixes.size() == 0) error->all(FLERR, "Need to define fix rheo/pressure to use pair rheo"); + if (fixes.size() > 1) error->all(FLERR, "Must have only one fix rheo/pressure defined"); fix_pressure = dynamic_cast(fixes[0]); compute_kernel = fix_rheo->compute_kernel; @@ -518,12 +522,10 @@ double PairRHEO::init_one(int i, int j) int PairRHEO::pack_reverse_comm(int n, int first, double *buf) { - int i, k, m, last; double **fp_store = compute_interface->fp_store; - - m = 0; - last = first + n; - for (i = first; i < last; i++) { + int m = 0; + int last = first + n; + for (int i = first; i < last; i++) { buf[m++] = fp_store[i][0]; buf[m++] = fp_store[i][1]; buf[m++] = fp_store[i][2]; @@ -536,12 +538,10 @@ int PairRHEO::pack_reverse_comm(int n, int first, double *buf) void PairRHEO::unpack_reverse_comm(int n, int *list, double *buf) { - int i, j, k, m; double **fp_store = compute_interface->fp_store; - - m = 0; - for (i = 0; i < n; i++) { - j = list[i]; + int m = 0; + for (int i = 0; i < n; i++) { + int j = list[i]; fp_store[j][0] += buf[m++]; fp_store[j][1] += buf[m++]; fp_store[j][2] += buf[m++]; diff --git a/src/compute_bond_local.cpp b/src/compute_bond_local.cpp index 9ed591f73f..e3c3e26a85 100644 --- a/src/compute_bond_local.cpp +++ b/src/compute_bond_local.cpp @@ -13,21 +13,23 @@ ------------------------------------------------------------------------- */ #include "compute_bond_local.h" -#include -#include + #include "atom.h" #include "atom_vec.h" -#include "molecule.h" -#include "update.h" -#include "domain.h" -#include "force.h" #include "bond.h" #include "comm.h" +#include "domain.h" +#include "error.h" +#include "force.h" #include "input.h" -#include "variable.h" #include "math_extra.h" #include "memory.h" -#include "error.h" +#include "molecule.h" +#include "update.h" +#include "variable.h" + +#include +#include using namespace LAMMPS_NS; @@ -80,8 +82,8 @@ ComputeBondLocal::ComputeBondLocal(LAMMPS *lmp, int narg, char **arg) : bstyle[nvalues++] = VARIABLE; vstr[nvar] = utils::strdup(&arg[iarg][2]); nvar++; - } else if (arg[iarg][0] == 'b') { - int n = atoi(&arg[iarg][1]); + } else if (utils::strmatch(arg[iarg], "^b\\d+$")) { // b1, b2, b3, ... bN + int n = std::stoi(&arg[iarg][1]); if (n <= 0) error->all(FLERR, "Invalid keyword {} in compute bond/local command", arg[iarg]); bstyle[nvalues] = BN; bindex[nvalues++] = n - 1; @@ -96,13 +98,13 @@ ComputeBondLocal::ComputeBondLocal(LAMMPS *lmp, int narg, char **arg) : while (iarg < narg) { if (strcmp(arg[iarg],"set") == 0) { setflag = 1; - if (iarg+3 > narg) error->all(FLERR,"Illegal compute bond/local command"); + if (iarg+3 > narg) utils::missing_cmd_args(FLERR,"compute bond/local set", error); if (strcmp(arg[iarg+1],"dist") == 0) { delete [] dstr; dstr = utils::strdup(arg[iarg+2]); - } else error->all(FLERR,"Illegal compute bond/local command"); + } else error->all(FLERR,"Unknown compute bond/local set keyword: {}", arg[iarg+2]); iarg += 3; - } else error->all(FLERR,"Illegal compute bond/local command"); + } else error->all(FLERR,"Unknown compute bond/local keyword: {}", arg[iarg]); } // error check @@ -113,9 +115,9 @@ ComputeBondLocal::ComputeBondLocal(LAMMPS *lmp, int narg, char **arg) : for (int i = 0; i < nvar; i++) { vvar[i] = input->variable->find(vstr[i]); if (vvar[i] < 0) - error->all(FLERR,"Variable name for copute bond/local does not exist"); + error->all(FLERR,"Variable name {} for copute bond/local does not exist", vstr[i]); if (!input->variable->equalstyle(vvar[i])) - error->all(FLERR,"Variable for compute bond/local is invalid style"); + error->all(FLERR,"Variable {} for compute bond/local is invalid style", vstr[i]); } if (dstr) { @@ -126,7 +128,7 @@ ComputeBondLocal::ComputeBondLocal(LAMMPS *lmp, int narg, char **arg) : error->all(FLERR,"Variable for compute bond/local is invalid style"); } } else if (setflag) - error->all(FLERR,"Compute bond/local set with no variable"); + error->all(FLERR,"Compute bond/local set used with without a variable"); // set singleflag if need to call bond->single() diff --git a/src/compute_pair_local.cpp b/src/compute_pair_local.cpp index 88991f7481..351b499468 100644 --- a/src/compute_pair_local.cpp +++ b/src/compute_pair_local.cpp @@ -66,15 +66,13 @@ ComputePairLocal::ComputePairLocal(LAMMPS *lmp, int narg, char **arg) : pstyle[nvalues++] = DY; else if (strcmp(arg[iarg], "dz") == 0) pstyle[nvalues++] = DZ; - else if (arg[iarg][0] == 'p') { - int n = atoi(&arg[iarg][1]); + else if (utils::strmatch(arg[iarg], "^p\\d+$")) { // p1, p2, p3, ... pN + int n = std::stoi(&arg[iarg][1]); if (n <= 0) error->all(FLERR, "Invalid keyword {} in compute pair/local command", arg[iarg]); pstyle[nvalues] = PN; pindex[nvalues++] = n - 1; - } else break; - iarg++; } @@ -84,22 +82,22 @@ ComputePairLocal::ComputePairLocal(LAMMPS *lmp, int narg, char **arg) : while (iarg < narg) { if (strcmp(arg[iarg], "cutoff") == 0) { - if (iarg + 2 > narg) error->all(FLERR, "Illegal compute pair/local command"); + if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "compute pair/local cutoff", error); if (strcmp(arg[iarg + 1], "type") == 0) cutstyle = TYPE; else if (strcmp(arg[iarg + 1], "radius") == 0) cutstyle = RADIUS; else - error->all(FLERR, "Illegal compute pair/local command"); + error->all(FLERR, "Unknown compute pair/local cutoff keyword: {}", arg[iarg + 1]); iarg += 2; } else - error->all(FLERR, "Illegal compute pair/local command"); + error->all(FLERR, "Unknown compute pair/local keyword: {}", arg[iarg]); } // error check if (cutstyle == RADIUS && !atom->radius_flag) - error->all(FLERR, "Compute pair/local requires atom attribute radius"); + error->all(FLERR, "This compute pair/local requires atom attribute radius"); // set singleflag if need to call pair->single() diff --git a/src/label_map.cpp b/src/label_map.cpp index a9ffd60877..174f16a7ec 100644 --- a/src/label_map.cpp +++ b/src/label_map.cpp @@ -369,35 +369,35 @@ void LabelMap::read_restart(FILE *fp) for (int i = 0; i < natomtypes; i++) { charlabel = read_string(fp); typelabel[i] = charlabel; - typelabel_map[charlabel] = i + 1; + if (strlen(charlabel) > 0) typelabel_map[charlabel] = i + 1; delete[] charlabel; } for (int i = 0; i < nbondtypes; i++) { charlabel = read_string(fp); btypelabel[i] = charlabel; - btypelabel_map[charlabel] = i + 1; + if (strlen(charlabel) > 0) btypelabel_map[charlabel] = i + 1; delete[] charlabel; } for (int i = 0; i < nangletypes; i++) { charlabel = read_string(fp); atypelabel[i] = charlabel; - atypelabel_map[charlabel] = i + 1; + if (strlen(charlabel) > 0) atypelabel_map[charlabel] = i + 1; delete[] charlabel; } for (int i = 0; i < ndihedraltypes; i++) { charlabel = read_string(fp); dtypelabel[i] = charlabel; - dtypelabel_map[charlabel] = i + 1; + if (strlen(charlabel) > 0) dtypelabel_map[charlabel] = i + 1; delete[] charlabel; } for (int i = 0; i < nimpropertypes; i++) { charlabel = read_string(fp); itypelabel[i] = charlabel; - itypelabel_map[charlabel] = i + 1; + if (strlen(charlabel) > 0) itypelabel_map[charlabel] = i + 1; delete[] charlabel; } } diff --git a/src/lammps.cpp b/src/lammps.cpp index b3659fdf50..8e81f785de 100644 --- a/src/lammps.cpp +++ b/src/lammps.cpp @@ -191,7 +191,7 @@ LAMMPS::LAMMPS(int narg, char **arg, MPI_Comm communicator) : int me,nprocs; MPI_Comm_rank(communicator,&me); MPI_Comm_size(communicator,&nprocs); - int color = atoi(arg[iarg+1]); + int color = std::stoi(arg[iarg+1]); MPI_Comm subcomm; MPI_Comm_split(communicator,color,me,&subcomm); external_comm = communicator; diff --git a/src/library.cpp b/src/library.cpp index 5af39ce910..32577e47b4 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -2532,7 +2532,7 @@ a char pointer and it should **not** be deallocated. Example: * * \param handle pointer to a previously created LAMMPS instance * \param name name of the variable - * \param group group-ID for atom style variable or ``NULL`` + * \param group group-ID for atom style variable or ``NULL`` or non-NULL to get vector length * \return pointer (cast to ``void *``) to the location of the * requested data or ``NULL`` if not found. */ diff --git a/src/lmptype.h b/src/lmptype.h index 4e62a1a7eb..d2181c9898 100644 --- a/src/lmptype.h +++ b/src/lmptype.h @@ -74,10 +74,8 @@ namespace LAMMPS_NS { #ifdef LAMMPS_LONGLONG_TO_LONG #define MPI_LL MPI_LONG -#define ATOLL atoll #else #define MPI_LL MPI_LONG_LONG -#define ATOLL atol #endif // for atomic problems that exceed 2 billion (2^31) atoms @@ -103,9 +101,6 @@ typedef int64_t bigint; #define TAGINT_FORMAT "%d" #define BIGINT_FORMAT "%" PRId64 -#define ATOTAGINT atoi -#define ATOBIGINT ATOLL - #define LAMMPS_TAGINT LAMMPS_INT #define LAMMPS_TAGINT_2D LAMMPS_INT_2D #define LAMMPS_BIGINT LAMMPS_INT64 @@ -141,9 +136,6 @@ typedef int64_t bigint; #define TAGINT_FORMAT "%" PRId64 #define BIGINT_FORMAT "%" PRId64 -#define ATOTAGINT ATOLL -#define ATOBIGINT ATOLL - #define LAMMPS_TAGINT LAMMPS_INT64 #define LAMMPS_TAGINT_2D LAMMPS_INT64_2D #define LAMMPS_BIGINT LAMMPS_INT64 @@ -178,9 +170,6 @@ typedef int bigint; #define TAGINT_FORMAT "%d" #define BIGINT_FORMAT "%d" -#define ATOTAGINT atoi -#define ATOBIGINT atoi - #define LAMMPS_TAGINT LAMMPS_INT #define LAMMPS_TAGINT_2D LAMMPS_INT_2D #define LAMMPS_BIGINT LAMMPS_INT diff --git a/src/lmpwindows.h b/src/lmpwindows.h deleted file mode 100644 index 34e90a140f..0000000000 --- a/src/lmpwindows.h +++ /dev/null @@ -1,60 +0,0 @@ -/* -*- c++ -*- ---------------------------------------------------------- - LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator - https://www.lammps.org/, Sandia National Laboratories - LAMMPS development team: developers@lammps.org - - Copyright (2003) Sandia Corporation. Under the terms of Contract - DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains - certain rights in this software. This software is distributed under - the GNU General Public License. - - See the README file in the top-level LAMMPS directory. -------------------------------------------------------------------------- */ - -#include -#include - -// some symbols have different names in Windows - -#undef ATOBIGINT -#define ATOBIGINT _atoi64 - -#define pclose _pclose -#define strdup _strdup - -// the following functions are defined to get rid of -// 'ambiguous call to overloaded function' error in VSS for mismatched type arguments -#if !defined(__MINGW32__) -inline double pow(int i, int j) -{ - return pow((double) i, j); -} -inline double fabs(int i) -{ - return fabs((double) i); -} -inline double sqrt(int i) -{ - return sqrt((double) i); -} -#endif - -inline double trunc(double x) -{ - return x > 0 ? floor(x) : ceil(x); -} - -// Windows version of mkdir function does not have permission flags -#ifndef S_IRWXU -#define S_IRWXU 0 -#endif -#ifndef S_IRGRP -#define S_IRGRP 0 -#endif -#ifndef S_IXGRP -#define S_IXGRP 0 -#endif -inline int mkdir(const char *path, int) -{ - return _mkdir(path); -} diff --git a/src/molecule.cpp b/src/molecule.cpp index 55dfb4d15d..617402d605 100644 --- a/src/molecule.cpp +++ b/src/molecule.cpp @@ -599,14 +599,14 @@ void Molecule::read(int flag) } else if (keyword == "Dihedrals") { if (ndihedrals == 0) error->all(FLERR, - "Found Dihedrals section" + "Found Dihedrals section " "but no ndihedrals setting in header"); dihedralflag = tag_require = 1; dihedrals(flag, line); } else if (keyword == "Impropers") { if (nimpropers == 0) error->all(FLERR, - "Found Impropers section" + "Found Impropers section " "but no nimpropers setting in header"); improperflag = tag_require = 1; impropers(flag, line); diff --git a/src/neighbor.h b/src/neighbor.h index cec6ae45c1..5ac40635b5 100644 --- a/src/neighbor.h +++ b/src/neighbor.h @@ -212,7 +212,6 @@ class Neighbor : protected Pointers { int npair_perpetual; // # of perpetual NeighPair classes int *slist; // indices of them in neigh_stencil int *plist; // indices of them in neigh_pair - int npair_occasional; // # of occasional NeighPair classes int *olist; // indices of them in neigh_pair diff --git a/src/read_data.cpp b/src/read_data.cpp index d83456e985..29696bd62d 100644 --- a/src/read_data.cpp +++ b/src/read_data.cpp @@ -535,7 +535,7 @@ void ReadData::command(int narg, char **arg) if (domain->dimension == 2) { if (triclinic_general == 0) { if (boxlo[2] >= 0.0 || boxhi[2] <= 0.0) - error->all(FLERR,"Read_data zlo/zhi for 2d simulation must straddle 0.0"); + error->all(FLERR, "Read_data zlo/zhi for 2d simulation must straddle 0.0"); } else if (triclinic_general == 1) { if (cvec[0] != 0.0 || cvec[1] != 0.0 || cvec[2] != 1.0 || abc_origin[2] != -0.5) error->all(FLERR,"Read_data cvec and/or abc_origin is invalid for " @@ -594,12 +594,12 @@ void ReadData::command(int narg, char **arg) domain->yz = yz; } - // general triclinic box - // define_general_triclinic() converts - // ABC edge vectors + abc_origin to restricted triclinic + // general triclinic box + // define_general_triclinic() converts + // ABC edge vectors + abc_origin to restricted triclinic } else if (triclinic_general) { - domain->define_general_triclinic(avec,bvec,cvec,abc_origin); + domain->define_general_triclinic(avec, bvec, cvec, abc_origin); } } @@ -615,7 +615,7 @@ void ReadData::command(int narg, char **arg) if (triclinic_general) { if (!domain->triclinic_general) - error->all(FLERR,"Read_data subsequent file cannot switch to general triclinic"); + error->all(FLERR, "Read_data subsequent file cannot switch to general triclinic"); int errflag = 0; if (avec[0] != domain->avec[0] || avec[1] != domain->avec[1] || avec[2] != domain->avec[2]) errflag = 1; @@ -627,22 +627,22 @@ void ReadData::command(int narg, char **arg) abc_origin[2] != domain->boxlo[2]) errflag = 1; if (errflag) - error->all(FLERR,"Read_data subsequent file ABC vectors must be same as first file"); + error->all(FLERR, "Read_data subsequent file ABC vectors must be same as first file"); if (shift[0] != 0.0 || shift[1] != 0.0 || shift[2] != 0.0) - error->all(FLERR,"Read_data subsequent file with ABC vectors cannot define shift"); + error->all(FLERR, "Read_data subsequent file with ABC vectors cannot define shift"); - // restricted triclinic - // tilt factors must match first data file + // restricted triclinic + // tilt factors must match first data file } else if (triclinic) { if (!domain->triclinic || domain->triclinic_general) - error->all(FLERR,"Read_data subsequent file cannot switch to restricted triclinic"); + error->all(FLERR, "Read_data subsequent file cannot switch to restricted triclinic"); if (xy != domain->xy || xz != domain->xz || yz != domain->yz) - error->all(FLERR,"Read_data subsequent file tilt factors must be same as first file"); + error->all(FLERR, "Read_data subsequent file tilt factors must be same as first file"); } else { if (domain->triclinic) - error->all(FLERR,"Read_data subsequent file cannot switch to orthogonal"); + error->all(FLERR, "Read_data subsequent file cannot switch to orthogonal"); } double oldboxlo[3] = {domain->boxlo[0], domain->boxlo[1], domain->boxlo[2]}; @@ -1068,8 +1068,7 @@ void ReadData::command(int narg, char **arg) // if general triclinic, perform general to restricted rotation operation // on any quantities read from data file which require it - if (triclinic_general) - atom->avec->read_data_general_to_restricted(nlocal_previous, atom->nlocal); + if (triclinic_general) atom->avec->read_data_general_to_restricted(nlocal_previous, atom->nlocal); // init per-atom fix/compute/variable values for created atoms @@ -1371,7 +1370,7 @@ void ReadData::header(int firstpass) else if (firstpass) atom->nimpropers += nimpropers; - // Atom class type settings are only set by first data file + // Atom class type settings are only set by first data file } else if (utils::strmatch(line, "^\\s*\\d+\\s+atom\\s+types\\s")) { ntypes = utils::inumeric(FLERR, words[0], false, lmp); @@ -1393,10 +1392,10 @@ void ReadData::header(int firstpass) nimpropertypes = utils::inumeric(FLERR, words[0], false, lmp); if (addflag == NONE) atom->nimpropertypes = nimpropertypes + extra_improper_types; - // these settings only used by first data file - // NOTE: these are now obsolete, we parse them to maintain backward compatibility - // the recommended way is to set them via command keywords in the input script - // if these flags are set both ways, the larger of the two values is used + // these settings only used by first data file + // NOTE: these are now obsolete, we parse them to maintain backward compatibility + // the recommended way is to set them via command keywords in the input script + // if these flags are set both ways, the larger of the two values is used } else if (strstr(line, "extra bond per atom")) { if (addflag == NONE) extra_flag_value = utils::inumeric(FLERR, words[0], false, lmp); @@ -1414,8 +1413,8 @@ void ReadData::header(int firstpass) if (addflag == NONE) extra_flag_value = utils::inumeric(FLERR, words[0], false, lmp); force->special_extra = MAX(force->special_extra, extra_flag_value); - // local copy of box info - // so can treat differently for first vs subsequent data files + // local copy of box info + // so can treat differently for first vs subsequent data files } else if (utils::strmatch(line, "^\\s*\\f+\\s+\\f+\\s+xlo\\s+xhi\\s")) { xloxhi_flag = 1; @@ -1527,8 +1526,8 @@ void ReadData::atoms() if (eof) error->all(FLERR, "Unexpected end of data file"); if (tlabelflag && !lmap->is_complete(Atom::ATOM)) error->all(FLERR, "Label map is incomplete: all types must be assigned a unique type label"); - atom->data_atoms(nchunk, buffer, id_offset, mol_offset, toffset, - shiftflag, shift, tlabelflag, lmap->lmap2lmap.atom, triclinic_general); + atom->data_atoms(nchunk, buffer, id_offset, mol_offset, toffset, shiftflag, shift, tlabelflag, + lmap->lmap2lmap.atom, triclinic_general); nread += nchunk; } diff --git a/src/read_restart.cpp b/src/read_restart.cpp index ee46d57fc3..ef30bc81d7 100644 --- a/src/read_restart.cpp +++ b/src/read_restart.cpp @@ -527,7 +527,8 @@ std::string ReadRestart::file_search(const std::string &inpfile) loc = pattern.find('*'); if (loc != std::string::npos) { // the regex matcher in utils::strmatch() only checks the first 256 characters. - if (loc > 256) + // a 64-bit integer timestep will consume 20 characters, so 236 chars is the cutoff. + if (loc > 236) error->one(FLERR, "Filename part before '*' is too long to find restart with largest step"); // convert pattern to equivalent regexp @@ -538,7 +539,7 @@ std::string ReadRestart::file_search(const std::string &inpfile) for (const auto &candidate : platform::list_directory(dirname)) { if (utils::strmatch(candidate,pattern)) { - bigint num = ATOBIGINT(utils::strfind(candidate.substr(loc),"\\d+").c_str()); + auto num = (bigint) std::stoll(utils::strfind(candidate.substr(loc),"\\d+")); if (num > maxnum) maxnum = num; } } diff --git a/src/reader_native.cpp b/src/reader_native.cpp index 4dac65e3cb..f252e01644 100644 --- a/src/reader_native.cpp +++ b/src/reader_native.cpp @@ -254,7 +254,7 @@ bigint ReaderNative::read_header(double box[3][3], int &boxinfo, int &triclinic, triclinic = 0; box[0][2] = box[1][2] = box[2][2] = 0.0; read_lines(1); - if (line[strlen("ITEM: BOX BOUNDS ")] == 'x') triclinic = 1; + if (utils::strmatch(line,"ITEM: BOX BOUNDS.*xy\\s+xz\\s+yz")) triclinic = 1; try { read_lines(1); @@ -484,7 +484,7 @@ void ReaderNative::read_atoms(int n, int nfield, double **fields) // convert selected fields to floats for (int m = 0; m < nfield; m++) - fields[i][m] = atof(words[fieldindex[m]].c_str()); + fields[i][m] = std::stod(words[fieldindex[m]]); } } } diff --git a/src/reader_xyz.cpp b/src/reader_xyz.cpp index 22b6df0cb9..385682533c 100644 --- a/src/reader_xyz.cpp +++ b/src/reader_xyz.cpp @@ -182,7 +182,7 @@ void ReaderXYZ::read_atoms(int n, int nfield, double **fields) // XXX: we could insert an element2type translation here // XXX: for now we flag unrecognized types as type 0, // XXX: which should trigger an error, if LAMMPS uses it. - mytype = atoi(line); + mytype = std::stoi(line); for (m = 0; m < nfield; m++) { switch (fieldindex[m]) { diff --git a/src/region_cone.cpp b/src/region_cone.cpp index c0cb2dbee6..dc37eeefe3 100644 --- a/src/region_cone.cpp +++ b/src/region_cone.cpp @@ -779,8 +779,7 @@ void RegCone::shape_update() if (lostyle == VARIABLE) lo = input->variable->compute_equal(lovar); if (histyle == VARIABLE) hi = input->variable->compute_equal(hivar); - if (radiuslo == 0.0 && radiushi == 0.0) - error->all(FLERR, "dtion in region gave bad value"); + if (radiuslo == 0.0 && radiushi == 0.0) error->all(FLERR, "dtion in region gave bad value"); if (axis == 'x') { if (c1style == VARIABLE) c1 *= yscale; diff --git a/src/thermo.cpp b/src/thermo.cpp index 892de240b2..b2e4f7f934 100644 --- a/src/thermo.cpp +++ b/src/thermo.cpp @@ -575,9 +575,9 @@ void Thermo::modify_params(int narg, char **arg) iarg += 2; - } else if (strcmp(arg[iarg],"triclinic/general") == 0) { - if (iarg + 2 > narg) error->all(FLERR,"Illegal thermo_modify command"); - triclinic_general = utils::logical(FLERR,arg[iarg+1],false,lmp); + } else if (strcmp(arg[iarg], "triclinic/general") == 0) { + if (iarg + 2 > narg) error->all(FLERR, "Illegal thermo_modify command"); + triclinic_general = utils::logical(FLERR, arg[iarg + 1], false, lmp); if (triclinic_general && !domain->triclinic_general) error->all(FLERR,"Thermo_modify triclinic/general cannot be used " "if simulation box is not general triclinic"); @@ -987,31 +987,35 @@ void Thermo::parse_fields(const std::string &str) } else if (word == "pxx") { if (triclinic_general) addfield("Pxx", &Thermo::compute_pxx_triclinic_general, FLOAT); - else addfield("Pxx", &Thermo::compute_pxx, FLOAT); + else + addfield("Pxx", &Thermo::compute_pxx, FLOAT); index_press_vector = add_compute(id_press, VECTOR); } else if (word == "pyy") { if (triclinic_general) addfield("Pyy", &Thermo::compute_pyy_triclinic_general, FLOAT); - else addfield("Pyy", &Thermo::compute_pyy, FLOAT); + else + addfield("Pyy", &Thermo::compute_pyy, FLOAT); index_press_vector = add_compute(id_press, VECTOR); } else if (word == "pzz") { if (triclinic_general) addfield("Pzz", &Thermo::compute_pzz_triclinic_general, FLOAT); - else addfield("Pzz", &Thermo::compute_pzz, FLOAT); + else + addfield("Pzz", &Thermo::compute_pzz, FLOAT); index_press_vector = add_compute(id_press, VECTOR); } else if (word == "pxy") { if (triclinic_general) addfield("Pxy", &Thermo::compute_pxy_triclinic_general, FLOAT); - else addfield("Pxy", &Thermo::compute_pxy, FLOAT); + else + addfield("Pxy", &Thermo::compute_pxy, FLOAT); index_press_vector = add_compute(id_press, VECTOR); } else if (word == "pxz") { if (triclinic_general) addfield("Pxz", &Thermo::compute_pxz_triclinic_general, FLOAT); - else addfield("Pxz", &Thermo::compute_pxz, FLOAT); + else + addfield("Pxz", &Thermo::compute_pxz, FLOAT); index_press_vector = add_compute(id_press, VECTOR); } else if (word == "pyz") { - if (triclinic_general) - addfield("Pyz", &Thermo::compute_pyz_triclinic_general, FLOAT); + if (triclinic_general) addfield("Pyz", &Thermo::compute_pyz_triclinic_general, FLOAT); addfield("Pyz", &Thermo::compute_pyz, FLOAT); index_press_vector = add_compute(id_press, VECTOR); @@ -1486,33 +1490,45 @@ int Thermo::evaluate_keyword(const std::string &word, double *answer) else if (word == "pxx") { check_press_vector(word); - if (triclinic_general) compute_pxx_triclinic_general(); - else compute_pxx(); + if (triclinic_general) + compute_pxx_triclinic_general(); + else + compute_pxx(); } else if (word == "pyy") { check_press_vector(word); - if (triclinic_general) compute_pyy_triclinic_general(); - else compute_pyy(); + if (triclinic_general) + compute_pyy_triclinic_general(); + else + compute_pyy(); } else if (word == "pzz") { check_press_vector(word); - if (triclinic_general) compute_pzz_triclinic_general(); - else compute_pzz(); + if (triclinic_general) + compute_pzz_triclinic_general(); + else + compute_pzz(); } else if (word == "pxy") { check_press_vector(word); - if (triclinic_general) compute_pxy_triclinic_general(); - else compute_pxy(); + if (triclinic_general) + compute_pxy_triclinic_general(); + else + compute_pxy(); } else if (word == "pxz") { check_press_vector(word); - if (triclinic_general) compute_pxz_triclinic_general(); - else compute_pxz(); + if (triclinic_general) + compute_pxz_triclinic_general(); + else + compute_pxz(); } else if (word == "pyz") { check_press_vector(word); - if (triclinic_general) compute_pyz_triclinic_general(); - else compute_pyz(); + if (triclinic_general) + compute_pyz_triclinic_general(); + else + compute_pyz(); } else if (word == "bonds") { compute_bonds(); @@ -1603,8 +1619,10 @@ void Thermo::compute_fix() // if index exceeds variable vector length, use a zero value // this can be useful if vector length is not known a priori - if (fix->size_vector_variable && argindex1[ifield] > fix->size_vector) dvalue = 0.0; - else dvalue = fix->compute_vector(argindex1[ifield] - 1); + if (fix->size_vector_variable && argindex1[ifield] > fix->size_vector) + dvalue = 0.0; + else + dvalue = fix->compute_vector(argindex1[ifield] - 1); if (normflag) { if (fix->extvector == 0) return; @@ -1618,8 +1636,10 @@ void Thermo::compute_fix() // if index exceeds variable array rows, use a zero value // this can be useful if array size is not known a priori - if (fix->size_array_rows_variable && argindex1[ifield] > fix->size_array_rows) dvalue = 0.0; - else dvalue = fix->compute_array(argindex1[ifield] - 1, argindex2[ifield] - 1); + if (fix->size_array_rows_variable && argindex1[ifield] > fix->size_array_rows) + dvalue = 0.0; + else + dvalue = fix->compute_array(argindex1[ifield] - 1, argindex2[ifield] - 1); if (normflag && fix->extarray) dvalue /= natoms; } } @@ -1639,8 +1659,10 @@ void Thermo::compute_variable() else { double *varvec; int nvec = input->variable->compute_vector(variables[field2index[ifield]], &varvec); - if (iarg > nvec) dvalue = 0.0; - else dvalue = varvec[iarg - 1]; + if (iarg > nvec) + dvalue = 0.0; + else + dvalue = varvec[iarg - 1]; } } @@ -2087,79 +2109,106 @@ void Thermo::compute_yz() void Thermo::compute_avecx() { - if (!domain->triclinic) dvalue = domain->xprd; - else if (triclinic_general) dvalue = domain->avec[0]; - else dvalue = domain->xprd; + if (!domain->triclinic) + dvalue = domain->xprd; + else if (triclinic_general) + dvalue = domain->avec[0]; + else + dvalue = domain->xprd; } /* ---------------------------------------------------------------------- */ void Thermo::compute_avecy() { - if (!domain->triclinic) dvalue = 0.0; - else if (triclinic_general) dvalue = domain->avec[1]; - else dvalue = 0.0; + if (!domain->triclinic) + dvalue = 0.0; + else if (triclinic_general) + dvalue = domain->avec[1]; + else + dvalue = 0.0; } /* ---------------------------------------------------------------------- */ void Thermo::compute_avecz() { - if (!domain->triclinic) dvalue = 0.0; - else if (triclinic_general) dvalue = domain->avec[2]; - else dvalue = 0.0; + if (!domain->triclinic) + dvalue = 0.0; + else if (triclinic_general) + dvalue = domain->avec[2]; + else + dvalue = 0.0; } /* ---------------------------------------------------------------------- */ void Thermo::compute_bvecx() { - if (!domain->triclinic) dvalue = 0.0; - else if (triclinic_general) dvalue = domain->bvec[0]; - else dvalue = domain->xy; + if (!domain->triclinic) + dvalue = 0.0; + else if (triclinic_general) + dvalue = domain->bvec[0]; + else + dvalue = domain->xy; } /* ---------------------------------------------------------------------- */ void Thermo::compute_bvecy() { - if (!domain->triclinic) dvalue = domain->yprd; - else if (triclinic_general) dvalue = domain->bvec[1]; - else dvalue = domain->yprd; + if (!domain->triclinic) + dvalue = domain->yprd; + else if (triclinic_general) + dvalue = domain->bvec[1]; + else + dvalue = domain->yprd; } /* ---------------------------------------------------------------------- */ void Thermo::compute_bvecz() { - if (!domain->triclinic) dvalue = 0.0; - else if (triclinic_general) dvalue = domain->bvec[2]; - else dvalue = 0.0; + if (!domain->triclinic) + dvalue = 0.0; + else if (triclinic_general) + dvalue = domain->bvec[2]; + else + dvalue = 0.0; } /* ---------------------------------------------------------------------- */ void Thermo::compute_cvecx() { - if (!domain->triclinic) dvalue = 0.0; - else if (triclinic_general) dvalue = domain->cvec[0]; - else dvalue = domain->xz; + if (!domain->triclinic) + dvalue = 0.0; + else if (triclinic_general) + dvalue = domain->cvec[0]; + else + dvalue = domain->xz; } /* ---------------------------------------------------------------------- */ void Thermo::compute_cvecy() { - if (!domain->triclinic) dvalue = 0.0; - else if (triclinic_general) dvalue = domain->cvec[1]; - else dvalue = domain->yz; + if (!domain->triclinic) + dvalue = 0.0; + else if (triclinic_general) + dvalue = domain->cvec[1]; + else + dvalue = domain->yz; } /* ---------------------------------------------------------------------- */ void Thermo::compute_cvecz() { - if (!domain->triclinic) dvalue = domain->zprd; - else if (triclinic_general) dvalue = domain->cvec[2]; - else dvalue = domain->zprd; + if (!domain->triclinic) + dvalue = domain->zprd; + else if (triclinic_general) + dvalue = domain->cvec[2]; + else + dvalue = domain->zprd; } /* ---------------------------------------------------------------------- */ @@ -2309,9 +2358,9 @@ void Thermo::compute_pyz() void Thermo::compute_pxx_triclinic_general() { - double middle[3][3],final[3][3]; - MathExtra::times3(domain->rotate_r2g,press_tensor,middle); - MathExtra::times3(middle,domain->rotate_g2r,final); + double middle[3][3], final[3][3]; + MathExtra::times3(domain->rotate_r2g, press_tensor, middle); + MathExtra::times3(middle, domain->rotate_g2r, final); dvalue = final[0][0]; } @@ -2319,9 +2368,9 @@ void Thermo::compute_pxx_triclinic_general() void Thermo::compute_pyy_triclinic_general() { - double middle[3][3],final[3][3]; - MathExtra::times3(domain->rotate_r2g,press_tensor,middle); - MathExtra::times3(middle,domain->rotate_g2r,final); + double middle[3][3], final[3][3]; + MathExtra::times3(domain->rotate_r2g, press_tensor, middle); + MathExtra::times3(middle, domain->rotate_g2r, final); dvalue = final[1][1]; } @@ -2329,9 +2378,9 @@ void Thermo::compute_pyy_triclinic_general() void Thermo::compute_pzz_triclinic_general() { - double middle[3][3],final[3][3]; - MathExtra::times3(domain->rotate_r2g,press_tensor,middle); - MathExtra::times3(middle,domain->rotate_g2r,final); + double middle[3][3], final[3][3]; + MathExtra::times3(domain->rotate_r2g, press_tensor, middle); + MathExtra::times3(middle, domain->rotate_g2r, final); dvalue = final[2][2]; } @@ -2339,9 +2388,9 @@ void Thermo::compute_pzz_triclinic_general() void Thermo::compute_pxy_triclinic_general() { - double middle[3][3],final[3][3]; - MathExtra::times3(domain->rotate_r2g,press_tensor,middle); - MathExtra::times3(middle,domain->rotate_g2r,final); + double middle[3][3], final[3][3]; + MathExtra::times3(domain->rotate_r2g, press_tensor, middle); + MathExtra::times3(middle, domain->rotate_g2r, final); dvalue = final[0][1]; } @@ -2349,9 +2398,9 @@ void Thermo::compute_pxy_triclinic_general() void Thermo::compute_pxz_triclinic_general() { - double middle[3][3],final[3][3]; - MathExtra::times3(domain->rotate_r2g,press_tensor,middle); - MathExtra::times3(middle,domain->rotate_g2r,final); + double middle[3][3], final[3][3]; + MathExtra::times3(domain->rotate_r2g, press_tensor, middle); + MathExtra::times3(middle, domain->rotate_g2r, final); dvalue = final[0][2]; } @@ -2359,9 +2408,9 @@ void Thermo::compute_pxz_triclinic_general() void Thermo::compute_pyz_triclinic_general() { - double middle[3][3],final[3][3]; - MathExtra::times3(domain->rotate_r2g,press_tensor,middle); - MathExtra::times3(middle,domain->rotate_g2r,final); + double middle[3][3], final[3][3]; + MathExtra::times3(domain->rotate_r2g, press_tensor, middle); + MathExtra::times3(middle, domain->rotate_g2r, final); dvalue = final[1][2]; } diff --git a/src/tokenizer.cpp b/src/tokenizer.cpp index 6b87f0c421..a2b95ab89e 100644 --- a/src/tokenizer.cpp +++ b/src/tokenizer.cpp @@ -287,8 +287,20 @@ std::string ValueTokenizer::next_string() int ValueTokenizer::next_int() { std::string current = tokens.next(); - if (!utils::is_integer(current)) { throw InvalidIntegerException(current); } - return atoi(current.c_str()); + try { + std::size_t end; + auto val = std::stoi(current, &end); + // only partially converted + if (current.size() != end) { throw InvalidIntegerException(current); } + return val; + + // rethrow exceptions from std::stoi() + } catch (std::out_of_range const &) { + throw InvalidIntegerException(current); + } catch (std::invalid_argument const &) { + throw InvalidIntegerException(current); + } + return 0; } /*! Retrieve next token and convert to bigint @@ -297,8 +309,22 @@ int ValueTokenizer::next_int() bigint ValueTokenizer::next_bigint() { std::string current = tokens.next(); - if (!utils::is_integer(current)) { throw InvalidIntegerException(current); } - return ATOBIGINT(current.c_str()); + try { + std::size_t end; + auto val = std::stoll(current, &end, 10); + // only partially converted + if (current.size() != end) { throw InvalidIntegerException(current); } + // out of range + if ((val < (-MAXBIGINT - 1) || (val > MAXBIGINT))) { throw InvalidIntegerException(current); }; + return (bigint) val; + + // rethrow exceptions from std::stoll() + } catch (std::out_of_range const &) { + throw InvalidIntegerException(current); + } catch (std::invalid_argument const &) { + throw InvalidIntegerException(current); + } + return 0; } /*! Retrieve next token and convert to tagint @@ -307,8 +333,22 @@ bigint ValueTokenizer::next_bigint() tagint ValueTokenizer::next_tagint() { std::string current = tokens.next(); - if (!utils::is_integer(current)) { throw InvalidIntegerException(current); } - return ATOTAGINT(current.c_str()); + try { + std::size_t end; + auto val = std::stoll(current, &end, 10); + // only partially converted + if (current.size() != end) { throw InvalidIntegerException(current); } + // out of range + if ((val < (-MAXTAGINT - 1) || (val > MAXTAGINT))) { throw InvalidIntegerException(current); } + return (tagint) val; + + // rethrow exceptions from std::stoll() + } catch (std::out_of_range const &) { + throw InvalidIntegerException(current); + } catch (std::invalid_argument const &) { + throw InvalidIntegerException(current); + } + return 0; } /*! Retrieve next token and convert to double @@ -317,8 +357,19 @@ tagint ValueTokenizer::next_tagint() double ValueTokenizer::next_double() { std::string current = tokens.next(); - if (!utils::is_double(current)) { throw InvalidFloatException(current); } - return atof(current.c_str()); + try { + std::size_t end; + auto val = std::stod(current, &end); + // only partially converted + if (current.size() != end) { throw InvalidFloatException(current); } + return val; + // rethrow exceptions from std::stod() + } catch (std::out_of_range const &) { + throw InvalidFloatException(current); + } catch (std::invalid_argument const &) { + throw InvalidFloatException(current); + } + return 0.0; } /*! Skip over a given number of tokens diff --git a/src/universe.cpp b/src/universe.cpp index edd5b01031..fb07e26759 100644 --- a/src/universe.cpp +++ b/src/universe.cpp @@ -181,10 +181,10 @@ void Universe::add_world(char *str) if ((found == 0) || (found == (part.size() - 1))) { valid = false; } else if (found == std::string::npos) { - nper = atoi(part.c_str()); + nper = std::stoi(part); } else { - n = atoi(part.substr(0,found).c_str()); - nper = atoi(part.substr(found+1).c_str()); + n = std::stoi(part.substr(0,found)); + nper = std::stoi(part.substr(found+1)); } } @@ -193,8 +193,7 @@ void Universe::add_world(char *str) if (n < 1 || nper < 1) valid = false; if (!valid) - error->universe_all(FLERR,fmt::format("Invalid partition string '{}'", - str)); + error->universe_all(FLERR, fmt::format("Invalid partition string '{}'", str)); } else nper = nprocs; memory->grow(procs_per_world,nworlds+n,"universe:procs_per_world"); diff --git a/src/update.cpp b/src/update.cpp index c9e57f9621..9494629221 100644 --- a/src/update.cpp +++ b/src/update.cpp @@ -404,7 +404,6 @@ void Update::create_minimize(int narg, char **arg, int trysuffix) minimize_style = utils::strdup(arg[0]); minimize = nullptr; - int sflag; new_minimize(arg[0], narg - 1, &arg[1], trysuffix, sflag); diff --git a/src/utils.cpp b/src/utils.cpp index 87c46c47ab..27d89f1bf5 100644 --- a/src/utils.cpp +++ b/src/utils.cpp @@ -398,17 +398,23 @@ double utils::numeric(const char *file, int line, const std::string &str, bool d } double rv = 0; + auto msg = fmt::format("Floating point number {} in input script or data file is invalid", buf); try { - rv = stod(buf); + std::size_t endpos; + rv = std::stod(buf, &endpos); + if (buf.size() != endpos) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); + } } catch (std::invalid_argument const &) { - auto msg = fmt::format("Floating point number {} in input script or data file is invalid", buf); if (do_abort) lmp->error->one(file, line, msg); else lmp->error->all(file, line, msg); } catch (std::out_of_range const &) { - auto msg = - fmt::format("Floating point number {} in input script or data file is out of range", buf); + msg = fmt::format("Floating point number {} in input script or data file is out of range", buf); if (do_abort) lmp->error->one(file, line, msg); else @@ -459,10 +465,23 @@ int utils::inumeric(const char *file, int line, const std::string &str, bool do_ } int rv = 0; + auto msg = fmt::format("Integer {} in input script or data file is invalid", buf); try { - rv = stoi(buf); + std::size_t endpos; + rv = std::stoi(buf, &endpos); + if (buf.size() != endpos) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); + } + } catch (std::invalid_argument const &) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); } catch (std::out_of_range const &) { - auto msg = fmt::format("Integer {} in input script or data file is out of range", buf); + msg = fmt::format("Integer {} in input script or data file is out of range", buf); if (do_abort) lmp->error->one(file, line, msg); else @@ -514,9 +533,22 @@ bigint utils::bnumeric(const char *file, int line, const std::string &str, bool } long long rv = 0; + auto msg = fmt::format("Integer {} in input script or data file is invalid", buf); try { - rv = stoll(buf); - if (rv > MAXBIGINT) throw std::out_of_range("64-bit"); + std::size_t endpos; + rv = std::stoll(buf, &endpos); + if (buf.size() != endpos) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); + } + if ((rv < (-MAXBIGINT - 1) || (rv > MAXBIGINT))) throw std::out_of_range("bigint"); + } catch (std::invalid_argument const &) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); } catch (std::out_of_range const &) { auto msg = fmt::format("Integer {} in input script or data file is out of range", buf); if (do_abort) @@ -570,9 +602,22 @@ tagint utils::tnumeric(const char *file, int line, const std::string &str, bool } long long rv = 0; + auto msg = fmt::format("Integer {} in input script or data file is invalid", buf); try { - rv = stoll(buf); - if (rv > MAXTAGINT) throw std::out_of_range("64-bit"); + std::size_t endpos; + rv = std::stoll(buf, &endpos); + if (buf.size() != endpos) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); + } + if ((rv < (-MAXTAGINT - 1) || (rv > MAXTAGINT))) throw std::out_of_range("tagint"); + } catch (std::invalid_argument const &) { + if (do_abort) + lmp->error->one(file, line, msg); + else + lmp->error->all(file, line, msg); } catch (std::out_of_range const &) { auto msg = fmt::format("Integer {} in input script or data file is out of range", buf); if (do_abort) @@ -614,19 +659,19 @@ void utils::bounds(const char *file, int line, const std::string &str, found = str.find_first_of('*'); if (found == std::string::npos) { // contains no '*' - nlo = nhi = strtol(str.c_str(), nullptr, 10); + nlo = nhi = std::stol(str, nullptr, 10); } else if (str.size() == 1) { // is only '*' nlo = nmin; nhi = nmax; } else if (found == 0) { // is '*j' nlo = nmin; - nhi = strtol(str.substr(1).c_str(), nullptr, 10); + nhi = std::stol(str.substr(1), nullptr, 10); } else if (str.size() == found + 1) { // is 'i*' - nlo = strtol(str.c_str(), nullptr, 10); + nlo = std::stol(str, nullptr, 10); nhi = nmax; } else { // is 'i*j' - nlo = strtol(str.c_str(), nullptr, 10); - nhi = strtol(str.substr(found + 1).c_str(), nullptr, 10); + nlo = std::stol(str, nullptr, 10); + nhi = std::stol(str.substr(found + 1), nullptr, 10); } if (error) { @@ -1697,10 +1742,10 @@ double utils::timespec2seconds(const std::string ×pec) int utils::date2num(const std::string &date) { std::size_t found = date.find_first_not_of("0123456789 "); - int num = strtol(date.substr(0, found).c_str(), nullptr, 10); + int num = std::stol(date.substr(0, found), nullptr, 10); auto month = date.substr(found); found = month.find_first_of("0123456789 "); - num += strtol(month.substr(found).c_str(), nullptr, 10) * 10000; + num += std::stol(month.substr(found), nullptr, 10) * 10000; if (num < 1000000) num += 20000000; if (strmatch(month, "^Jan")) diff --git a/src/variable.cpp b/src/variable.cpp index f308ed6efc..12ae516983 100644 --- a/src/variable.cpp +++ b/src/variable.cpp @@ -46,6 +46,7 @@ #include #include #include +#include #include using namespace LAMMPS_NS; @@ -799,8 +800,8 @@ int Variable::next(int narg, char **arg) fclose(fp); if (strlen(buf) > 0) { - nextindex = atoi(buf); - break; + nextindex = std::stoi(buf); + break; } delay = (int) (1000000*random->uniform()); platform::usleep(delay); @@ -1013,6 +1014,8 @@ char *Variable::retrieve(const char *name) } else if (style[ivar] == EQUAL) { double answer = evaluate(data[ivar][0],nullptr,ivar); + // round to zero on underflow + if (fabs(answer) < std::numeric_limits::min()) answer = 0.0; delete[] data[ivar][1]; data[ivar][1] = utils::strdup(fmt::format("{:.15g}",answer)); str = data[ivar][1]; @@ -1115,9 +1118,15 @@ double Variable::compute_equal(int ivar) if (ifunc < 0) print_var_error(FLERR,fmt::format("cannot find python function {}",data[ivar][0]),ivar); python->invoke_function(ifunc,data[ivar][1]); - value = atof(data[ivar][1]); + try { + value = std::stod(data[ivar][1]); + } catch (std::exception &e) { + print_var_error(FLERR,"has an invalid value", ivar); + } } + // round to zero on underflow + if (fabs(value) < std::numeric_limits::min()) value = 0.0; eval_in_progress[ivar] = 0; return value; } @@ -1132,6 +1141,7 @@ double Variable::compute_equal(const std::string &str) { char *ptr = utils::strdup(str); double val = evaluate(ptr,nullptr,-1); + if (fabs(val) < std::numeric_limits::min()) val = 0.0; delete[] ptr; return val; } @@ -1461,9 +1471,9 @@ double Variable::evaluate(char *str, Tree **tree, int ivar) if (tree) { auto newtree = new Tree(); newtree->type = VALUE; - newtree->value = atof(number); + newtree->value = std::stod(number); treestack[ntreestack++] = newtree; - } else argstack[nargstack++] = atof(number); + } else argstack[nargstack++] = std::stod(number); delete[] number; @@ -2066,9 +2076,9 @@ double Variable::evaluate(char *str, Tree **tree, int ivar) if (tree) { auto newtree = new Tree(); newtree->type = VALUE; - newtree->value = atof(var); + newtree->value = std::stod(var); treestack[ntreestack++] = newtree; - } else argstack[nargstack++] = atof(var); + } else argstack[nargstack++] = std::stod(var); // vector from vector-style variable // evaluate the vector-style variable, put result in newtree @@ -3555,7 +3565,7 @@ tagint Variable::int_between_brackets(char *&ptr, int varallow) } else { varflag = 0; while (*ptr && *ptr != ']') { - if (!isdigit(*ptr)) + if (!(isdigit(*ptr) || (*ptr == '-') || (*ptr == '+') || (*ptr == '*') || (*ptr == '/'))) error->all(FLERR,"Non digit character between brackets in variable"); ptr++; } @@ -3566,21 +3576,23 @@ tagint Variable::int_between_brackets(char *&ptr, int varallow) *ptr = '\0'; - // evaluate index as floating point variable or as tagint via ATOTAGINT() + // evaluate index as floating point variable or as tagint via stoll() - if (varflag) { - char *id = start+2; - int ivar = find(id); - if (ivar < 0) - error->all(FLERR,"Invalid variable name {} in variable formula", id); - - char *var = retrieve(id); - if (var == nullptr) - error->all(FLERR,"Invalid variable evaluation for variable {} in variable formula", id); - index = static_cast (atof(var)); - - } else index = ATOTAGINT(start); + try { + if (varflag) { + char *id = start+2; + int ivar = find(id); + if (ivar < 0) + error->all(FLERR,"Invalid variable name {} in variable formula", id); + char *var = retrieve(id); + if (var == nullptr) + error->all(FLERR,"Invalid variable evaluation for variable {} in variable formula", id); + index = static_cast(std::stod(var)); + } else index = static_cast(std::stoll(start)); + } catch (std::exception &e) { + error->all(FLERR,"Illegal value in brackets: {}({})", e.what(), start); + } *ptr = ']'; if (index <= 0) @@ -4685,7 +4697,7 @@ int Variable::special_function(const std::string &word, char *contents, Tree **t // save value in tree or on argstack if (style[ivar] == SCALARFILE) { - double value = atof(data[ivar][0]); + double value = std::stod(data[ivar][0]); int done = reader[ivar]->read_scalar(data[ivar][0]); if (done) remove(ivar); @@ -5278,7 +5290,7 @@ double Variable::evaluate_boolean(char *str) onechar = str[i]; str[i] = '\0'; - argstack[nargstack].value = atof(&str[istart]); + argstack[nargstack].value = std::stod(&str[istart]); str[i] = onechar; argstack[nargstack++].flag = 0; diff --git a/src/write_coeff.cpp b/src/write_coeff.cpp index 11987d4b37..5a7177f62f 100644 --- a/src/write_coeff.cpp +++ b/src/write_coeff.cpp @@ -152,7 +152,7 @@ void WriteCoeff::command(int narg, char **arg) } // parse type number and skip over it - int type = atoi(str); + int type = std::stoi(str); char *p = str; while ((*p != '\0') && (*p == ' ')) ++p; while ((*p != '\0') && isdigit(*p)) ++p; diff --git a/tools/lammps-gui/CMakeLists.txt b/tools/lammps-gui/CMakeLists.txt index 435516a521..7abf8f1e67 100644 --- a/tools/lammps-gui/CMakeLists.txt +++ b/tools/lammps-gui/CMakeLists.txt @@ -1,6 +1,6 @@ cmake_minimum_required(VERSION 3.16) -project(lammps-gui VERSION 1.6.3 LANGUAGES CXX) +project(lammps-gui VERSION 1.6.5 LANGUAGES CXX) set(CMAKE_AUTOUIC ON) set(CMAKE_AUTOMOC ON) @@ -22,7 +22,13 @@ function(check_omp_h_include) set(CMAKE_REQUIRED_INCLUDES ${OpenMP_CXX_INCLUDE_DIRS}) set(CMAKE_REQUIRED_LINK_OPTIONS ${OpenMP_CXX_FLAGS}) set(CMAKE_REQUIRED_LIBRARIES ${OpenMP_CXX_LIBRARIES}) - check_include_file_cxx(omp.h _have_omp_h) + # there are all kinds of problems with finding omp.h + # for Clang and derived compilers so we pretend it is there. + if(CMAKE_CXX_COMPILER_ID STREQUAL "Clang") + set(_have_omp_h TRUE) + else() + check_include_file_cxx(omp.h _have_omp_h) + endif() else() set(_have_omp_h FALSE) endif() @@ -199,7 +205,7 @@ if(APPLE) MACOSX_BUNDLE_BUNDLE_VERSION ${PROJECT_VERSION} MACOSX_BUNDLE_SHORT_VERSION_STRING ${PROJECT_VERSION_MAJOR}.${PROJECT_VERSION_MINOR} MACOSX_BUNDLE_ICON_FILE lammps.icns - MACOSX_BUNDLE_COPYRIGHT "(c) 2003 - 2023, The LAMMPS Developers" + MACOSX_BUNDLE_COPYRIGHT "(c) 2003 - 2024, The LAMMPS Developers" MACOSX_BUNDLE TRUE ) # additional targets to populate the bundle tree and create the .dmg image file diff --git a/tools/lammps-gui/TODO.md b/tools/lammps-gui/TODO.md index d0b78b5591..bef30c1309 100644 --- a/tools/lammps-gui/TODO.md +++ b/tools/lammps-gui/TODO.md @@ -2,8 +2,9 @@ LAMMPS-GUI TODO list: # Short term goals (v1.x) -- add a "file viewer", also view file option in editor context menu if word under cursor is a file -- add "export to YAML" to chart viewer. +- figure out how widgets can be resized to fraction of available screen size. +- figure out stacking order of frames and whether it can be more flexible +- bundle LAMMPS tutorial input files - implement indenting regions for (nested) loops? - implement data file manager GUI with the following features: diff --git a/tools/lammps-gui/chartviewer.cpp b/tools/lammps-gui/chartviewer.cpp index d176dd98b7..2bbe029ebf 100644 --- a/tools/lammps-gui/chartviewer.cpp +++ b/tools/lammps-gui/chartviewer.cpp @@ -17,6 +17,10 @@ #include #include +#include +#include +#include +#include #include #include #include @@ -30,13 +34,17 @@ #include #include #include +#include +#include #include using namespace QtCharts; ChartWindow::ChartWindow(const QString &_filename, QWidget *parent) : - QWidget(parent), menu(new QMenuBar), file(new QMenu("&File")), filename(_filename) + QWidget(parent), menu(new QMenuBar), file(new QMenu("&File")), saveAsAct(nullptr), + exportCsvAct(nullptr), exportDatAct(nullptr), exportYamlAct(nullptr), closeAct(nullptr), + stopAct(nullptr), quitAct(nullptr), filename(_filename) { auto *top = new QHBoxLayout; menu->addMenu(file); @@ -61,6 +69,8 @@ ChartWindow::ChartWindow(const QString &_filename, QWidget *parent) : exportCsvAct->setIcon(QIcon(":/icons/application-calc.png")); exportDatAct = file->addAction("Export data to &Gnuplot...", this, &ChartWindow::exportDat); exportDatAct->setIcon(QIcon(":/icons/application-plot.png")); + exportYamlAct = file->addAction("Export data to &YAML...", this, &ChartWindow::exportYaml); + exportYamlAct->setIcon(QIcon(":/icons/yaml-file-icon.png")); file->addSeparator(); stopAct = file->addAction("Stop &Run", this, &ChartWindow::stop_run); stopAct->setIcon(QIcon(":/icons/process-stop.png")); @@ -227,6 +237,40 @@ void ChartWindow::exportCsv() } } } +void ChartWindow::exportYaml() +{ + if (charts.empty()) return; + QString defaultname = filename + ".yaml"; + if (filename.isEmpty()) defaultname = "lammpsdata.yaml"; + QString fileName = QFileDialog::getSaveFileName(this, "Save Chart as YAML data", defaultname, + "Image Files (*.yaml, *.yml)"); + if (!fileName.isEmpty()) { + QFile file(fileName); + if (file.open(QIODevice::WriteOnly | QIODevice::Text)) { + QTextStream out(&file); + out.setRealNumberPrecision(8); + out << "---\n"; + + out << "keywords: ['Step'"; + for (auto &c : charts) + out << ", " << c->get_title(); + out << "]\n"; + + out << "data: \n"; + int lines = charts[0]->get_count(); + for (int i = 0; i < lines; ++i) { + // timestep + out << " - [" << charts[0]->get_step(i); + // data + for (auto &c : charts) + out << ", " << c->get_data(i); + out << "]\n"; + } + out << "...\n"; + file.close(); + } + } +} void ChartWindow::change_chart(int) { diff --git a/tools/lammps-gui/chartviewer.h b/tools/lammps-gui/chartviewer.h index da0468eaf8..37a26c8b57 100644 --- a/tools/lammps-gui/chartviewer.h +++ b/tools/lammps-gui/chartviewer.h @@ -20,6 +20,8 @@ #include class QAction; +class QCloseEvent; +class QEvent; class QMenuBar; class QMenu; namespace QtCharts { @@ -50,6 +52,7 @@ private slots: void saveAs(); void exportDat(); void exportCsv(); + void exportYaml(); void change_chart(int index); @@ -61,7 +64,7 @@ private: QMenuBar *menu; QMenu *file; QComboBox *columns; - QAction *saveAsAct, *exportCsvAct, *exportDatAct; + QAction *saveAsAct, *exportCsvAct, *exportDatAct, *exportYamlAct; QAction *closeAct, *stopAct, *quitAct; QString filename; @@ -70,10 +73,10 @@ private: /* -------------------------------------------------------------------- */ -#include #include #include #include +class QChart; namespace QtCharts { class ChartViewer : public QChartView { diff --git a/tools/lammps-gui/codeeditor.cpp b/tools/lammps-gui/codeeditor.cpp index 2d349e28ab..36a811b92d 100644 --- a/tools/lammps-gui/codeeditor.cpp +++ b/tools/lammps-gui/codeeditor.cpp @@ -24,21 +24,27 @@ #include #include #include +#include #include #include #include #include #include +#include #include #include #include #include #include -#include #include +#include #include #include +#include +#include +#include +#include #include #include @@ -775,6 +781,11 @@ void CodeEditor::contextMenuEvent(QContextMenuEvent *event) action->setData(QString()); connect(action, &QAction::triggered, this, &CodeEditor::open_help); + action = menu->addAction(QString("LAMMPS Tutorial")); + action->setIcon(QIcon(":/icons/help-tutorial.png")); + action->setData(QString("https://lammpstutorials.github.io/")); + connect(action, &QAction::triggered, this, &CodeEditor::open_url); + menu->exec(event->globalPos()); delete menu; } @@ -1208,6 +1219,12 @@ void CodeEditor::open_help() QUrl(QString("https://docs.lammps.org/%1").arg(act->data().toString()))); } +void CodeEditor::open_url() +{ + auto *act = qobject_cast(sender()); + QDesktopServices::openUrl(QUrl(act->data().toString())); +} + void CodeEditor::view_file() { auto *act = qobject_cast(sender()); diff --git a/tools/lammps-gui/codeeditor.h b/tools/lammps-gui/codeeditor.h index b36f0d8de0..873b13e302 100644 --- a/tools/lammps-gui/codeeditor.h +++ b/tools/lammps-gui/codeeditor.h @@ -14,15 +14,23 @@ #ifndef CODEEDITOR_H #define CODEEDITOR_H -#include #include #include #include #include class QCompleter; -class QStringListModel; +class QContextMenuEvent; +class QDragEnterEvent; +class QDropEvent; +class QFont; +class QKeyEvent; +class QMimeData; +class QPaintEvent; +class QRect; +class QResizeEvent; class QShortcut; +class QWidget; class CodeEditor : public QPlainTextEdit { Q_OBJECT @@ -78,6 +86,7 @@ private slots: void get_help(); void find_help(QString &page, QString &help); void open_help(); + void open_url(); void view_file(); void reformatCurrentLine(); void runCompletion(); diff --git a/tools/lammps-gui/fileviewer.cpp b/tools/lammps-gui/fileviewer.cpp index f93039d809..69bda7d91b 100644 --- a/tools/lammps-gui/fileviewer.cpp +++ b/tools/lammps-gui/fileviewer.cpp @@ -15,15 +15,16 @@ #include "lammpsgui.h" -#include #include -#include #include +#include #include #include +#include #include #include #include +#include #include #include @@ -39,24 +40,65 @@ FileViewer::FileViewer(const QString &_filename, QWidget *parent) : // open and read file. Set editor to read-only. QFile file(fileName); - if (file.open(QIODevice::Text | QIODevice::ReadOnly)) { - QTextStream in(&file); - QString content = in.readAll(); - file.close(); + QFileInfo finfo(file); + QString command; + QString content; + QProcess decomp; + QStringList args = {"-cdf", fileName}; + bool compressed = false; - QFont text_font; - QSettings settings; - text_font.fromString(settings.value("textfont", text_font.toString()).toString()); - document()->setDefaultFont(text_font); - - document()->setPlainText(content); - moveCursor(QTextCursor::Start, QTextCursor::MoveAnchor); - setReadOnly(true); - setLineWrapMode(NoWrap); - setMinimumSize(800, 500); - setWindowIcon(QIcon(":/icons/lammps-icon-128x128.png")); - setWindowTitle("LAMMPS-GUI - Viewer - " + fileName); + // match suffix with decompression program + if (finfo.suffix() == "gz") { + command = "gzip"; + compressed = true; + } else if (finfo.suffix() == "bz2") { + command = "bzip2"; + compressed = true; + } else if (finfo.suffix() == "zst") { + command = "zstd"; + compressed = true; + } else if (finfo.suffix() == "xz") { + command = "xz"; + compressed = true; + } else if (finfo.suffix() == "lzma") { + command = "xz"; + args.insert(1, "--format=lzma"); + compressed = true; + } else if (finfo.suffix() == "lz4") { + command = "lz4"; + compressed = true; } + + // read compressed file from pipe + if (compressed) { + decomp.start(command, args, QIODevice::ReadOnly); + if (decomp.waitForStarted()) { + while (decomp.waitForReadyRead()) + content += decomp.readAll(); + } else { + content = "\nCould not open compressed file %1 with decompression program %2\n"; + content = content.arg(fileName).arg(command); + } + decomp.close(); + } else if (file.open(QIODevice::Text | QIODevice::ReadOnly)) { + // read plain text + QTextStream in(&file); + content = in.readAll(); + file.close(); + } + + QFont text_font; + QSettings settings; + text_font.fromString(settings.value("textfont", text_font.toString()).toString()); + document()->setDefaultFont(text_font); + + document()->setPlainText(content); + moveCursor(QTextCursor::Start, QTextCursor::MoveAnchor); + setReadOnly(true); + setLineWrapMode(NoWrap); + setMinimumSize(800, 500); + setWindowIcon(QIcon(":/icons/lammps-icon-128x128.png")); + setWindowTitle("LAMMPS-GUI - Viewer - " + fileName); } void FileViewer::quit() diff --git a/tools/lammps-gui/icons/application-yaml.png b/tools/lammps-gui/icons/application-yaml.png new file mode 100644 index 0000000000..d1457147d4 Binary files /dev/null and b/tools/lammps-gui/icons/application-yaml.png differ diff --git a/tools/lammps-gui/icons/help-browser.png b/tools/lammps-gui/icons/help-browser.png index f3edf2c430..f61fca4573 100644 Binary files a/tools/lammps-gui/icons/help-browser.png and b/tools/lammps-gui/icons/help-browser.png differ diff --git a/tools/lammps-gui/icons/help-tutorial.png b/tools/lammps-gui/icons/help-tutorial.png new file mode 100644 index 0000000000..94df016852 Binary files /dev/null and b/tools/lammps-gui/icons/help-tutorial.png differ diff --git a/tools/lammps-gui/icons/move-recenter.png b/tools/lammps-gui/icons/move-recenter.png new file mode 100644 index 0000000000..6de5429306 Binary files /dev/null and b/tools/lammps-gui/icons/move-recenter.png differ diff --git a/tools/lammps-gui/imageviewer.cpp b/tools/lammps-gui/imageviewer.cpp index 388f4251ae..0894855a02 100644 --- a/tools/lammps-gui/imageviewer.cpp +++ b/tools/lammps-gui/imageviewer.cpp @@ -20,27 +20,27 @@ #include #include #include +#include #include +#include #include +#include #include #include #include #include -#include +#include #include -#include #include -#include +#include #include -#include #include #include #include +#include #include -#include #include -#include -#include +#include #include @@ -152,6 +152,7 @@ ImageViewer::ImageViewer(const QString &fileName, LammpsWrapper *_lammps, QWidge vdwfactor = 0.5; auto pix = QPixmap(":/icons/emblem-photos.png"); + xcenter = ycenter = zcenter = 0.5; auto *renderstatus = new QLabel(QString()); renderstatus->setPixmap(pix.scaled(22, 22, Qt::KeepAspectRatio)); @@ -211,6 +212,8 @@ ImageViewer::ImageViewer(const QString &fileName, LammpsWrapper *_lammps, QWidge rotup->setToolTip("Rotate up by 15 degrees"); auto *rotdown = new QPushButton(QIcon(":/icons/gtk-go-down.png"), ""); rotdown->setToolTip("Rotate down by 15 degrees"); + auto *recenter = new QPushButton(QIcon(":/icons/move-recenter.png"), ""); + recenter->setToolTip("Recenter on group"); auto *reset = new QPushButton(QIcon(":/icons/gtk-zoom-fit.png"), ""); reset->setToolTip("Reset view to defaults"); auto *combo = new QComboBox; @@ -245,6 +248,7 @@ ImageViewer::ImageViewer(const QString &fileName, LammpsWrapper *_lammps, QWidge menuLayout->addWidget(rotright); menuLayout->addWidget(rotup); menuLayout->addWidget(rotdown); + menuLayout->addWidget(recenter); menuLayout->addWidget(reset); menuLayout->addWidget(new QLabel(" Group: ")); menuLayout->addWidget(combo); @@ -260,6 +264,7 @@ ImageViewer::ImageViewer(const QString &fileName, LammpsWrapper *_lammps, QWidge connect(rotright, &QPushButton::released, this, &ImageViewer::do_rot_right); connect(rotup, &QPushButton::released, this, &ImageViewer::do_rot_up); connect(rotdown, &QPushButton::released, this, &ImageViewer::do_rot_down); + connect(recenter, &QPushButton::released, this, &ImageViewer::do_recenter); connect(reset, &QPushButton::released, this, &ImageViewer::reset_view); connect(combo, SIGNAL(currentIndexChanged(int)), this, SLOT(change_group(int))); @@ -301,6 +306,7 @@ void ImageViewer::reset_view() showaxes = settings.value("axes", false).toBool(); usessao = settings.value("ssao", false).toBool(); antialias = settings.value("antialias", false).toBool(); + xcenter = ycenter = zcenter = 0.5; settings.endGroup(); // reset state of checkable push buttons and combo box (if accessible) @@ -421,6 +427,24 @@ void ImageViewer::do_rot_up() createImage(); } +void ImageViewer::do_recenter() +{ + QString commands = QString("variable LAMMPSGUI_CX delete\n" + "variable LAMMPSGUI_CY delete\n" + "variable LAMMPSGUI_CZ delete\n" + "variable LAMMPSGUI_CX equal (xcm(%1,x)-xlo)/lx\n" + "variable LAMMPSGUI_CY equal (xcm(%1,y)-ylo)/ly\n" + "variable LAMMPSGUI_CZ equal (xcm(%1,z)-zlo)/lz\n").arg(group); + lammps->commands_string(commands.toLocal8Bit()); + xcenter = lammps->extract_variable("LAMMPSGUI_CX"); + ycenter = lammps->extract_variable("LAMMPSGUI_CZ"); + zcenter = lammps->extract_variable("LAMMPSGUI_CZ"); + lammps->commands_string("variable LAMMPSGUI_CX delete\n" + "variable LAMMPSGUI_CY delete\n" + "variable LAMMPSGUI_CZ delete\n"); + createImage(); +} + void ImageViewer::cmd_to_clipboard() { auto words = last_dump_cmd.split(" "); @@ -534,6 +558,7 @@ void ImageViewer::createImage() else dumpcmd += " axes no 0.0 0.0"; + dumpcmd += QString(" center s %1 %2 %3").arg(xcenter).arg(ycenter).arg(zcenter); dumpcmd += " modify boxcolor " + settings.value("boxcolor", "yellow").toString(); dumpcmd += " backcolor " + settings.value("background", "black").toString(); if (useelements) dumpcmd += blank + elements + blank + adiams + blank; diff --git a/tools/lammps-gui/imageviewer.h b/tools/lammps-gui/imageviewer.h index 9612a0e5d9..f9b935640e 100644 --- a/tools/lammps-gui/imageviewer.h +++ b/tools/lammps-gui/imageviewer.h @@ -55,6 +55,7 @@ private slots: void do_rot_right(); void do_rot_up(); void do_rot_down(); + void do_recenter(); void cmd_to_clipboard(); void change_group(int); @@ -90,6 +91,7 @@ private: int xsize, ysize; int hrot, vrot; double zoom, vdwfactor; + double xcenter, ycenter, zcenter; bool showbox, showaxes, antialias, usessao, useelements, usediameter, usesigma; }; #endif diff --git a/tools/lammps-gui/lammpsgui.cpp b/tools/lammps-gui/lammpsgui.cpp index 99fc75d869..d42d59bca2 100644 --- a/tools/lammps-gui/lammpsgui.cpp +++ b/tools/lammps-gui/lammpsgui.cpp @@ -36,8 +36,6 @@ #include #include #include -#include -#include #include #include #include @@ -46,15 +44,15 @@ #include #include #include -#include #include #include +#include +#include #include #include #if defined(_OPENMP) -#include #include #endif @@ -210,6 +208,7 @@ LammpsGui::LammpsGui(QWidget *parent, const char *filename) : connect(ui->action_Help, &QAction::triggered, this, &LammpsGui::help); connect(ui->actionLAMMPS_GUI_Howto, &QAction::triggered, this, &LammpsGui::howto); connect(ui->actionLAMMPS_Manual, &QAction::triggered, this, &LammpsGui::manual); + connect(ui->actionLAMMPS_Tutorial, &QAction::triggered, this, &LammpsGui::tutorial); connect(ui->actionPreferences, &QAction::triggered, this, &LammpsGui::preferences); connect(ui->actionDefaults, &QAction::triggered, this, &LammpsGui::defaults); connect(ui->actionView_in_OVITO, &QAction::triggered, this, &LammpsGui::start_exe); @@ -1332,6 +1331,11 @@ void LammpsGui::manual() QDesktopServices::openUrl(QUrl("https://docs.lammps.org/")); } +void LammpsGui::tutorial() +{ + QDesktopServices::openUrl(QUrl("https://lammpstutorials.github.io/")); +} + void LammpsGui::howto() { QDesktopServices::openUrl(QUrl("https://docs.lammps.org/Howto_lammps_gui.html")); diff --git a/tools/lammps-gui/lammpsgui.h b/tools/lammps-gui/lammpsgui.h index f41266b485..29372efc1c 100644 --- a/tools/lammps-gui/lammpsgui.h +++ b/tools/lammps-gui/lammpsgui.h @@ -21,6 +21,7 @@ #include #include #include +#include #include #include "lammpswrapper.h" @@ -99,6 +100,7 @@ private slots: void about(); void help(); void manual(); + void tutorial(); void howto(); void logupdate(); void modified(); diff --git a/tools/lammps-gui/lammpsgui.qrc b/tools/lammps-gui/lammpsgui.qrc index b4bf6b578d..eb7c1fbc19 100644 --- a/tools/lammps-gui/lammpsgui.qrc +++ b/tools/lammps-gui/lammpsgui.qrc @@ -9,6 +9,7 @@ icons/application-calc.png icons/application-exit.png icons/application-plot.png + icons/application-yaml.png icons/axes-img.png icons/document-new.png icons/document-open-recent.png @@ -40,9 +41,11 @@ icons/help-about.png icons/help-browser.png icons/help-faq.png + icons/help-tutorial.png icons/image-x-generic.png icons/media-playback-start-2.png icons/media-playlist-repeat.png + icons/move-recenter.png icons/object-rotate-left.png icons/object-rotate-right.png icons/ovito.png diff --git a/tools/lammps-gui/lammpsgui.ui b/tools/lammps-gui/lammpsgui.ui index d0354feb6c..2607e715b9 100644 --- a/tools/lammps-gui/lammpsgui.ui +++ b/tools/lammps-gui/lammpsgui.ui @@ -86,6 +86,7 @@ + @@ -313,6 +314,17 @@ Ctrl+Shift+M + + + + + + LAMMPS &Tutorial + + + Ctrl+Shift+T + + @@ -454,7 +466,7 @@ - LAMMPS-GUI Howto + LAMMPS-&GUI Howto Ctrl+Shift+G diff --git a/tools/lammps-gui/lammpswrapper.cpp b/tools/lammps-gui/lammpswrapper.cpp index 1ccb225a0a..ed2bde1c9f 100644 --- a/tools/lammps-gui/lammpswrapper.cpp +++ b/tools/lammps-gui/lammpswrapper.cpp @@ -102,6 +102,26 @@ void *LammpsWrapper::extract_atom(const char *keyword) return val; } +// note: equal style and compatible variables only +double LammpsWrapper::extract_variable(const char *keyword) +{ + void *ptr = nullptr; + if (lammps_handle) { +#if defined(LAMMPS_GUI_USE_PLUGIN) + ptr = ((liblammpsplugin_t *)plugin_handle)->extract_variable(lammps_handle, keyword, nullptr); +#else + ptr = lammps_extract_variable(lammps_handle, keyword, nullptr); +#endif + } + double val = *((double *)ptr); +#if defined(LAMMPS_GUI_USE_PLUGIN) + ptr = ((liblammpsplugin_t *)plugin_handle)->free(ptr); +#else + lammps_free(ptr); +#endif + return val; +} + int LammpsWrapper::id_count(const char *keyword) { int val = 0; diff --git a/tools/lammps-gui/lammpswrapper.h b/tools/lammps-gui/lammpswrapper.h index 8719ef4491..9157bf77b5 100644 --- a/tools/lammps-gui/lammpswrapper.h +++ b/tools/lammps-gui/lammpswrapper.h @@ -34,6 +34,7 @@ public: void *extract_global(const char *keyword); void *extract_pair(const char *keyword); void *extract_atom(const char *keyword); + double extract_variable(const char *keyword); int id_count(const char *idtype); int id_name(const char *idtype, int idx, char *buf, int buflen); diff --git a/tools/lammps-gui/logwindow.cpp b/tools/lammps-gui/logwindow.cpp index 56dce35179..4527bcb0dd 100644 --- a/tools/lammps-gui/logwindow.cpp +++ b/tools/lammps-gui/logwindow.cpp @@ -17,7 +17,6 @@ #include #include -#include #include #include #include diff --git a/tools/lammps-gui/main.cpp b/tools/lammps-gui/main.cpp index d70e9d3e46..4820d48ebb 100644 --- a/tools/lammps-gui/main.cpp +++ b/tools/lammps-gui/main.cpp @@ -14,7 +14,6 @@ #include "lammpsgui.h" #include -#include #include #include diff --git a/tools/lammps-gui/preferences.cpp b/tools/lammps-gui/preferences.cpp index 207b68cb66..b601a9995f 100644 --- a/tools/lammps-gui/preferences.cpp +++ b/tools/lammps-gui/preferences.cpp @@ -23,7 +23,6 @@ #include #include #include -#include #include #include #include @@ -40,7 +39,9 @@ #include #include #include +#if defined(_OPENMP) #include +#endif #include #if defined(_OPENMP) diff --git a/tools/lammps-gui/setvariables.cpp b/tools/lammps-gui/setvariables.cpp index db5eb1cdea..1b2a54df8a 100644 --- a/tools/lammps-gui/setvariables.cpp +++ b/tools/lammps-gui/setvariables.cpp @@ -14,12 +14,12 @@ #include "setvariables.h" #include -#include #include #include #include #include #include +#include SetVariables::SetVariables(QList> &_vars, QWidget *parent) : QDialog(parent), vars(_vars), layout(new QVBoxLayout) diff --git a/tools/lammps-gui/slideshow.cpp b/tools/lammps-gui/slideshow.cpp index dcc85fc34a..08b854becd 100644 --- a/tools/lammps-gui/slideshow.cpp +++ b/tools/lammps-gui/slideshow.cpp @@ -32,7 +32,6 @@ #include #include #include -#include #include #include #include diff --git a/unittest/c-library/CMakeLists.txt b/unittest/c-library/CMakeLists.txt index f5793a804e..62ec750a3f 100644 --- a/unittest/c-library/CMakeLists.txt +++ b/unittest/c-library/CMakeLists.txt @@ -3,6 +3,7 @@ add_executable(test_library_open test_library_open.cpp test_main.cpp) target_link_libraries(test_library_open PRIVATE lammps GTest::GMock) add_test(NAME LibraryOpen COMMAND test_library_open) +set_tests_properties(LibraryOpen PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=4;OMP_PROC_BIND=false") add_executable(test_library_commands test_library_commands.cpp test_main.cpp) target_link_libraries(test_library_commands PRIVATE lammps GTest::GMock) @@ -16,7 +17,7 @@ add_executable(test_library_properties test_library_properties.cpp test_main.cpp target_link_libraries(test_library_properties PRIVATE lammps GTest::GMock) target_compile_definitions(test_library_properties PRIVATE -DTEST_INPUT_FOLDER=${CMAKE_CURRENT_SOURCE_DIR}) add_test(NAME LibraryProperties COMMAND test_library_properties) -set_tests_properties(LibraryProperties PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR}") +set_tests_properties(LibraryProperties PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};OMP_NUM_THREADS=4;OMP_PROC_BIND=false") add_executable(test_library_objects test_library_objects.cpp test_main.cpp) target_link_libraries(test_library_objects PRIVATE lammps GTest::GMock) diff --git a/unittest/c-library/test_library_open.cpp b/unittest/c-library/test_library_open.cpp index 426b2adaa7..6f2ea9ac33 100644 --- a/unittest/c-library/test_library_open.cpp +++ b/unittest/c-library/test_library_open.cpp @@ -2,10 +2,12 @@ #include "lammps.h" #define LAMMPS_LIB_MPI 1 +#include "info.h" #include "library.h" #include // for stdin, stdout #include #include +#include #include "gmock/gmock.h" #include "gtest/gtest.h" @@ -78,9 +80,38 @@ TEST(lammps_open, with_args) TEST(lammps_open, with_kokkos) { if (!LAMMPS_NS::LAMMPS::is_installed_pkg("KOKKOS")) GTEST_SKIP(); - const char *args[] = {"liblammps", "-k", "on", "t", "2", "-sf", "kk", "-log", "none", nullptr}; - char **argv = (char **)args; - int argc = (sizeof(args) / sizeof(char *)) - 1; + std::vector args = {(char *)"lammps", (char *)"-log", (char *)"none", (char *)"-echo", + (char *)"screen", (char *)"-sf", (char *)"kk"}; + + char *one = (char *)"1"; + char *four = (char *)"4"; + char *tee = (char *)"t"; + char *gee = (char *)"g"; + char *kay = (char *)"-k"; + char *yes = (char *)"on"; + + args.push_back(kay); + args.push_back(yes); + + // when GPU support is enabled in KOKKOS, it *must* be used + if (lammps_config_accelerator("KOKKOS", "api", "hip") || + lammps_config_accelerator("KOKKOS", "api", "cuda") || + lammps_config_accelerator("KOKKOS", "api", "sycl")) { + args.push_back(gee); + args.push_back(one); + } + + // use threads or serial + args.push_back(tee); + if (lammps_config_accelerator("KOKKOS", "api", "openmp")) { + args.push_back(four); + } else if (lammps_config_accelerator("KOKKOS", "api", "pthreads")) { + args.push_back(four); + } else { + args.push_back(one); + } + int argc = args.size(); + char **argv = args.data(); ::testing::internal::CaptureStdout(); void *alt_ptr; diff --git a/unittest/c-library/test_library_properties.cpp b/unittest/c-library/test_library_properties.cpp index 3900265f8c..a2df55d120 100644 --- a/unittest/c-library/test_library_properties.cpp +++ b/unittest/c-library/test_library_properties.cpp @@ -7,6 +7,7 @@ #include "lmptype.h" #include "platform.h" #include +#include #include "gmock/gmock.h" #include "gtest/gtest.h" @@ -742,14 +743,45 @@ TEST_F(AtomProperties, position) TEST(SystemSettings, kokkos) { if (!lammps_config_has_package("KOKKOS")) GTEST_SKIP(); - if (!lammps_config_accelerator("KOKKOS", "api", "openmp")) GTEST_SKIP(); + std::vector args = {(char *)"lammps", (char *)"-log", (char *)"none", + (char *)"-echo", (char *)"screen", (char *)"-nocite", + (char *)"-sf", (char *)"kk"}; - // clang-format off - const char *args[] = {"SystemSettings", "-log", "none", "-echo", "screen", "-nocite", - "-k", "on", "t", "4", "-sf", "kk", nullptr}; - // clang-format on - char **argv = (char **)args; - int argc = (sizeof(args) / sizeof(char *)) - 1; + char *one = (char *)"1"; + char *four = (char *)"4"; + char *tee = (char *)"t"; + char *gee = (char *)"g"; + char *kay = (char *)"-k"; + char *yes = (char *)"on"; + + args.push_back(kay); + args.push_back(yes); + + bool has_gpu = false; + bool has_threads = false; + + // when GPU support is enabled in KOKKOS, it *must* be used + if (lammps_config_accelerator("KOKKOS", "api", "hip") || + lammps_config_accelerator("KOKKOS", "api", "cuda") || + lammps_config_accelerator("KOKKOS", "api", "sycl")) { + has_gpu = true; + args.push_back(gee); + args.push_back(one); + } + + // use threads or serial + args.push_back(tee); + if (lammps_config_accelerator("KOKKOS", "api", "openmp")) { + has_threads = true; + args.push_back(four); + } else if (lammps_config_accelerator("KOKKOS", "api", "pthreads")) { + has_threads = true; + args.push_back(four); + } else { + args.push_back(one); + } + int argc = args.size(); + char **argv = args.data(); ::testing::internal::CaptureStdout(); void *lmp = lammps_open_no_mpi(argc, argv, nullptr); @@ -758,7 +790,13 @@ TEST(SystemSettings, kokkos) EXPECT_THAT(output, StartsWith("LAMMPS (")); EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_active"), 1); - EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_nthreads"), 4); - EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_ngpus"), 0); + if (has_threads) + EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_nthreads"), 4); + else + EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_nthreads"), 1); + if (has_gpu) + EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_ngpus"), 1); + else + EXPECT_EQ(lammps_extract_setting(lmp, "kokkos_ngpus"), 0); lammps_close(lmp); } diff --git a/unittest/commands/test_groups.cpp b/unittest/commands/test_groups.cpp index efeb00f685..4035a4a4d4 100644 --- a/unittest/commands/test_groups.cpp +++ b/unittest/commands/test_groups.cpp @@ -314,7 +314,7 @@ TEST_F(GroupTest, Dynamic) command("group ramp variable grow");); } -constexpr double EPSILON = 1.0e-13; +static constexpr double EPSILON = 1.0e-13; TEST_F(GroupTest, VariableFunctions) { diff --git a/unittest/commands/test_variables.cpp b/unittest/commands/test_variables.cpp index 8f2f4ba9d9..2390b1b675 100644 --- a/unittest/commands/test_variables.cpp +++ b/unittest/commands/test_variables.cpp @@ -287,6 +287,7 @@ TEST_F(VariableTest, AtomicSystem) ASSERT_DOUBLE_EQ(variable->compute_equal("f_press[1]"), 0.0); ASSERT_DOUBLE_EQ(variable->compute_equal("c_press"), 0.0); ASSERT_DOUBLE_EQ(variable->compute_equal("c_press[2]"), 0.0); + ASSERT_DOUBLE_EQ(variable->compute_equal("c_press[1+1]"), 0.0); ASSERT_DOUBLE_EQ(variable->compute_equal("1.5+3.25"), 4.75); ASSERT_DOUBLE_EQ(variable->compute_equal("-2.5*1.5"), -3.75); @@ -302,8 +303,18 @@ TEST_F(VariableTest, AtomicSystem) variable->compute_equal("v_self");); TEST_FAILURE(".*ERROR: Variable sum2: Inconsistent lengths in vector-style variable.*", variable->compute_equal("max(v_sum2)");); - TEST_FAILURE("ERROR: Mismatched fix in variable formula.*", + TEST_FAILURE(".*ERROR: Mismatched fix in variable formula.*", variable->compute_equal("f_press");); + TEST_FAILURE(".*ERROR .*Variable formula compute vector is accessed out-of-range.*", + variable->compute_equal("c_press[10]");); + TEST_FAILURE(".*ERROR: Non digit character between brackets in variable.*", + variable->compute_equal("c_press[axy]");); + TEST_FAILURE(".*ERROR: Illegal value in brackets: stoll.*", + variable->compute_equal("c_press[73786976294838206464]");); + TEST_FAILURE(".*ERROR: Index between variable brackets must be positive.*", + variable->compute_equal("c_press[-2]");); + TEST_FAILURE(".*ERROR: Index between variable brackets must be positive.*", + variable->compute_equal("c_press[0]");); } TEST_F(VariableTest, Expressions) diff --git a/unittest/cplusplus/CMakeLists.txt b/unittest/cplusplus/CMakeLists.txt index 445e0fffeb..dac4cf3e10 100644 --- a/unittest/cplusplus/CMakeLists.txt +++ b/unittest/cplusplus/CMakeLists.txt @@ -3,7 +3,7 @@ add_executable(test_lammps_class test_lammps_class.cpp) target_link_libraries(test_lammps_class PRIVATE lammps GTest::GMockMain) add_test(NAME LammpsClass COMMAND test_lammps_class) -set_tests_properties(LammpsClass PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=1") +set_tests_properties(LammpsClass PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=2;OMP_PROC_BIND=false") add_executable(test_input_class test_input_class.cpp) target_link_libraries(test_input_class PRIVATE lammps GTest::GTestMain) diff --git a/unittest/cplusplus/test_lammps_class.cpp b/unittest/cplusplus/test_lammps_class.cpp index 23d83c4ecb..464ecca925 100644 --- a/unittest/cplusplus/test_lammps_class.cpp +++ b/unittest/cplusplus/test_lammps_class.cpp @@ -253,6 +253,15 @@ protected: { LAMMPS::argv args = {"LAMMPS_test", "-log", "none", "-echo", "none", "-screen", "none", "-k", "on", "t", "1", "-sf", "kk"}; + + // when GPU support is enabled in KOKKOS, it *must* be used + if (Info::has_accelerator_feature("KOKKOS", "api", "hip") || + Info::has_accelerator_feature("KOKKOS", "api", "cuda") || + Info::has_accelerator_feature("KOKKOS", "api", "sycl")) { + args = {"LAMMPS_test", "-log", "none", "-echo", "none", "-screen", "none", "-k", + "on", "t", "1", "g", "1", "-sf", "kk"}; + } + if (Info::has_accelerator_feature("KOKKOS", "api", "openmp")) args[10] = "2"; if (LAMMPS::is_installed_pkg("KOKKOS")) { diff --git a/unittest/force-styles/CMakeLists.txt b/unittest/force-styles/CMakeLists.txt index 7d29395e39..f47bb75305 100644 --- a/unittest/force-styles/CMakeLists.txt +++ b/unittest/force-styles/CMakeLists.txt @@ -204,11 +204,7 @@ foreach(TEST ${FIX_TIMESTEP_TESTS}) continue() endif() add_test(NAME ${TNAME} COMMAND test_fix_timestep ${TEST}) -if(WIN32) - set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}\\\;${LAMMPS_PYTHON_DIR}") -else() - set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:${LAMMPS_PYTHON_DIR}:$ENV{PYTHONPATH};PYTHONDONTWRITEBYTECODE=1") -endif() + set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "${FORCE_TEST_ENVIRONMENT}") set_tests_properties(${TNAME} PROPERTIES LABELS "${TEST_TAGS}") endforeach() @@ -225,7 +221,7 @@ foreach(TEST ${DIHEDRAL_TESTS}) continue() endif() add_test(NAME ${TNAME} COMMAND test_dihedral_style ${TEST}) - set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH};PYTHONDONTWRITEBYTECODE=1") + set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "${FORCE_TEST_ENVIRONMENT}") set_tests_properties(${TNAME} PROPERTIES LABELS "${TEST_TAGS}") endforeach() @@ -242,7 +238,7 @@ foreach(TEST ${IMPROPER_TESTS}) continue() endif() add_test(NAME ${TNAME} COMMAND test_improper_style ${TEST}) - set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH};PYTHONDONTWRITEBYTECODE=1") + set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "${FORCE_TEST_ENVIRONMENT}") set_tests_properties(${TNAME} PROPERTIES LABELS "${TEST_TAGS}") endforeach() @@ -250,7 +246,7 @@ if(MLIAP_ENABLE_PYTHON AND (NOT WIN32)) add_executable(test_mliappy_unified test_mliappy_unified.cpp) target_link_libraries(test_mliappy_unified PRIVATE lammps GTest::GMockMain) add_test(NAME TestMliapPyUnified COMMAND test_mliappy_unified) - set_tests_properties(TestMliapPyUnified PROPERTIES ENVIRONMENT "PYTHONPATH=${LAMMPS_PYTHON_DIR};PYTHONDONTWRITEBYTECODE=1") + set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "${FORCE_TEST_ENVIRONMENT}") endif() add_executable(test_pair_list test_pair_list.cpp) diff --git a/unittest/force-styles/test_angle_style.cpp b/unittest/force-styles/test_angle_style.cpp index e2b76cc7b6..e706cc11ac 100644 --- a/unittest/force-styles/test_angle_style.cpp +++ b/unittest/force-styles/test_angle_style.cpp @@ -530,6 +530,126 @@ TEST(AngleStyle, omp) if (!verbose) ::testing::internal::GetCapturedStdout(); }; +TEST(AngleStyle, kokkos_omp) +{ + if (!LAMMPS::is_installed_pkg("KOKKOS")) GTEST_SKIP(); + if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP(); + if (!Info::has_accelerator_feature("KOKKOS", "api", "openmp")) GTEST_SKIP(); + + LAMMPS::argv args = {"AngleStyle", "-log", "none", "-echo", "screen", "-nocite", + "-k", "on", "t", "4", "-sf", "kk"}; + + ::testing::internal::CaptureStdout(); + LAMMPS *lmp = init_lammps(args, test_config, true); + + std::string output = ::testing::internal::GetCapturedStdout(); + if (verbose) std::cout << output; + + if (!lmp) { + std::cerr << "One or more prerequisite styles with /kk suffix\n" + "are not available in this LAMMPS configuration:\n"; + for (auto &prerequisite : test_config.prerequisites) { + std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n"; + } + GTEST_SKIP(); + } + + EXPECT_THAT(output, StartsWith("LAMMPS (")); + EXPECT_THAT(output, HasSubstr("Loop time")); + + // abort if running in parallel and not all atoms are local + const int nlocal = lmp->atom->nlocal; + ASSERT_EQ(lmp->atom->natoms, nlocal); + + // relax error a bit for KOKKOS package + double epsilon = 5.0 * test_config.epsilon; + + ErrorStats stats; + auto angle = lmp->force->angle; + + EXPECT_FORCES("init_forces (newton on)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton on)", angle->virial, test_config.init_stress, epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(angle->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton on)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton on)", angle->virial, test_config.run_stress, epsilon); + + stats.reset(); + int id = lmp->modify->find_compute("sum"); + double energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(angle->energy, test_config.run_energy, epsilon); + EXPECT_FP_LE_WITH_EPS(angle->energy, energy, epsilon); + if (print_stats) std::cerr << "run_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + lmp = init_lammps(args, test_config, false); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + // skip over these tests if newton bond is forced to be on + if (lmp->force->newton_bond == 0) { + angle = lmp->force->angle; + + EXPECT_FORCES("init_forces (newton off)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton off)", angle->virial, test_config.init_stress, + 2 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(angle->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton off)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton off)", angle->virial, test_config.run_stress, epsilon); + + stats.reset(); + id = lmp->modify->find_compute("sum"); + energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(angle->energy, test_config.run_energy, epsilon); + EXPECT_FP_LE_WITH_EPS(angle->energy, energy, epsilon); + if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl; + } + + if (!verbose) ::testing::internal::CaptureStdout(); + restart_lammps(lmp, test_config); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + angle = lmp->force->angle; + EXPECT_FORCES("restart_forces", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("restart_stress", angle->virial, test_config.init_stress, epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(angle->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "restart_energy stats:" << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + data_lammps(lmp, test_config); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + angle = lmp->force->angle; + EXPECT_FORCES("data_forces", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("data_stress", angle->virial, test_config.init_stress, epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(angle->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "data_energy stats:" << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + if (!verbose) ::testing::internal::GetCapturedStdout(); +}; + + TEST(AngleStyle, numdiff) { if (!LAMMPS::is_installed_pkg("EXTRA-FIX")) GTEST_SKIP(); diff --git a/unittest/force-styles/test_bond_style.cpp b/unittest/force-styles/test_bond_style.cpp index d2523ff51d..81b105259d 100644 --- a/unittest/force-styles/test_bond_style.cpp +++ b/unittest/force-styles/test_bond_style.cpp @@ -532,6 +532,112 @@ TEST(BondStyle, omp) if (!verbose) ::testing::internal::GetCapturedStdout(); }; +TEST(BondStyle, kokkos_omp) +{ + if (!LAMMPS::is_installed_pkg("KOKKOS")) GTEST_SKIP(); + if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP(); + if (!Info::has_accelerator_feature("KOKKOS", "api", "openmp")) GTEST_SKIP(); + + LAMMPS::argv args = {"BondStyle", "-log", "none", "-echo", "screen", "-nocite", + "-k", "on", "t", "4", "-sf", "kk"}; + + ::testing::internal::CaptureStdout(); + LAMMPS *lmp = init_lammps(args, test_config, true); + + std::string output = ::testing::internal::GetCapturedStdout(); + if (verbose) std::cout << output; + + if (!lmp) { + std::cerr << "One or more prerequisite styles with /kk suffix\n" + "are not available in this LAMMPS configuration:\n"; + for (auto &prerequisite : test_config.prerequisites) { + std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n"; + } + GTEST_SKIP(); + } + + EXPECT_THAT(output, StartsWith("LAMMPS (")); + EXPECT_THAT(output, HasSubstr("Loop time")); + + // abort if running in parallel and not all atoms are local + const int nlocal = lmp->atom->nlocal; + ASSERT_EQ(lmp->atom->natoms, nlocal); + + // relax error a bit for KOKKOS package + double epsilon = 5.0 * test_config.epsilon; + + ErrorStats stats; + auto bond = lmp->force->bond; + + EXPECT_FORCES("init_forces (newton on)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton on)", bond->virial, test_config.init_stress, 10 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(bond->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton on)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton on)", bond->virial, test_config.run_stress, 10 * epsilon); + + stats.reset(); + int id = lmp->modify->find_compute("sum"); + double energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(bond->energy, test_config.run_energy, epsilon); + + // FIXME: this is currently broken ??? for KOKKOS with bond style hybrid + // needs to be fixed in the main code somewhere. Not sure where, though. + //if (test_config.bond_style.substr(0, 6) != "hybrid") + // EXPECT_FP_LE_WITH_EPS(bond->energy, energy, epsilon); + + if (print_stats) std::cerr << "run_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + lmp = init_lammps(args, test_config, false); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + // skip over these tests if newton bond is forced to be on + if (lmp->force->newton_bond == 0) { + bond = lmp->force->bond; + + EXPECT_FORCES("init_forces (newton off)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton off)", bond->virial, test_config.init_stress, + 10 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(bond->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton off)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton off)", bond->virial, test_config.run_stress, + 10 * epsilon); + + stats.reset(); + id = lmp->modify->find_compute("sum"); + energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(bond->energy, test_config.run_energy, epsilon); + + // FIXME: this is currently broken ??? for KOKKOS with bond style hybrid + // needs to be fixed in the main code somewhere. Not sure where, though. + //if (test_config.bond_style.substr(0, 6) != "hybrid") + // EXPECT_FP_LE_WITH_EPS(bond->energy, energy, epsilon); + + if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl; + } + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + if (!verbose) ::testing::internal::GetCapturedStdout(); +}; + TEST(BondStyle, numdiff) { if (!LAMMPS::is_installed_pkg("EXTRA-FIX")) GTEST_SKIP(); diff --git a/unittest/force-styles/test_dihedral_style.cpp b/unittest/force-styles/test_dihedral_style.cpp index e6c34843c3..324745cc41 100644 --- a/unittest/force-styles/test_dihedral_style.cpp +++ b/unittest/force-styles/test_dihedral_style.cpp @@ -534,6 +534,112 @@ TEST(DihedralStyle, omp) if (!verbose) ::testing::internal::GetCapturedStdout(); }; +TEST(DihedralStyle, kokkos_omp) +{ + if (!LAMMPS::is_installed_pkg("KOKKOS")) GTEST_SKIP(); + if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP(); + if (!Info::has_accelerator_feature("KOKKOS", "api", "openmp")) GTEST_SKIP(); + + LAMMPS::argv args = {"DihedralStyle", "-log", "none", "-echo", "screen", "-nocite", + "-k", "on", "t", "4", "-sf", "kk"}; + + ::testing::internal::CaptureStdout(); + LAMMPS *lmp = init_lammps(args, test_config, true); + + std::string output = ::testing::internal::GetCapturedStdout(); + if (verbose) std::cout << output; + + if (!lmp) { + std::cerr << "One or more prerequisite styles with /kk suffix\n" + "are not available in this LAMMPS configuration:\n"; + for (auto &prerequisite : test_config.prerequisites) { + std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n"; + } + GTEST_SKIP(); + } + + EXPECT_THAT(output, StartsWith("LAMMPS (")); + EXPECT_THAT(output, HasSubstr("Loop time")); + + // abort if running in parallel and not all atoms are local + const int nlocal = lmp->atom->nlocal; + ASSERT_EQ(lmp->atom->natoms, nlocal); + + // relax error a bit for KOKKOS package + double epsilon = 5.0 * test_config.epsilon; + + ErrorStats stats; + auto dihedral = lmp->force->dihedral; + + EXPECT_FORCES("init_forces (newton on)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton on)", dihedral->virial, test_config.init_stress, + 10 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton on)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton on)", dihedral->virial, test_config.run_stress, 10 * epsilon); + + stats.reset(); + int id = lmp->modify->find_compute("sum"); + double energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.run_energy, epsilon); + + // FIXME: this is currently broken ??? for KOKKOS with dihedral style hybrid + // needs to be fixed in the main code somewhere. Not sure where, though. + //if (test_config.dihedral_style.substr(0, 6) != "hybrid") + // EXPECT_FP_LE_WITH_EPS(dihedral->energy, energy, epsilon); + //if (print_stats) std::cerr << "run_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + lmp = init_lammps(args, test_config, false); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + // skip over these tests if newton bond is forced to be on + if (lmp->force->newton_bond == 0) { + dihedral = lmp->force->dihedral; + + EXPECT_FORCES("init_forces (newton off)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton off)", dihedral->virial, test_config.init_stress, + 10 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton off)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton off)", dihedral->virial, test_config.run_stress, + 10 * epsilon); + + stats.reset(); + id = lmp->modify->find_compute("sum"); + energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.run_energy, epsilon); + + // FIXME: this is currently broken ??? for KOKKOS with dihedral style hybrid + // needs to be fixed in the main code somewhere. Not sure where, though. + //if (test_config.dihedral_style.substr(0, 6) != "hybrid") + // EXPECT_FP_LE_WITH_EPS(dihedral->energy, energy, epsilon); + + if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl; + } + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + if (!verbose) ::testing::internal::GetCapturedStdout(); +}; + TEST(DihedralStyle, numdiff) { if (!LAMMPS::is_installed_pkg("EXTRA-FIX")) GTEST_SKIP(); diff --git a/unittest/force-styles/test_improper_style.cpp b/unittest/force-styles/test_improper_style.cpp index 59028fb1ff..a65bf438df 100644 --- a/unittest/force-styles/test_improper_style.cpp +++ b/unittest/force-styles/test_improper_style.cpp @@ -527,6 +527,108 @@ TEST(ImproperStyle, omp) if (!verbose) ::testing::internal::GetCapturedStdout(); }; +TEST(ImproperStyle, kokkos_omp) +{ + if (!LAMMPS::is_installed_pkg("KOKKOS")) GTEST_SKIP(); + if (test_config.skip_tests.count(test_info_->name())) GTEST_SKIP(); + + LAMMPS::argv args = {"ImproperStyle", "-log", "none", "-echo", "screen", "-nocite", + "-k", "on", "t", "4", "-sf", "kk"}; + + ::testing::internal::CaptureStdout(); + LAMMPS *lmp = init_lammps(args, test_config, true); + + std::string output = ::testing::internal::GetCapturedStdout(); + if (verbose) std::cout << output; + + if (!lmp) { + std::cerr << "One or more prerequisite styles with /kk suffix\n" + "are not available in this LAMMPS configuration:\n"; + for (auto &prerequisite : test_config.prerequisites) { + std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n"; + } + GTEST_SKIP(); + } + + EXPECT_THAT(output, StartsWith("LAMMPS (")); + EXPECT_THAT(output, HasSubstr("Loop time")); + + // abort if running in parallel and not all atoms are local + const int nlocal = lmp->atom->nlocal; + ASSERT_EQ(lmp->atom->natoms, nlocal); + + // relax error a bit for KOKKOS package + double epsilon = 5.0 * test_config.epsilon; + + ErrorStats stats; + auto improper = lmp->force->improper; + + EXPECT_FORCES("init_forces (newton on)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton on)", improper->virial, test_config.init_stress, + 10 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(improper->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton on)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton on)", improper->virial, test_config.run_stress, 10 * epsilon); + + stats.reset(); + int id = lmp->modify->find_compute("sum"); + double energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(improper->energy, test_config.run_energy, epsilon); + // FIXME: this is currently broken ??? for KOKKOS with improper style hybrid + // needs to be fixed in the main code somewhere. Not sure where, though. + //if (test_config.improper_style.substr(0, 6) != "hybrid") + // EXPECT_FP_LE_WITH_EPS(improper->energy, energy, epsilon); + if (print_stats) std::cerr << "run_energy stats, newton on: " << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + lmp = init_lammps(args, test_config, false); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + // skip over these tests if newton bond is forced to be on + if (lmp->force->newton_bond == 0) { + improper = lmp->force->improper; + + EXPECT_FORCES("init_forces (newton off)", lmp->atom, test_config.init_forces, epsilon); + EXPECT_STRESS("init_stress (newton off)", improper->virial, test_config.init_stress, + 10 * epsilon); + + stats.reset(); + EXPECT_FP_LE_WITH_EPS(improper->energy, test_config.init_energy, epsilon); + if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl; + + if (!verbose) ::testing::internal::CaptureStdout(); + run_lammps(lmp); + if (!verbose) ::testing::internal::GetCapturedStdout(); + + EXPECT_FORCES("run_forces (newton off)", lmp->atom, test_config.run_forces, 10 * epsilon); + EXPECT_STRESS("run_stress (newton off)", improper->virial, test_config.run_stress, + 10 * epsilon); + + stats.reset(); + id = lmp->modify->find_compute("sum"); + energy = lmp->modify->compute[id]->compute_scalar(); + EXPECT_FP_LE_WITH_EPS(improper->energy, test_config.run_energy, epsilon); + // FIXME: this is currently broken ??? for KOKKOS with improper style hybrid + // needs to be fixed in the main code somewhere. Not sure where, though. + //if (test_config.improper_style.substr(0, 6) != "hybrid") + // EXPECT_FP_LE_WITH_EPS(improper->energy, energy, epsilon); + if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl; + } + + if (!verbose) ::testing::internal::CaptureStdout(); + cleanup_lammps(lmp, test_config); + if (!verbose) ::testing::internal::GetCapturedStdout(); +}; + TEST(ImproperStyle, numdiff) { if (!LAMMPS::is_installed_pkg("EXTRA-FIX")) GTEST_SKIP(); diff --git a/unittest/force-styles/tests/improper-class2.yaml b/unittest/force-styles/tests/improper-class2.yaml index 38e8b9120b..4350441eeb 100644 --- a/unittest/force-styles/tests/improper-class2.yaml +++ b/unittest/force-styles/tests/improper-class2.yaml @@ -2,7 +2,7 @@ lammps_version: 17 Feb 2022 date_generated: Fri Mar 18 22:18:02 2022 epsilon: 2.5e-13 -skip_tests: +skip_tests: kokkos_omp prerequisites: ! | atom full improper class2 diff --git a/unittest/formats/test_atom_styles.cpp b/unittest/formats/test_atom_styles.cpp index 68bc0a4437..921d469e31 100644 --- a/unittest/formats/test_atom_styles.cpp +++ b/unittest/formats/test_atom_styles.cpp @@ -88,7 +88,7 @@ static void create_molecule_files(const std::string &h2o_filename, const std::st // whether to print verbose output (i.e. not capturing LAMMPS screen output). bool verbose = false; -const double EPSILON = 5.0e-14; +static const double EPSILON = 5.0e-14; namespace LAMMPS_NS { using ::testing::Eq; diff --git a/unittest/formats/test_molecule_file.cpp b/unittest/formats/test_molecule_file.cpp index 743a8fbbfa..e11a8ea4f1 100644 --- a/unittest/formats/test_molecule_file.cpp +++ b/unittest/formats/test_molecule_file.cpp @@ -32,7 +32,7 @@ using testing::StrEq; using utils::split_words; -const double EPSILON = 5.0e-14; +static constexpr double EPSILON = 5.0e-14; #define test_name test_info_->name() diff --git a/unittest/utils/test_tokenizer.cpp b/unittest/utils/test_tokenizer.cpp index 6db4a5fe5e..d5d6955b32 100644 --- a/unittest/utils/test_tokenizer.cpp +++ b/unittest/utils/test_tokenizer.cpp @@ -11,11 +11,14 @@ See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ +#include "fmt/format.h" #include "lmptype.h" #include "tokenizer.h" #include "gmock/gmock.h" #include "gtest/gtest.h" +#include + using namespace LAMMPS_NS; using ::testing::Eq; @@ -334,46 +337,73 @@ TEST(ValueTokenizer, move_assignment) TEST(ValueTokenizer, bad_integer) { - ValueTokenizer values("f10 f11 f12"); + ValueTokenizer values("f10 f11 f12 0xff 109951162777 " + "36893488147419103232 36893488147419103232"); ASSERT_THROW(values.next_int(), InvalidIntegerException); ASSERT_THROW(values.next_bigint(), InvalidIntegerException); ASSERT_THROW(values.next_tagint(), InvalidIntegerException); + ASSERT_THROW(values.next_int(), InvalidIntegerException); + ASSERT_THROW(values.next_int(), InvalidIntegerException); + ASSERT_THROW(values.next_tagint(), InvalidIntegerException); + ASSERT_THROW(values.next_bigint(), InvalidIntegerException); } TEST(ValueTokenizer, bad_double) { - ValueTokenizer values("1a.0"); + ValueTokenizer values("1a.0 --2.0 2.4d3 -1e20000 1.0e-1.0"); + ASSERT_THROW(values.next_double(), InvalidFloatException); + ASSERT_THROW(values.next_double(), InvalidFloatException); + ASSERT_THROW(values.next_double(), InvalidFloatException); + ASSERT_THROW(values.next_double(), InvalidFloatException); ASSERT_THROW(values.next_double(), InvalidFloatException); } TEST(ValueTokenizer, valid_int) { - ValueTokenizer values("10"); + ValueTokenizer values(fmt::format("10 {} {}", -MAXSMALLINT - 1, MAXSMALLINT)); ASSERT_EQ(values.next_int(), 10); + ASSERT_EQ(values.next_int(), -MAXSMALLINT - 1); + ASSERT_EQ(values.next_int(), MAXSMALLINT); } TEST(ValueTokenizer, valid_tagint) { - ValueTokenizer values("42"); + ValueTokenizer values( + fmt::format("42 {} {} {} {}", -MAXSMALLINT - 1, MAXSMALLINT, -MAXTAGINT - 1, MAXTAGINT)); ASSERT_EQ(values.next_tagint(), 42); + ASSERT_EQ(values.next_tagint(), -MAXSMALLINT - 1); + ASSERT_EQ(values.next_tagint(), MAXSMALLINT); + ASSERT_EQ(values.next_tagint(), -MAXTAGINT - 1); + ASSERT_EQ(values.next_tagint(), MAXTAGINT); } TEST(ValueTokenizer, valid_bigint) { - ValueTokenizer values("42"); + ValueTokenizer values( + fmt::format("42 {} {} {} {}", -MAXSMALLINT - 1, MAXSMALLINT, -MAXBIGINT - 1, MAXBIGINT)); ASSERT_EQ(values.next_bigint(), 42); + ASSERT_EQ(values.next_bigint(), -MAXSMALLINT - 1); + ASSERT_EQ(values.next_bigint(), MAXSMALLINT); + ASSERT_EQ(values.next_bigint(), -MAXBIGINT - 1); + ASSERT_EQ(values.next_bigint(), MAXBIGINT); } TEST(ValueTokenizer, valid_double) { - ValueTokenizer values("3.14"); + ValueTokenizer values("3.14 -0.00002 .1 0xff " + std::to_string(MAXBIGINT)); ASSERT_DOUBLE_EQ(values.next_double(), 3.14); + ASSERT_DOUBLE_EQ(values.next_double(), -0.00002); + ASSERT_DOUBLE_EQ(values.next_double(), 0.1); + ASSERT_DOUBLE_EQ(values.next_double(), 255); + ASSERT_DOUBLE_EQ(values.next_double(), MAXBIGINT); } TEST(ValueTokenizer, valid_double_with_exponential) { - ValueTokenizer values("3.14e22"); + ValueTokenizer values(fmt::format("3.14e22 {} {}", DBL_MAX, DBL_MIN)); ASSERT_DOUBLE_EQ(values.next_double(), 3.14e22); + ASSERT_DOUBLE_EQ(values.next_double(), DBL_MAX); + ASSERT_DOUBLE_EQ(values.next_double(), DBL_MIN); } TEST(ValueTokenizer, contains)