Merge release-2019 into master
authorPaul Bauer <paul.bauer.q@gmail.com>
Wed, 31 Jul 2019 13:28:34 +0000 (15:28 +0200)
committerBerk Hess <hess@kth.se>
Thu, 1 Aug 2019 07:53:51 +0000 (09:53 +0200)
Resolved Conflicts:
cmake/gmxVersionInfo.cmake
docs/CMakeLists.txt
docs/release-notes/2019/2019.2.rst
src/external/build-fftw/CMakeLists.txt
src/gromacs/awh/read-params.cpp
src/gromacs/awh/read-params.h
src/gromacs/domdec/domdec_internal.h
src/gromacs/ewald/pme-gather.clh
src/gromacs/fileio/gmx_internal_xdr.cpp
src/gromacs/gmxana/gmx_editconf.cpp
src/gromacs/gmxpreprocess/readir.cpp
src/gromacs/listed-forces/manage-threading.cpp
src/gromacs/mdlib/broadcaststructs.cpp
src/gromacs/mdlib/membed.cpp
src/gromacs/mdlib/nbnxn_atomdata.cpp
src/gromacs/mdlib/nbnxn_atomdata.h
src/gromacs/mdlib/nbnxn_search.cpp
src/gromacs/mdlib/ns.cpp
src/gromacs/mdrun/runner.cpp
src/programs/mdrun/tests/spc-and-methanol.top

Change-Id: I2c0853732620246cade31f3a66ae23b7ba5b902b

61 files changed:
cmake/gmxDetectAvx512FmaUnits.cmake
cmake/gmxManageNvccConfig.cmake
docs/CMakeLists.txt
docs/install-guide/index.rst
docs/reference-manual/algorithms/molecular-dynamics.rst
docs/reference-manual/functions/restraints.rst
docs/reference-manual/references.rst
docs/release-notes/2018/2018.6.rst
docs/release-notes/2018/2018.7.rst [new file with mode: 0644]
docs/release-notes/2019/2019.1.rst
docs/release-notes/2019/2019.2.rst
docs/release-notes/2019/2019.3.rst [new file with mode: 0644]
docs/release-notes/2019/2019.4.rst [new file with mode: 0644]
docs/release-notes/index.rst
docs/user-guide/mdp-options.rst
docs/user-guide/mdrun-performance.rst
docs/user-guide/run-time-errors.rst
docs/user-guide/system-preparation.rst
src/api/cpp/tests/testingconfiguration.h
src/external/build-fftw/CMakeLists.txt
src/external/tinyxml2/tinyxml2.cpp
src/external/tng_io/external/zlib/inflate.c
src/gromacs/awh/awh.cpp
src/gromacs/awh/read_params.cpp
src/gromacs/awh/read_params.h
src/gromacs/domdec/domdec.cpp
src/gromacs/domdec/domdec_internal.h
src/gromacs/domdec/domdec_struct.h
src/gromacs/ewald/pme.cpp
src/gromacs/ewald/pme_gather.clh
src/gromacs/fileio/gmx_internal_xdr.cpp
src/gromacs/fileio/pdbio.cpp
src/gromacs/gmxana/gmx_anaeig.cpp
src/gromacs/gmxana/gmx_cluster.cpp
src/gromacs/gmxana/gmx_hbond.cpp
src/gromacs/gmxana/gmx_wham.cpp
src/gromacs/gmxana/gmx_xpm2ps.cpp
src/gromacs/gmxpreprocess/editconf.cpp
src/gromacs/gmxpreprocess/gmxcpp.cpp
src/gromacs/gmxpreprocess/gmxcpp.h
src/gromacs/gmxpreprocess/grompp.cpp
src/gromacs/gmxpreprocess/toppush.cpp
src/gromacs/hardware/tests/hardwaretopology.cpp
src/gromacs/listed_forces/manage_threading.cpp
src/gromacs/mdlib/broadcaststructs.cpp
src/gromacs/mdlib/membed.cpp
src/gromacs/mdlib/shake.cpp
src/gromacs/mdrun/legacymdrunoptions.h
src/gromacs/mdrun/minimize.cpp
src/gromacs/mdrun/runner.cpp
src/gromacs/nbnxm/pairlist.cpp
src/gromacs/pulling/pull.cpp
src/gromacs/pulling/pullutil.cpp
src/gromacs/selection/indexutil.cpp
src/gromacs/selection/sm_simple.cpp
src/gromacs/simd/impl_reference/impl_reference_definitions.h
src/programs/mdrun/tests/refdata/MdrunTest_WritesHelp.xml
src/programs/mdrun/tests/refdata/MinimizersWork_EnergyMinimizationTest_WithinTolerances_5.xml
src/testutils/simulationdatabase/spc_and_methane.gro [new file with mode: 0644]
src/testutils/simulationdatabase/spc_and_methane.ndx [new file with mode: 0644]
src/testutils/simulationdatabase/spc_and_methane.top [new file with mode: 0644]

index 5f278b3bff4461e09aa218dadec7d7ef58908f08..7836c4b2a99192442adb121e097692f3bc352675 100644 (file)
@@ -82,7 +82,7 @@ function(gmx_detect_avx_512_fma_units RESULT)
                         ERROR_QUIET)
                     if (RESULT_VAR EQUAL 0)
                         string(STRIP "${OUTPUT_VAR_TEMP}" OUTPUT_VAR)
-                        set(${RESULT} ${OUTPUT_VAR_TEMP} CACHE INTERNAL "Result of test for number of AVX-512 FMA units")
+                        set(${RESULT} ${OUTPUT_VAR} CACHE INTERNAL "Result of test for number of AVX-512 FMA units")
                     else()
                         message(STATUS "Could not identify number of AVX-512 units - detection program did run successfully")
                         set(${RESULT} -1 CACHE INTERNAL "Result of test for number of AVX-512 FMA units")
index f88297e960f2c1a52cecaf9f11a96cefe63bf250..63fe528599d4b3a43f55177cf9be19bd82c8c575 100644 (file)
@@ -172,7 +172,7 @@ gmx_check_if_changed(_cuda_nvcc_executable_or_flags_changed CUDA_NVCC_EXECUTABLE
 # code. Set the CMake variable GMX_NVCC_WORKS on if you want to
 # bypass this check.
 if((_cuda_nvcc_executable_or_flags_changed OR CUDA_HOST_COMPILER_CHANGED OR NOT GMX_NVCC_WORKS) AND NOT WIN32)
-    message(STATUS "Check for working NVCC/C compiler combination")
+    message(STATUS "Check for working NVCC/C++ compiler combination with nvcc '${CUDA_NVCC_EXECUTABLE}'")
     execute_process(COMMAND ${CUDA_NVCC_EXECUTABLE} -ccbin ${CUDA_HOST_COMPILER} -c ${CUDA_NVCC_FLAGS} ${CUDA_NVCC_FLAGS_${_build_type}} ${CMAKE_SOURCE_DIR}/cmake/TestCUDA.cu
         RESULT_VARIABLE _cuda_test_res
         OUTPUT_VARIABLE _cuda_test_out
@@ -184,12 +184,12 @@ if((_cuda_nvcc_executable_or_flags_changed OR CUDA_HOST_COMPILER_CHANGED OR NOT
         message(STATUS "${CUDA_NVCC_EXECUTABLE} standard output: '${_cuda_test_out}'")
         message(STATUS "${CUDA_NVCC_EXECUTABLE} standard error:  '${_cuda_test_err}'")
         if(${_cuda_test_err} MATCHES "nsupported")
-            message(FATAL_ERROR "NVCC/C compiler combination does not seem to be supported. CUDA frequently does not support the latest versions of the host compiler, so you might want to try an earlier C/C++ compiler version and make sure your CUDA compiler and driver are as recent as possible.")
+            message(FATAL_ERROR "NVCC/C++ compiler combination does not seem to be supported. CUDA frequently does not support the latest versions of the host compiler, so you might want to try an earlier C++ compiler version and make sure your CUDA compiler and driver are as recent as possible.")
         else()
             message(FATAL_ERROR "CUDA compiler does not seem to be functional.")
         endif()
     elseif(NOT GMX_CUDA_TEST_COMPILER_QUIETLY)
-        message(STATUS "Check for working NVCC/C compiler combination - works")
+        message(STATUS "Check for working NVCC/C++ compiler combination - works")
         set(GMX_NVCC_WORKS TRUE CACHE INTERNAL "Nvcc can compile a trivial test program")
     endif()
 endif() # GMX_CHECK_NVCC
index 9580f88d4a68449f0242e19e4b831ab2f41b3e02..07d96196939a410427989511b14a35caba20c152 100644 (file)
@@ -369,6 +369,8 @@ if (SPHINX_FOUND)
         release-notes/2020/major/deprecated-functionality.rst
         release-notes/2020/major/portability.rst
         release-notes/2020/major/miscellaneous.rst
+        release-notes/2019/2019.4.rst
+        release-notes/2019/2019.3.rst
         release-notes/2019/2019.2.rst
         release-notes/2019/2019.1.rst
         release-notes/2019/major/highlights.rst
@@ -380,6 +382,7 @@ if (SPHINX_FOUND)
         release-notes/2019/major/deprecated-functionality.rst
         release-notes/2019/major/portability.rst
         release-notes/2019/major/miscellaneous.rst
+        release-notes/2018/2018.7.rst
         release-notes/2018/2018.6.rst
         release-notes/2018/2018.5.rst
         release-notes/2018/2018.4.rst
@@ -599,24 +602,27 @@ if (SPHINX_FOUND)
 
 else()
     set(MANUAL_BUILD_IS_POSSIBLE OFF)
-    set(MANUAL_BUILD_NOT_POSSIBLE_REASON "Sphinx version ${EXPECTED_SPHINX_VERSION} is not available")
+    set(MANUAL_BUILD_NOT_POSSIBLE_REASON "Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is not available")
 
     add_custom_target(webpage-sphinx
         COMMAND ${CMAKE_COMMAND} -E echo
-            "HTML pages cannot be built because Sphinx version ${EXPECTED_SPHINX_VERSION} is not available"
+            "HTML pages cannot be built because Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is not available"
         VERBATIM)
     add_custom_target(install-guide
         COMMAND ${CMAKE_COMMAND} -E echo
-            "INSTALL cannot be built because Sphinx version ${EXPECTED_SPHINX_VERSION} is not available"
+            "INSTALL cannot be built because Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is not available"
         VERBATIM)
     add_custom_target(man
         COMMAND ${CMAKE_COMMAND} -E echo
-            "man pages cannot be built because Sphinx version ${EXPECTED_SPHINX_VERSION} is not available"
+            "man pages cannot be built because Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is not available"
         VERBATIM)
     add_custom_target(sphinx-create-texman
         COMMAND ${CMAKE_COMMAND} -E echo
-            "Cannot prepare LaTeX input files because Sphinx version ${EXPECTED_SPHINX_VERSION} is not available"
+            "Cannot prepare LaTeX input files because Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is not available"
         VERBATIM)
+    add_custom_target(manual
+        COMMAND ${CMAKE_COMMAND} -E echo
+            "manual cannot be built because Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is not available")
 endif()
 
 if (MAN_PAGE_DIR)
@@ -645,7 +651,7 @@ if (NOT PythonInterp_FOUND)
 elseif (NOT SPHINX_FOUND)
     # Hardly anything gets built if Sphinx is not available, so don't bother.
     set(HTML_BUILD_IS_POSSIBLE OFF)
-    set(HTML_BUILD_NOT_POSSIBLE_REASON "Sphinx version ${EXPECTED_SPHINX_VERSION} is required")
+    set(HTML_BUILD_NOT_POSSIBLE_REASON "Sphinx expected minimum version ${EXPECTED_SPHINX_VERSION} is required")
 endif()
 if (NOT MANUAL_BUILD_IS_POSSIBLE)
     list(APPEND HTML_BUILD_WARNINGS
index 2b2c78eaf341e5497cc6bbba02f9f0a1f523a193..48da20d21efdf4b9ab6668490bf0f7ec5e0a0142 100644 (file)
@@ -236,6 +236,10 @@ networks might depend on accelerations only available in the vendor's
 library. LAM-MPI_ might work, but since it has
 been deprecated for years, it is not supported.
 
+For example, depending on your actual MPI library, use ``cmake
+-DCMAKE_C_COMPILER=mpicc -DCMAKE_CXX_COMPILER=mpicxx -DGMX_MPI=on``.
+
+
 CMake
 ^^^^^
 
index e881a0f2a2c8cd18a21b25f26299abc37a296706..abbf71ad1204f19e5f355a87ae5ca3327d519698 100644 (file)
@@ -1641,7 +1641,7 @@ Trotter expansion
 .. math:: \begin{aligned}
           \exp(iL{\Delta t}) &\approx& \exp\left(iL_{\mathrm{NHC}}{\Delta t}/2\right) \exp\left(iL_2 {\Delta t}/2\right) \nonumber \\
           &&\exp\left(iL_1 {\Delta t}\right) \exp\left(iL_2 {\Delta t}/2\right) \exp\left(iL_{\mathrm{NHC}}{\Delta t}/2\right).\end{aligned}
-          :label:eqnVVNHTrotter
+          :label: eqnVVNHTrotter
 
 If the Nosé-Hoover chain is sufficiently slow with respect to the
 motions of the system, we can write an alternate integrator over
index 07df0ca9933a4f122cf8735b92517fc20c5f5650..34946c73a8dfaeef24a4e16c53f652bac6128608 100644 (file)
@@ -62,6 +62,15 @@ harmonically restrained to a plane or a line. Position restraints are
 applied to a special fixed list of atoms. Such a list is usually
 generated by the :ref:`pdb2gmx <gmx pdb2gmx>` program.
 
+.. _reference-manual-position-restraints:
+
+Note that position restraints make the potential dependent on absolute
+coordinates in space. Therefore, in general the pressure (and virial)
+is not well defined, as the pressure is the derivative of the free-energy
+of the system with respect to the volume. When the reference coordinates
+are scaled along with the system, which can be selected with the mdp option
+:mdp-value:`refcoord-scaling=all`, the pressure and virial are well defined.
+
 Flat-bottomed position restraints
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
@@ -190,9 +199,9 @@ as:
 
 .. math:: V_{dihr}(\phi') ~=~ \left\{
           \begin{array}{lcllll}
-          {\frac{1}{2}}k_{dihr}(\phi'-\phi_0-\Delta\phi)^2      
-                          &\mbox{for}&     \phi' & >   & \Delta\phi       \\[1.5ex]
-          0               &\mbox{for}&     \phi' & \le & \Delta\phi       \\[1.5ex]
+          {\frac{1}{2}}k_{dihr}(\phi'-\Delta\phi)^2      
+                          &\mbox{for}&     \|\phi'\| & >   & \Delta\phi       \\[1.5ex]
+          0               &\mbox{for}&     \|\phi'\| & \le & \Delta\phi       \\[1.5ex]
           \end{array}\right.
           :label: eqndihre
 
index ec95099ea587961130cdc416f1ad7b48bdf8bf4b..85b14b943a40ff712203a0cc0a670e8549a388df 100644 (file)
@@ -2538,6 +2538,21 @@ molecular dynamics simulations*, (2002).
 
    </div>
 
+.. raw:: html
+
+   <div id="ref-GroenhofEwaldArtefact">
+
+.. _refGroenhofEwaldArtefact:
+
+:sup:`181` Hub, J. S., de Groot, B. L., Grubmüller, H., Groenhof, G.,
+"Quantifying artifacts in Ewald simulations of inhomogeneous systems with a net charge,"
+*J. Chem. Theory Comput.*, **10**, 381–390 (2014).
+
+.. raw:: html
+
+   </div>
+
+
 .. raw:: html
 
    </div>
index a2ee52c24c23bfc886406bd6273ca160b0ec723e..5ec2d1c00c501946fb0fc9b18728025dc67d6814 100644 (file)
@@ -1,19 +1,54 @@
 GROMACS 2018.6 release notes
 ----------------------------
 
-This version was released on TODO, 2019. These release notes document
+This version was released on February 22, 2019. These release notes document
 the changes that have taken place in GROMACS since version 2018.5, to fix known
-issues. It also incorporates all fixes made in version TODO and
-earlier, which you can find described in the :ref:`release-notes`.
+issues. It also incorporates all fixes made in previous versions,
+which you can find described in the :ref:`release-notes`.
 
 Fixes where mdrun could behave incorrectly
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
+Correct free-energy Delta H output with mass lambda's
+"""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+When separate lambda parameters were used for perturbed mass
+free-energy contributions, these contributions were double counted
+in the Delta H output used for BAR calculations. Note that dH/dlambda
+was always correct
+
+:issue:`2703`
+:issue:`2849`
+
+.. _release-notes-2018-6-gpu:
+
+Fix incorrect LJ repulsion force switching on GPUs
+""""""""""""""""""""""""""""""""""""""""""""""""""
+
+When using a CUDA or OpenCL GPU, the coefficient for the second order
+term for the LJ repulsion in the force switching function, called 'A'
+in the manual, had the wrong sign. This lead to very small errors in
+the forces and the pressure. Note that the dispersion force switching
+was correct. The effects of this bug on any physical results seems to
+be negligible. Note that force switching is usually only used in
+combination with the CHARMM force field.
+
+:issue:`2845`
+
 Fixes for ``gmx`` tools
 ^^^^^^^^^^^^^^^^^^^^^^^
 
 Fixes to improve portability
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
+Fix compiler flags for Power8
+""""""""""""""""""""""""""""""""""""""""""""""""""
+
+A compiler flag for Power8 processors lead to errors in the code and was removed.
+
+:issue:`2747`
+:issue:`2746`
+:issue:`2734`
+
 Miscellaneous
 ^^^^^^^^^^^^^
diff --git a/docs/release-notes/2018/2018.7.rst b/docs/release-notes/2018/2018.7.rst
new file mode 100644 (file)
index 0000000..1545234
--- /dev/null
@@ -0,0 +1,44 @@
+GROMACS 2018.7 release notes
+----------------------------
+
+This version was released on May 29, 2019. These release notes document
+the changes that have taken place in GROMACS since version 2018.6, to fix known
+issues. It also incorporates all fixes made in previous versions,
+which you can find described in the :ref:`release-notes`.
+
+Fixes where mdrun could behave incorrectly
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Reverted broken change made in 2018.1
+"""""""""""""""""""""""""""""""""""""""""""""""""
+
+Reverted a change made in 2018.1 that broke simulations that used the
+SHAKE constraint algorithm.
+
+:issue:`2879`
+
+Work around gcc 7 AVX512 compiler bug
+"""""""""""""""""""""""""""""""""""""""
+
+With gcc version 7 a compiler bug caused a large part of non-bonded
+interactions to be ignored when compiled with AVX512 and running on more
+than 16 OpenMP threads.
+
+:issue:`2762`
+
+Fixes for ``gmx`` tools
+^^^^^^^^^^^^^^^^^^^^^^^
+
+Fixes to improve portability
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Miscellaneous
+^^^^^^^^^^^^^
+
+Updated release notes for 2018.6
+""""""""""""""""""""""""""""""""
+
+A :ref:`fix <release-notes-2018-6-gpu>` made to GPU kernels in 2018.6 was
+thought to resolve :issue:`2845` but further investigation suggests that
+the real cause is not yet known.
+
index b6972117e832f3635ea1dc44f4eb5e803b91bbc0..1f105605ff89886d97fb38e458852d4bf748501f 100644 (file)
@@ -19,6 +19,8 @@ has shifted too much .." when a cell size was limited.
 
 :issue:`2830`
 
+.. _release-notes-2019-1-gpu:
+
 Fix incorrect LJ repulsion force switching on GPUs
 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
 
@@ -26,9 +28,9 @@ When using a CUDA or OpenCL GPU, the coefficient for the second order
 term for the LJ repulsion in the force switching function, called 'A'
 in the manual, had the wrong sign. This lead to very small errors in
 the forces and the pressure. Note that the dispersion force switching
-was correct. Although the effect on individual atoms pairs was negligible,
-their combined effect on the pressure could lead to deformation of
-CHARMM membrane systems, where LJ force switching is regularly applied.
+was correct. The effects of this bug on any physical results seems to
+be negligible. Note that force switching is usually only used in
+combination with the CHARMM force field.
 
 :issue:`2845`
 
index 03f8e029f797b20d9aa415a7ede3b7ffdbd178d2..0f867ce1c061ca3a00ef45876e5ae26780416084 100644 (file)
@@ -1,10 +1,10 @@
 GROMACS 2019.2 release notes
 ----------------------------
 
-This version was released on TODO, 2019. These release notes
+This version was released on April 16th, 2019. These release notes
 document the changes that have taken place in GROMACS since the
-initial version 2019.1, to fix known issues. It also incorporates all
-fixes made in version 2018.5 and earlier, which you can find described
+previous 2019.1 version, to fix known issues. It also incorporates all
+fixes made in version 2018.6 and earlier, which you can find described
 in the :ref:`release-notes`.
 
 .. Note to developers!
@@ -16,12 +16,72 @@ in the :ref:`release-notes`.
 Fixes where mdrun could behave incorrectly
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
+Fix L-BGFS minimizer
+""""""""""""""""""""""""""""""""""""""""""""""""
+
+The minimizer could fail on a number of systems.
+
+:issue:`2641`
+
+Disallow pull geometry direction-periodic with AWH
+""""""""""""""""""""""""""""""""""""""""""""""""""
+
+This could lead to incorrect behavior or a cryptic error message.
+
+:issue:`2923`
+
+Fixed mdrun -nsteps option
+""""""""""""""""""""""""""
+
+Fixed that the, deprecated, mdrun option -nsteps only allowed extension
+of the simulation under certain conditions.
+
+:issue:`2881`
+
 Fixes for ``gmx`` tools
 ^^^^^^^^^^^^^^^^^^^^^^^
 
+gmx cluster -clndx indices now correct
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+The reported indices of trajectory frames in clusters were
+too small by one.
+
+:issue:`2926`
+
+gmx editconf -f in.pdb -o out.pdb again preserves chain IDs
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+This had been inadvertently broken and is now fixed.
+
+:issue:`2900`
+
+
+Tools again accept .tpr files as input
+"""""""""""""""""""""""""""""""""""""""
+
+The pdb2gmx, solvate, and insert-molecules tools could no longer
+accept input configurations contained in .tpr format files. This
+is now fixed.
+
+:issue:`2900`
+
+Fix segmentation fault when preparing simulated annealing inputs
+""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+grompp was unable to prepare tpr files for inputs containing simulated annealing
+procedures. The code has been fixed to allow the generation of those files again.
+
+:issue:`2871`
+       
 Fixes that affect portability
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
+Fix error in AVX 512 detection code
+"""""""""""""""""""""""""""""""""""
+
+The CMake detection code had a typo that could lead to wrong detection results.
+
 Miscellaneous
 ^^^^^^^^^^^^^
 
@@ -35,3 +95,25 @@ physical properties, such as the density, might be off from the intended values.
 
 :issue:`2884`
 
+Prevented internal build of FFTW with clang and AVX-512 SIMD
+""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+Prevented the internal build of FFTW with clang from attempting to
+configure FFTW to compile with AVX-512 support. That SIMD level is not
+supported by FFTW with the clang compiler, and compilation fails.
+
+:issue:`2892`
+
+Updated performance guide for recent Intel processors with AVX512 instruction support
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+Noted the tradeoffs between CPU frequency and SIMD throughput and advising users to
+prefer AVX2 over AVX512 in GPU-offload or highly parallel MPI cases.
+
+Updated release notes for 2019.1
+""""""""""""""""""""""""""""""""
+
+A :ref:`fix <release-notes-2019-1-gpu>` made to GPU kernels in 2019.1 was
+thought to resolve :issue:`2845` but further investigation suggests that
+the real cause is not yet known. 
+
diff --git a/docs/release-notes/2019/2019.3.rst b/docs/release-notes/2019/2019.3.rst
new file mode 100644 (file)
index 0000000..c5e7b7c
--- /dev/null
@@ -0,0 +1,118 @@
+GROMACS 2019.3 release notes
+----------------------------
+
+This version was released on June 14, 2019. These release notes
+document the changes that have taken place in GROMACS since the
+previous 2019.2 version, to fix known issues. It also incorporates all
+fixes made in version 2018.7 and earlier, which you can find described
+in the :ref:`release-notes`.
+
+.. Note to developers!
+   Please use """"""" to underline the individual entries for fixed issues in the subfolders,
+   otherwise the formatting on the webpage is messed up.
+   Also, please use the syntax :issue:`number` to reference issues on redmine, without the
+   a space between the colon and number!
+
+Fixes where mdrun could behave incorrectly
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Fix missing interactions with domain decomposition
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+When running with domain decomposition, any interactions described by
+the rarely-used topology file section
+``[ intermolecular_interactions ]`` were ignored. This did not
+affect normal non-bonded or intra-molecular interactions.
+
+:issue:`2953`
+
+Fix possible floating point exception during minimization.
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+It was possible that very small forces during minimization could lead to
+a crash due to a divide by zero error. Fixed by introducing a check.
+
+:issue:`2917`
+
+Fix segmentation fault when using membrane embedding
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+:issue:`2947`
+
+Allow AWH with pull-geometry 'direction' to be periodic
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+When applying AWH to a pull coordinate with geometry 'direction'
+with a AWH interval length of more than 95% of the box size,
+the dimension is now made periodic.
+
+:issue:`2946`
+       
+Fixes for ``gmx`` tools
+^^^^^^^^^^^^^^^^^^^^^^^
+
+Fixed residue and molecule indexing in selections
+"""""""""""""""""""""""""""""""""""""""""""""""""
+
+:issue:`2951`
+
+Fix PQR formatting
+""""""""""""""""""""
+
+The formatting was incorrect for some tools that use PQR files.
+
+:issue:`2955`
+
+Fix gmx wham with angle geometries
+""""""""""""""""""""""""""""""""""
+
+gmx wham would mix up degree and radian units leading to no overlap
+or not-a-number output.
+
+:issue:`2609`
+
+Add some information for grompp error with wrong line endings
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+Give meaningful error with too large grid in hbond
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+When using a grid that is too large :ref:`hbond <gmx hbond>` could try to
+allocate enough memory to cause a crash.
+
+:issue:`2962`
+
+Add some information for syntax errors with include delimiters in grompp
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+:issue:`2911`
+
+Fixes that affect portability
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Fixed wider reference SIMD setups
+"""""""""""""""""""""""""""""""""
+
+The reference SIMD builds could use a too small memory alignment,
+leading to mdrun exiting with an alignment error
+
+:issue:`2952`
+
+Fixed build failure with Apple Clang
+""""""""""""""""""""""""""""""""""""
+
+Builds would fail because of qsort being undefined.
+
+Miscellaneous
+^^^^^^^^^^^^^
+
+Removed non-existent mdp option awh1-dim1-period from user guide
+""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+:issue:`2940`
+
+Add checks for too many interactions during memory allocation
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+:issue:`2932`
+
diff --git a/docs/release-notes/2019/2019.4.rst b/docs/release-notes/2019/2019.4.rst
new file mode 100644 (file)
index 0000000..caa2a0a
--- /dev/null
@@ -0,0 +1,58 @@
+GROMACS 2019.4 release notes
+----------------------------
+
+This version was released on TODO, 2019. These release notes
+document the changes that have taken place in GROMACS since the
+previous 2019.3 version, to fix known issues. It also incorporates all
+fixes made in version 2018.7 and earlier, which you can find described
+in the :ref:`release-notes`.
+
+.. Note to developers!
+   Please use """"""" to underline the individual entries for fixed issues in the subfolders,
+   otherwise the formatting on the webpage is messed up.
+   Also, please use the syntax :issue:`number` to reference issues on redmine, without the
+   a space between the colon and number!
+
+Fixes where mdrun could behave incorrectly
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Fix segmentation fault in grompp and mdrun with cosine COM pulling
+""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+:issue:`3023`
+
+
+Fixes for ``gmx`` tools
+^^^^^^^^^^^^^^^^^^^^^^^
+
+Fix bug in gmx xpm2ps
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+The tool would fail when not being provided with a library file to read in.
+
+:issue:`3012`
+
+
+Fix bug in gmx anaeig
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+An issue was noted when reading a second set
+set of eigenvectors that could lead to problems when the number
+of eigenvectors was less than the three times the number of atoms.
+
+:issue:`2972`
+
+Fixes that affect portability
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Disable PME OpenCL on Apple
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+The Apple OpenCL compilers fail to produce a functional clFFT build.
+The OpenCL PME support is therefore disabled on Apple platforms.
+
+:issue:`2941`
+
+Miscellaneous
+^^^^^^^^^^^^^
+
index 4d7f1ec9ef4dc79412bca1b344c3261eb0dccea9..3ed1cf3c45354ee65c9d0b47e1ccdec45b6bec3b 100644 (file)
@@ -49,6 +49,8 @@ Patch releases
 .. toctree::
    :maxdepth: 1
 
+   2019/2019.4
+   2019/2019.3
    2019/2019.2
    2019/2019.1
 
@@ -78,6 +80,7 @@ Patch releases
 .. toctree::
    :maxdepth: 1
 
+   2018/2018.7
    2018/2018.6
    2018/2018.5
    2018/2018.4
index 15b8b8aaa988f2b0abd9c4f6a71797f6effbdc59..8ae3921743ec42be9b1014437929645d4614fefb 100644 (file)
@@ -1217,8 +1217,8 @@ Pressure coupling
 
       The reference coordinates for position restraints are not
       modified. Note that with this option the virial and pressure
-      will depend on the absolute positions of the reference
-      coordinates.
+      might be ill defined, see :ref:`here <reference-manual-position-restraints>`
+      for more details.
 
    .. mdp-value:: all
 
@@ -1233,7 +1233,9 @@ Pressure coupling
       one COM is used, even when there are multiple molecules with
       position restraints. For calculating the COM of the reference
       coordinates in the starting configuration, periodic boundary
-      conditions are not taken into account.
+      conditions are not taken into account. Note that with this option
+      the virial and pressure might be ill defined, see
+      :ref:`here <reference-manual-position-restraints>` for more details.
 
 
 Simulated annealing
@@ -1907,7 +1909,7 @@ AWH adaptive biasing
       multidimensional and is defined by mapping each dimension to a pull coordinate index.
       This is only allowed if :mdp-value:`pull-coord1-type=external-potential` and
       :mdp:`pull-coord1-potential-provider` = ``awh`` for the concerned pull coordinate
-      indices.
+      indices. Pull geometry 'direction-periodic' is not supported by AWH.
 
 .. mdp:: awh-potential
 
@@ -2131,19 +2133,18 @@ AWH adaptive biasing
    (0.0) [nm] or [rad]
    Start value of the sampling interval along this dimension. The range of allowed
    values depends on the relevant pull geometry (see :mdp:`pull-coord1-geometry`).
-   For periodic geometries :mdp:`awh1-dim1-start` greater than :mdp:`awh1-dim1-end`
+   For dihedral geometries :mdp:`awh1-dim1-start` greater than :mdp:`awh1-dim1-end`
    is allowed. The interval will then wrap around from +period/2 to -period/2.
+   For the direction geometry, the dimension is made periodic when
+   the direction is along a box vector and covers more than 95%
+   of the box length. Note that one should not apply pressure coupling
+   along a periodic dimension.
 
 .. mdp:: awh1-dim1-end
 
    (0.0) [nm] or [rad]
    End value defining the sampling interval together with :mdp:`awh1-dim1-start`.
 
-.. mdp:: awh1-dim1-period
-
-   (0.0) [nm] or [rad]
-   The period of this reaction coordinate, use 0 when the coordinate is not periodic.
-
 .. mdp:: awh1-dim1-diffusion
 
    (10\ :sup:`-5`) [nm\ :sup:`2`/ps] or [rad\ :sup:`2`/ps]
index 20311c0476d264e604358c02b49660e7ca06303b..6ee02260f8e693bc94881d7f85b7c01cfced3be8 100644 (file)
@@ -173,7 +173,9 @@ Parallelization schemes
 There are multiple parallelization schemes available, therefore a simulation can be run on a
 given hardware with different choices of run configuration.
 
-Core level parallelization via SIMD: SSE, AVX, etc.
+.. _intra-core-parallelization:
+
+Intra-core parallelization via SIMD: SSE, AVX, etc.
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
 One level of performance improvement available in |Gromacs| is through the use of
@@ -196,13 +198,20 @@ Thus, you need to configure and compile |Gromacs| for the SIMD capabilities of t
 By default, the build system will detect the highest supported
 acceleration of the host where the compilation is carried out. For cross-compiling for
 a machine with a different highest SIMD instructions set, in order to set the target acceleration,
-the ``-DGMX_SIMD`` CMake option can be used. For best performance always pick the highest
-(latest) SIMD instruction set supported by the target architecture (and |Gromacs|). To use a single
+the ``-DGMX_SIMD`` CMake option can be used.
+To use a single
 installation on multiple different machines, it is convenient to compile the analysis tools with
 the lowest common SIMD instruction set (as these rely little on SIMD acceleration), but for best
-performance :ref:`mdrun <gmx mdrun>` should be compiled separately for each machine.
+performance :ref:`mdrun <gmx mdrun>` should be compiled be compiled separately with the
+highest (latest) ``native`` SIMD instruction set of the target architecture (supported by |Gromacs|).
 
-.. TODO add a note on AVX throttle and its impact on MPI-parallel and GPU accelerated runs
+Recent Intel CPU architectures bring tradeoffs between the maximum clock frequency of the
+CPU (ie. its speed), and the width of the SIMD instructions it executes (ie its throughput
+at a given speed). In particular, the Intel ``Skylake`` and ``Cascade Lake`` processors
+(e.g. Xeon SP Gold/Platinum), can offer better throughput when using narrower SIMD because
+of the better clock frequency available. Consider building :ref:`mdrun <gmx mdrun>`
+configured with ``GMX_SIMD=AVX2_256`` instead of ``GMX_SIMD=AVX512`` for better
+performance in GPU accelerated or highly parallel MPI runs.
 
 Process(-or) level parallelization via OpenMP
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
@@ -439,6 +448,10 @@ behavior.
     default, 0, will start one thread on each available core.
     Alternatively, :ref:`mdrun <gmx mdrun>` will honor the appropriate system
     environment variable (e.g. ``OMP_NUM_THREADS``) if set.
+    Note that the maximum number of OpenMP threads (per rank) is,
+    for efficiency reasons, limited to 64. While it is rarely beneficial to use
+    a number of threads higher than this, the GMX_OPENMP_MAX_THREADS CMake variable
+    can be used to increase the limit.
 
 ``-npme``
     The total number of ranks to dedicate to the long-ranged
@@ -1060,7 +1073,7 @@ Known limitations
 
 **Please note again the limitations outlined below!**
 
-- Only compilation with CUDA is supported.
+- PME GPU offload is supported on NVIDIA hardware with CUDA and AMD hardware with OpenCL.
 
 - Only a PME order of 4 is supported on GPUs.
 
@@ -1295,9 +1308,13 @@ of 2. So it can be useful go through the checklist.
 * If you have GPUs that support either CUDA or OpenCL, use them.
 
   * Configure with ``-DGMX_GPU=ON`` (add ``-DGMX_USE_OPENCL=ON`` for OpenCL).
-  * For CUDA, use the newest CUDA availabe for your GPU to take advantage of the
+  * For CUDA, use the newest CUDA available for your GPU to take advantage of the
     latest performance enhancements.
   * Use a recent GPU driver.
+  * Make sure you use an :ref:`gmx mdrun` with ``GMX_SIMD`` appropriate for the CPU
+    architecture; the log file will contain a warning note if suboptimal setting is used.
+    However, prefer ``AVX2` over ``AVX512`` in GPU or highly parallel MPI runs (for more
+    information see the :ref:`intra-core parallelization information <intra-core-parallelization>`).
   * If compiling on a cluster head node, make sure that ``GMX_SIMD``
     is appropriate for the compute nodes.
 
index 58f5a0dd8deb7573ee5b141cd6a4bec6b30cd6ec..b47083a4a72028f27c671ad9dca98358492392f8 100644 (file)
@@ -272,11 +272,8 @@ If the charge is already close to an integer, then the difference is caused by
 Note for PME users: It is possible to use a uniform neutralizing background
 charge in PME to compensate for a system with a net background charge.
 This may however, especially for non-homogeneous systems, lead to unwanted artifacts, as
-shown in `Hub, J. S., de Groot, B. L., Grubmüller, H. & Groenhof, G. Quantifying
-artifacts in Ewald simulations of inhomogeneous systems with a net charge.
-*J. Chem. Theory Comput.* **10**, 381–390 (2014) <http://pubs.acs.org/doi/abs/10.1021/ct400626b>`.
-Nevertheless, it is a standard
-practice to actually add counter-ions to make the system net neutral.
+shown in \ :ref:`181 <refGroenhofEwaldArtefact>` (http://pubs.acs.org/doi/abs/10.1021/ct400626b).
+Nevertheless, it is a standard practice to actually add counter-ions to make the system net neutral.
 
 Incorrect number of parameters
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
index 5999da484d0fdda6d89391de3d93c405f3d82157..a1c3343f10538d9c486b1866b1ce43f9b8f72607 100644 (file)
@@ -62,7 +62,7 @@ simulations. Some stages are optional for some kinds of simulations.
    favourite text editor in concert with chapter 5 of the |Gromacs|
    `Reference Manual`_. For the AMBER force fields, `antechamber
    <http://amber.scripps.edu/antechamber/antechamber.html>`__ or
-   `acpype <https://github.com/choderalab/mmtools/blob/master/converters/acpype.py>`__
+   `acpype <https://github.com/alanwilter/acpype>`__
    might be appropriate.
 
 6. Describe a simulation box (e.g. using :ref:`gmx editconf`) whose
index ddc727001d7f23fc4e6f28f92911c9b98f6f3435..70b05dfa19d52b6acdc34b18aa9ca4863fa08441 100644 (file)
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2018, by the GROMACS development team, led by
+ * Copyright (c) 2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -89,7 +89,7 @@ class GmxApiTest : public gmx::test::MdrunTestFixture
          */
         void makeTprFile(int steps)
         {
-            runner_.useTopGroAndNdxFromDatabase("spc-and-methanol");
+            runner_.useTopGroAndNdxFromDatabase("spc_and_methane");
             runner_.useStringAsMdpFile(gmx::formatString("integrator = md\n"
                                                          "cutoff-scheme = Verlet\n"
                                                          "nsteps = %d\n"
index aada5438d667487825ea941c7010c8909824b678..a0a1379a0124f0e41e0f42a6bb10c4f95e7b0fe9 100644 (file)
@@ -80,9 +80,11 @@ elseif(${GMX_SIMD_ACTIVE} MATCHES "^(AVX)" AND NOT ${GMX_SIMD_ACTIVE} MATCHES "^
     set(_fftw_simd_support_level --enable-sse2;--enable-avx;--enable-avx2)
 elseif(${GMX_SIMD_ACTIVE} MATCHES "^(AVX_512)")
     # MSVC, GCC < 4.9, Clang < 3.9 do not support AVX-512, so
-    # we should not enable it.
+    # we should not enable it there. FFTW does not support clang with
+    # AVX-512, so we should not enable that either.
     if(MSVC OR (CMAKE_COMPILER_IS_GNUCC AND CMAKE_C_COMPILER_VERSION VERSION_LESS 4.9.0) OR
-       (CMAKE_C_COMPILER_ID MATCHES "Clang" AND CMAKE_C_COMPILER_VERSION VERSION_LESS 3.9.0))
+        (CMAKE_C_COMPILER_ID MATCHES "Clang" AND CMAKE_C_COMPILER_VERSION VERSION_LESS 3.9.0) OR
+        (CMAKE_C_COMPILER_ID MATCHES "Clang" AND ${GMX_SIMD_ACTIVE} MATCHES "^(AVX_512)"))
         set(_fftw_simd_support_level --enable-sse2;--enable-avx;--enable-avx2)
     else()
         set(_fftw_simd_support_level --enable-sse2;--enable-avx;--enable-avx2;--enable-avx512)
index e53ac5691c6b9a360993f76568fec8457bd46982..aa67b48ddd48bbaf81349b0cbe035107dc535034 100755 (executable)
@@ -30,7 +30,7 @@ distribution.
 #   include <cstddef>
 #endif
 
-#if defined(__GNUC__) && __GNUC__>=7
+#if defined(__GNUC__) && !defined(__INTEL_COMPILER) && __GNUC__>=7
 #pragma GCC diagnostic ignored "-Wimplicit-fallthrough"
 #endif
 
index 09ee0bcadd4b8763fb24ec770dff93ed639321a2..d51f4646e18325b611c3fc373264832d2fae20d4 100644 (file)
@@ -91,7 +91,7 @@
 #  endif
 #endif
 
-#if defined(__GNUC__) && __GNUC__>=7
+#if defined(__GNUC__) && !defined(__INTEL_COMPILER) && __GNUC__>=7
 #pragma GCC diagnostic ignored "-Wimplicit-fallthrough"
 #endif
 
index cf3a1576793ec2d810a92bd29274a0ee93e27c03..3d035c3ff8ad6eab3ed8f9bfef30e72b79b43093 100644 (file)
@@ -163,6 +163,7 @@ Awh::Awh(FILE                 *fplog,
             const AwhDimParams &awhDimParams      = awhBiasParams.dimParams[d];
             GMX_RELEASE_ASSERT(awhDimParams.eCoordProvider == eawhcoordproviderPULL, "Currently only the pull code is supported as coordinate provider");
             const t_pull_coord &pullCoord         = inputRecord.pull->coord[awhDimParams.coordIndex];
+            GMX_RELEASE_ASSERT(pullCoord.eGeom != epullgDIRPBC, "Pull geometry 'direction-periodic' is not supported by AWH");
             double              conversionFactor  = pull_coordinate_is_angletype(&pullCoord) ? DEG2RAD : 1;
             dimParams.emplace_back(conversionFactor, awhDimParams.forceConstant, beta);
 
index a02b54fc4ac00165174f2c3d5dca001ac7771847..ea4e6e219875aad3218aa75f4ca43a56dda40831 100644 (file)
@@ -553,77 +553,45 @@ AwhParams *readAndCheckAwhParams(std::vector<t_inpfile> *inp, const t_inputrec *
 /*! \brief
  * Gets the period of a pull coordinate.
  *
- * \param[in] pull_params      Pull parameters.
- * \param[in] coord_ind        Pull coordinate index.
- * \param[in] box              Box vectors.
+ * \param[in] pullCoordParams   The parameters for the pull coordinate.
+ * \param[in] pbc               The PBC setup
+ * \param[in] intervalLength    The length of the AWH interval for this pull coordinate
  * \returns the period (or 0 if not periodic).
  */
-static double get_pull_coord_period(const pull_params_t *pull_params,
-                                    int                  coord_ind,
-                                    const matrix         box)
+static double get_pull_coord_period(const t_pull_coord &pullCoordParams,
+                                    const t_pbc        &pbc,
+                                    const real          intervalLength)
 {
-    double        period;
-    t_pull_coord *pcrd_params = &pull_params->coord[coord_ind];
-
-    if (pcrd_params->eGeom == epullgDIRPBC)
-    {
-        /* For direction periodic, we need the pull vector to be one of the box vectors
-           (or more generally I guess it could be an integer combination of boxvectors).
-           This boxvector should to be orthogonal to the (periodic) plane spanned by the other two box vectors.
-           Here we assume that the pull vector is either x, y or z.
-         * E.g. for pull vec = (1, 0, 0) the box vector tensor should look like:
-         * | x 0 0 |
-         * | 0 a c |
-         * | 0 b d |
-         *
-           The period is then given by the box length x.
-
-           Note: we make these checks here for AWH and not in pull because we allow pull to be more general.
-         */
-        int m_pullvec = -1, count_nonzeros = 0;
-
-        /* Check that pull vec has only one component and which component it is. This component gives the relevant box vector */
-        for (int m = 0; m < DIM; m++)
-        {
-            if (pcrd_params->vec[m] != 0)
-            {
-                m_pullvec = m;
-                count_nonzeros++;
-            }
-        }
-        if (count_nonzeros != 1)
-        {
-            gmx_fatal(FARGS, "For AWH biasing pull coordinate %d with pull geometry %s, the pull vector needs to be parallel to "
-                      "a box vector that is parallel to either the x, y or z axis and is orthogonal to the other box vectors.",
-                      coord_ind + 1, EPULLGEOM(epullgDIRPBC));
-        }
+    double period = 0;
 
-        /* Check that there is a box vec parallel to pull vec and that this boxvec is orthogonal to the other box vectors */
-        for (int m = 0; m < DIM; m++)
+    if (pullCoordParams.eGeom == epullgDIR)
+    {
+        const real margin           = 0.001;
+        // Make dims periodic when the interval covers > 95%
+        const real periodicFraction = 0.95;
+
+        // Check if the pull direction is along a box vector
+        for (int dim = 0; dim < pbc.ndim_ePBC; dim++)
         {
-            for (int n = 0; n < DIM; n++)
+            const real boxLength    = norm(pbc.box[dim]);
+            const real innerProduct = iprod(pullCoordParams.vec, pbc.box[dim]);
+            if (innerProduct >= (1 - margin)*boxLength &&
+                innerProduct <= (1 + margin)*boxLength)
             {
-                if ((n != m) && (n == m_pullvec || m == m_pullvec) && box[m][n] > 0)
+                GMX_RELEASE_ASSERT(intervalLength < (1 + margin)*boxLength,
+                                   "We have checked before that interval <= period");
+                if (intervalLength > periodicFraction*boxLength)
                 {
-                    gmx_fatal(FARGS, "For AWH biasing pull coordinate %d with pull geometry %s, there needs to be a box vector parallel to the pull vector that is "
-                              "orthogonal to the other box vectors.",
-                              coord_ind + 1, EPULLGEOM(epullgDIRPBC));
+                    period = boxLength;
                 }
             }
         }
-
-        /* If this box vector only has one component as we assumed the norm should be equal to the absolute value of that component */
-        period = static_cast<double>(norm(box[m_pullvec]));
     }
-    else if (pcrd_params->eGeom == epullgDIHEDRAL)
+    else if (pullCoordParams.eGeom == epullgDIHEDRAL)
     {
         /* The dihedral angle is periodic in -180 to 180 deg */
         period = 360;
     }
-    else
-    {
-        period = 0;
-    }
 
     return period;
 }
@@ -733,7 +701,7 @@ static void checkInputConsistencyInterval(const AwhParams *awhParams, warninp_t
 
 void setStateDependentAwhParams(AwhParams *awhParams,
                                 const pull_params_t *pull_params, pull_t *pull_work,
-                                const matrix box,  int ePBC,
+                                const matrix box, int ePBC, const tensor &compressibility,
                                 const t_grpopts *inputrecGroupOptions, warninp_t wi)
 {
     /* The temperature is not really state depenendent but is not known
@@ -761,17 +729,46 @@ void setStateDependentAwhParams(AwhParams *awhParams,
         AwhBiasParams *awhBiasParams = &awhParams->awhBiasParams[k];
         for (int d = 0; d < awhBiasParams->ndim; d++)
         {
-            AwhDimParams *dimParams = &awhBiasParams->dimParams[d];
+            AwhDimParams       *dimParams       = &awhBiasParams->dimParams[d];
+            const t_pull_coord &pullCoordParams = pull_params->coord[dimParams->coordIndex];
 
-            /* The periodiciy of the AWH grid in certain cases depends on the simulation box */
-            dimParams->period = get_pull_coord_period(pull_params, dimParams->coordIndex, box);
+            if (pullCoordParams.eGeom == epullgDIRPBC)
+            {
+                gmx_fatal(FARGS, "AWH does not support pull geometry '%s'. "
+                          "If the maximum distance between the groups is always less than half the box size, "
+                          "you can use geometry '%s' instead.",
+                          EPULLGEOM(epullgDIRPBC),
+                          EPULLGEOM(epullgDIR));
+
+            }
+
+            dimParams->period = get_pull_coord_period(pullCoordParams, pbc, dimParams->end - dimParams->origin);
+            // We would like to check for scaling, but we don't have the full inputrec available here
+            if (dimParams->period > 0 && !(pullCoordParams.eGeom == epullgANGLE ||
+                                           pullCoordParams.eGeom == epullgDIHEDRAL))
+            {
+                bool coordIsScaled = false;
+                for (int d2 = 0; d2 < DIM; d2++)
+                {
+                    if (pullCoordParams.vec[d2] != 0 && norm2(compressibility[d2]) != 0)
+                    {
+                        coordIsScaled = true;
+                    }
+                }
+                if (coordIsScaled)
+                {
+                    std::string mesg = gmx::formatString("AWH dimension %d of bias %d is periodic with pull geometry '%s', "
+                                                         "while you should are applying pressure scaling to the corresponding box vector, this is not supported.",
+                                                         d + 1, k + 1, EPULLGEOM(pullCoordParams.eGeom));
+                    warning(wi, mesg.c_str());
+                }
+            }
 
             /* The initial coordinate value, converted to external user units. */
             dimParams->coordValueInit =
                 get_pull_coord_value(pull_work, dimParams->coordIndex, &pbc);
 
-            t_pull_coord *pullCoord = &pull_params->coord[dimParams->coordIndex];
-            dimParams->coordValueInit *= pull_conversion_factor_internal2userinput(pullCoord);
+            dimParams->coordValueInit *= pull_conversion_factor_internal2userinput(&pullCoordParams);
         }
     }
     checkInputConsistencyInterval(awhParams, wi);
index fa953974ea39ff96843b45a8c5556134d6190857..ef2f8eab06df98fe93567e300711818a65c37c0b 100644 (file)
@@ -78,6 +78,7 @@ AwhParams *readAndCheckAwhParams(std::vector<t_inpfile> *inp,
  * \param[in,out] pull_work             Pull working struct to register AWH bias in.
  * \param[in]     box                   Box vectors.
  * \param[in]     ePBC                  Periodic boundary conditions enum.
+ * \param[in]     compressibility       Compressibility matrix for pressure coupling, pass all 0 without pressure coupling
  * \param[in]     inputrecGroupOptions  Parameters for atom groups.
  * \param[in,out] wi                    Struct for bookeeping warnings.
  *
@@ -88,6 +89,7 @@ void setStateDependentAwhParams(AwhParams           *awhParams,
                                 pull_t              *pull_work,
                                 const matrix         box,
                                 int                  ePBC,
+                                const tensor        &compressibility,
                                 const t_grpopts     *inputrecGroupOptions,
                                 warninp_t            wi);
 
index 4ba6cc34846d2cf8ea4dd16d9ce07fe844799759..40488725e73536f38ca5c281fc9b237310276768 100644 (file)
@@ -1240,7 +1240,7 @@ static void make_load_communicators(gmx_domdec_t gmx_unused *dd)
         fprintf(debug, "Making load communicators\n");
     }
 
-    snew(dd->comm->load,          std::max(dd->ndim, 1));
+    dd->comm->load = new domdec_load_t[std::max(dd->ndim, 1)];
     snew(dd->comm->mpi_comm_load, std::max(dd->ndim, 1));
 
     if (dd->ndim == 0)
index 46a460ee1be0259af908a35019190739e905d13c..bd7b57f9aa91198bcb6942ef5333eb7dfbc2e154 100644 (file)
@@ -59,21 +59,27 @@ struct t_commrec;
 
 struct BalanceRegion;
 
+//! Indices to communicate in a dimension
 struct gmx_domdec_ind_t
 {
-    /* The numbers of charge groups to send and receive for each cell
-     * that requires communication, the last entry contains the total
+    //! @{
+    /*! \brief The numbers of charge groups to send and receive for each
+     * cell that requires communication, the last entry contains the total
      * number of atoms that needs to be communicated.
      */
-    int              nsend[DD_MAXIZONE+2];
-    int              nrecv[DD_MAXIZONE+2];
-    /* The charge groups to send */
+    int              nsend[DD_MAXIZONE+2] = {};
+    int              nrecv[DD_MAXIZONE+2] = {};
+    //! @}
+    //! The charge groups to send
     std::vector<int> index;
+    //! @{
     /* The atom range for non-in-place communication */
-    int              cell2at0[DD_MAXIZONE];
-    int              cell2at1[DD_MAXIZONE];
+    int              cell2at0[DD_MAXIZONE] = {};
+    int              cell2at1[DD_MAXIZONE] = {};
+    //! @}
 };
 
+//! Things relating to index communication
 struct gmx_domdec_comm_dim_t
 {
     /* Returns the number of grid pulses (the number of domains in the halo along this dimension) */
@@ -82,9 +88,12 @@ struct gmx_domdec_comm_dim_t
         return ind.size();
     }
 
-    int                           np_dlb;         /* For dlb, for use with edlbAUTO          */
-    std::vector<gmx_domdec_ind_t> ind;            /* The indices to communicate, size np     */
-    bool                          receiveInPlace; /* Can we receive data in place?            */
+    /**< For dlb, for use with edlbAUTO          */
+    int                           np_dlb = 0;
+    /**< The indices to communicate, size np     */
+    std::vector<gmx_domdec_ind_t> ind;
+    /**< Can we receive data in place?            */
+    bool                          receiveInPlace = false;
 };
 
 /*! \brief Load balancing data along a dim used on the master rank of that dim */
@@ -92,30 +101,45 @@ struct RowMaster
 {
     struct Bounds
     {
-        real cellFracLowerMax; /**< State var.: max lower bound., incl. neighbors */
-        real cellFracUpperMin; /**< State var.: min upper bound., incl. neighbors */
-        real boundMin;         /**< Temp. var.: lower limit for cell boundary */
-        real boundMax;         /**< Temp. var.: upper limit for cell boundary */
+        /**< State var.: max lower bound., incl. neighbors */
+        real cellFracLowerMax = 0;
+        /**< State var.: min upper bound., incl. neighbors */
+        real cellFracUpperMin = 0;
+        /**< Temp. var.: lower limit for cell boundary */
+        real boundMin = 0;
+        /**< Temp. var.: upper limit for cell boundary */
+        real boundMax = 0;
     };
 
-    std::vector<bool>   isCellMin;    /**< Temp. var.: is this cell size at the limit */
-    std::vector<real>   cellFrac;     /**< State var.: cell boundaries, box relative */
-    std::vector<real>   oldCellFrac;  /**< Temp. var.: old cell size */
-    std::vector<Bounds> bounds;       /**< Cell bounds */
-    bool                dlbIsLimited; /**< State var.: is DLB limited in this row */
-    std::vector<real>   buf_ncd;      /**< Temp. var.  */
+    /**< Temp. var.: is this cell size at the limit */
+    std::vector<bool>   isCellMin;
+    /**< State var.: cell boundaries, box relative */
+    std::vector<real>   cellFrac;
+    /**< Temp. var.: old cell size */
+    std::vector<real>   oldCellFrac;
+    /**< Cell bounds */
+    std::vector<Bounds> bounds;
+    /**< State var.: is DLB limited in this row */
+    bool                dlbIsLimited = false;
+    /**< Temp. var.  */
+    std::vector<real>   buf_ncd;
 };
 
 /*! \brief Struct for managing cell sizes with DLB along a dimension */
 struct DDCellsizesWithDlb
 {
-    /* Cell sizes for dynamic load balancing */
-    std::unique_ptr<RowMaster> rowMaster;    /**< Cell row root struct, only available on the first rank in a row */
-    std::vector<real>          fracRow;      /**< The cell sizes, in fractions, along a row, not available on the first rank in a row */
-    real                       fracLower;    /**< The lower corner, in fractions, in triclinic space */
-    real                       fracUpper;    /**< The upper corner, in fractions, in triclinic space */
-    real                       fracLowerMax; /**< The maximum lower corner among all our neighbors */
-    real                       fracUpperMin; /**< The minimum upper corner among all our neighbors */
+    /**< Cell row root struct, only available on the first rank in a row */
+    std::unique_ptr<RowMaster> rowMaster;
+    /**< The cell sizes, in fractions, along a row, not available on the first rank in a row */
+    std::vector<real>          fracRow;
+    /**< The lower corner, in fractions, in triclinic space */
+    real                       fracLower = 0;
+    /**< The upper corner, in fractions, in triclinic space */
+    real                       fracUpper = 0;
+    /**< The maximum lower corner among all our neighbors */
+    real                       fracLowerMax = 0;
+    /**< The minimum upper corner among all our neighbors */
+    real                       fracUpperMin = 0;
 };
 
 /*! \brief Struct for compute load commuication
@@ -125,32 +149,48 @@ struct DDCellsizesWithDlb
  */
 typedef struct
 {
-    int    nload;     /**< The number of load recordings */
-    float *load;      /**< Scan of the sum of load over dimensions */
-    float  sum;       /**< The sum of the load over the ranks up to our current dimension */
-    float  max;       /**< The maximum over the ranks contributing to \p sum */
-    float  sum_m;     /**< Like \p sum, but takes the maximum when the load balancing is limited */
-    float  cvol_min;  /**< Minimum cell volume, relative to the box */
-    float  mdf;       /**< The PP time during which PME can overlap */
-    float  pme;       /**< The PME-only rank load */
-    int    flags;     /**< Bit flags that tell if DLB was limited, per dimension */
+    /**< The number of load recordings */
+    int    nload = 0;
+    /**< Scan of the sum of load over dimensions */
+    float *load = nullptr;
+    /**< The sum of the load over the ranks up to our current dimension */
+    float  sum = 0;
+    /**< The maximum over the ranks contributing to \p sum */
+    float  max = 0;
+    /**< Like \p sum, but takes the maximum when the load balancing is limited */
+    float  sum_m = 0;
+    /**< Minimum cell volume, relative to the box */
+    float  cvol_min = 0;
+    /**< The PP time during which PME can overlap */
+    float  mdf = 0;
+    /**< The PME-only rank load */
+    float  pme = 0;
+    /**< Bit flags that tell if DLB was limited, per dimension */
+    int    flags = 0;
 } domdec_load_t;
 
 /*! \brief Data needed to sort an atom to the desired location in the local state */
 typedef struct
 {
-    int  nsc;     /**< Neighborsearch grid cell index */
-    int  ind_gl;  /**< Global atom/charge group index */
-    int  ind;     /**< Local atom/charge group index */
+    /**< Neighborsearch grid cell index */
+    int  nsc = 0;
+    /**< Global atom/charge group index */
+    int  ind_gl = 0;
+    /**< Local atom/charge group index */
+    int  ind = 0;
 } gmx_cgsort_t;
 
 /*! \brief Temporary buffers for sorting atoms */
 typedef struct
 {
-    std::vector<gmx_cgsort_t> sorted;     /**< Sorted array of indices */
-    std::vector<gmx_cgsort_t> stationary; /**< Array of stationary atom/charge group indices */
-    std::vector<gmx_cgsort_t> moved;      /**< Array of moved atom/charge group indices */
-    std::vector<int>          intBuffer;  /**< Integer buffer for sorting */
+    /**< Sorted array of indices */
+    std::vector<gmx_cgsort_t> sorted;
+    /**< Array of stationary atom/charge group indices */
+    std::vector<gmx_cgsort_t> stationary;
+    /**< Array of moved atom/charge group indices */
+    std::vector<gmx_cgsort_t> moved;
+    /**< Integer buffer for sorting */
+    std::vector<int>          intBuffer;
 } gmx_domdec_sort_t;
 
 /*! \brief Manages atom ranges and order for the local state atom vectors */
@@ -257,25 +297,40 @@ enum class DlbState
 /*! \brief The PME domain decomposition for one dimension */
 typedef struct
 {
-    int      dim;       /**< The dimension */
-    gmx_bool dim_match; /**< Tells if DD and PME dims match */
-    int      nslab;     /**< The number of PME ranks/domains in this dimension */
-    real    *slb_dim_f; /**< Cell sizes for determining the PME comm. with SLB */
-    int     *pp_min;    /**< The minimum pp node location, size nslab */
-    int     *pp_max;    /**< The maximum pp node location, size nslab */
-    int      maxshift;  /**< The maximum shift for coordinate redistribution in PME */
+    /**< The dimension */
+    int      dim = 0;
+    /**< Tells if DD and PME dims match */
+    gmx_bool dim_match = false;
+    /**< The number of PME ranks/domains in this dimension */
+    int      nslab = 0;
+    /**< Cell sizes for determining the PME comm. with SLB */
+    real    *slb_dim_f = nullptr;
+    /**< The minimum pp node location, size nslab */
+    int     *pp_min = nullptr;
+    /**< The maximum pp node location, size nslab */
+    int     *pp_max = nullptr;
+    /**< The maximum shift for coordinate redistribution in PME */
+    int      maxshift = 0;
 } gmx_ddpme_t;
 
 struct gmx_ddzone_t
 {
-    real min0;    /* The minimum bottom of this zone                        */
-    real max1;    /* The maximum top of this zone                           */
-    real min1;    /* The minimum top of this zone                           */
-    real mch0;    /* The maximum bottom communicaton height for this zone   */
-    real mch1;    /* The maximum top communicaton height for this zone      */
-    real p1_0;    /* The bottom value of the first cell in this zone        */
-    real p1_1;    /* The top value of the first cell in this zone           */
-    real dataSet; /* Bool disguised as a real, 1 when the above data has been set. 0 otherwise */
+    /**< The minimum bottom of this zone                        */
+    real min0 = 0;
+    /**< The maximum top of this zone                           */
+    real max1 = 0;
+    /**< The minimum top of this zone                           */
+    real min1 = 0;
+    /**< The maximum bottom communicaton height for this zone   */
+    real mch0 = 0;
+    /**< The maximum top communicaton height for this zone      */
+    real mch1 = 0;
+    /**< The bottom value of the first cell in this zone        */
+    real p1_0 = 0;
+    /**< The top value of the first cell in this zone           */
+    real p1_1 = 0;
+    /**< Bool disguised as a real, 1 when the above data has been set. 0 otherwise */
+    real dataSet = 0;
 };
 
 /*! \brief The number of reals in gmx_ddzone_t */
@@ -365,14 +420,19 @@ class DDBufferAccess
         gmx::ArrayRef<T>  buffer;    /**< The access to the memory buffer */
 };
 
-/*! brief Temporary buffer for setting up communiation over one pulse and all zones in the halo */
+/*! \brief Temporary buffer for setting up communiation over one pulse and all zones in the halo */
 struct dd_comm_setup_work_t
 {
-    std::vector<int>       localAtomGroupBuffer; /**< The local atom group indices to send */
-    std::vector<int>       atomGroupBuffer;      /**< Buffer for collecting the global atom group indices to send */
-    std::vector<gmx::RVec> positionBuffer;       /**< Buffer for collecting the atom group positions to send */
-    int                    nat;                  /**< The number of atoms contained in the atom groups to send */
-    int                    nsend_zone;           /**< The number of atom groups to send for the last zone */
+    /**< The local atom group indices to send */
+    std::vector<int>       localAtomGroupBuffer;
+    /**< Buffer for collecting the global atom group indices to send */
+    std::vector<int>       atomGroupBuffer;
+    /**< Buffer for collecting the atom group positions to send */
+    std::vector<gmx::RVec> positionBuffer;
+    /**< The number of atoms contained in the atom groups to send */
+    int                    nat = 0;
+    /**< The number of atom groups to send for the last zone */
+    int                    nsend_zone = 0;
 };
 
 /*! \brief Struct for domain decomposition communication
@@ -388,77 +448,107 @@ struct dd_comm_setup_work_t
 struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
 {
     /* PME and Cartesian communicator stuff */
-    int         npmedecompdim;     /**< The number of decomposition dimensions for PME, 0: no PME */
-    int         npmenodes;         /**< The number of ranks doing PME (PP/PME or only PME) */
-    int         npmenodes_x;       /**< The number of PME ranks/domains along x */
-    int         npmenodes_y;       /**< The number of PME ranks/domains along y */
-    gmx_bool    bCartesianPP_PME;  /**< Use Cartesian communication between PP and PME ranks */
-    ivec        ntot;              /**< Cartesian grid for combinted PP+PME ranks */
-    int         cartpmedim;        /**< The number of dimensions for the PME setup that are Cartesian */
-    int        *pmenodes;          /**< The PME ranks, size npmenodes */
-    int        *ddindex2simnodeid; /**< The Cartesian index to sim rank conversion, used with bCartesianPP_PME */
-    gmx_ddpme_t ddpme[2];          /**< The 1D or 2D PME domain decomposition setup */
+    /**< The number of decomposition dimensions for PME, 0: no PME */
+    int         npmedecompdim = 0;
+    /**< The number of ranks doing PME (PP/PME or only PME) */
+    int         npmenodes = 0;
+    /**< The number of PME ranks/domains along x */
+    int         npmenodes_x = 0;
+    /**< The number of PME ranks/domains along y */
+    int         npmenodes_y = 0;
+    /**< Use Cartesian communication between PP and PME ranks */
+    gmx_bool    bCartesianPP_PME = false;
+    /**< Cartesian grid for combinted PP+PME ranks */
+    ivec        ntot = { };
+    /**< The number of dimensions for the PME setup that are Cartesian */
+    int         cartpmedim = 0;
+    /**< The PME ranks, size npmenodes */
+    int        *pmenodes = nullptr;
+    /**< The Cartesian index to sim rank conversion, used with bCartesianPP_PME */
+    int        *ddindex2simnodeid = nullptr;
+    /**< The 1D or 2D PME domain decomposition setup */
+    gmx_ddpme_t ddpme[2];
 
     /* The DD particle-particle nodes only */
-    gmx_bool bCartesianPP;        /**< Use a Cartesian communicator for PP */
-    int     *ddindex2ddnodeid;    /**< The Cartesian index to DD rank conversion, used with bCartesianPP */
+    /**< Use a Cartesian communicator for PP */
+    gmx_bool bCartesianPP = false;
+    /**< The Cartesian index to DD rank conversion, used with bCartesianPP */
+    int     *ddindex2ddnodeid = nullptr;
 
     /* The DLB state, used for reloading old states, during e.g. EM */
-    t_block cgs_gl;               /**< The global charge groups, this defined the DD state (except for the DLB state) */
+    /**< The global charge groups, this defined the DD state (except for the DLB state) */
+    t_block cgs_gl = { };
 
     /* Charge group / atom sorting */
-    std::unique_ptr<gmx_domdec_sort_t> sort; /**< Data structure for cg/atom sorting */
+    /**< Data structure for cg/atom sorting */
+    std::unique_ptr<gmx_domdec_sort_t> sort;
 
     //! True when update groups are used
-    bool                                  useUpdateGroups;
+    bool                                  useUpdateGroups = false;
     //!  Update atom grouping for each molecule type
     std::vector<gmx::RangePartitioning>   updateGroupingPerMoleculetype;
     //! Centers of mass of local update groups
     std::unique_ptr<gmx::UpdateGroupsCog> updateGroupsCog;
 
     /* Are there charge groups? */
-    bool haveInterDomainBondeds;          /**< Are there inter-domain bonded interactions? */
-    bool haveInterDomainMultiBodyBondeds; /**< Are there inter-domain multi-body interactions? */
+    bool haveInterDomainBondeds          = false; /**< Are there inter-domain bonded interactions? */
+    bool haveInterDomainMultiBodyBondeds = false; /**< Are there inter-domain multi-body interactions? */
 
     /* Data for the optional bonded interaction atom communication range */
-    gmx_bool  bBondComm;          /**< Only communicate atoms beyond the non-bonded cut-off when they are involved in bonded interactions with non-local atoms */
-    t_blocka *cglink;             /**< Links between cg's through bonded interactions */
-    char     *bLocalCG;           /**< Local cg availability, TODO: remove when group scheme is removed */
+    /**< Only communicate atoms beyond the non-bonded cut-off when they are involved in bonded interactions with non-local atoms */
+    gmx_bool  bBondComm = false;
+    /**< Links between cg's through bonded interactions */
+    t_blocka *cglink = nullptr;
+    /**< Local cg availability, TODO: remove when group scheme is removed */
+    char     *bLocalCG = nullptr;
 
     /* The DLB state, possible values are defined above */
     DlbState dlbState;
     /* With dlbState=DlbState::offCanTurnOn, should we check if to DLB on at the next DD? */
-    gmx_bool bCheckWhetherToTurnDlbOn;
+    gmx_bool bCheckWhetherToTurnDlbOn = false;
     /* The first DD count since we are running without DLB */
     int      ddPartioningCountFirstDlbOff = 0;
 
     /* Cell sizes for static load balancing, first index cartesian */
-    real **slb_frac;
+    real **slb_frac = nullptr;
 
     /* The width of the communicated boundaries */
-    real     cutoff_mbody;        /**< Cut-off for multi-body interactions, also 2-body bonded when \p cutoff_mody > \p cutoff */
-    real     cutoff;              /**< Cut-off for non-bonded/2-body interactions */
-    rvec     cellsize_min;        /**< The minimum guaranteed cell-size, Cartesian indexing */
-    rvec     cellsize_min_dlb;    /**< The minimum guaranteed cell-size with dlb=auto */
-    real     cellsize_limit;      /**< The lower limit for the DD cell size with DLB */
-    gmx_bool bVacDLBNoLimit;      /**< Effectively no NB cut-off limit with DLB for systems without PBC? */
+    /**< Cut-off for multi-body interactions, also 2-body bonded when \p cutoff_mody > \p cutoff */
+    real     cutoff_mbody = 0;
+    /**< Cut-off for non-bonded/2-body interactions */
+    real     cutoff = 0;
+    /**< The minimum guaranteed cell-size, Cartesian indexing */
+    rvec     cellsize_min = { };
+    /**< The minimum guaranteed cell-size with dlb=auto */
+    rvec     cellsize_min_dlb = { };
+    /**< The lower limit for the DD cell size with DLB */
+    real     cellsize_limit = 0;
+    /**< Effectively no NB cut-off limit with DLB for systems without PBC? */
+    gmx_bool bVacDLBNoLimit = false;
 
     /** With PME load balancing we set limits on DLB */
-    gmx_bool bPMELoadBalDLBLimits;
+    gmx_bool bPMELoadBalDLBLimits = false;
     /** DLB needs to take into account that we want to allow this maximum
      *  cut-off (for PME load balancing), this could limit cell boundaries.
      */
-    real PMELoadBal_max_cutoff;
+    real PMELoadBal_max_cutoff = 0;
 
-    ivec tric_dir;                /**< tric_dir from \p gmx_ddbox_t is only stored here because dd_get_ns_ranges needs it */
-    rvec box0;                    /**< box lower corner, required with dim's without pbc when avoiding communication */
-    rvec box_size;                /**< box size, required with dim's without pbc when avoiding communication */
+    /**< tric_dir from \p gmx_ddbox_t is only stored here because dd_get_ns_ranges needs it */
+    ivec tric_dir = { };
+    /**< box lower corner, required with dim's without pbc and -gcom */
+    rvec box0 = { };
+    /**< box size, required with dim's without pbc and -gcom */
+    rvec box_size = { };
 
-    rvec cell_x0;                 /**< The DD cell lower corner, in triclinic space */
-    rvec cell_x1;                 /**< The DD cell upper corner, in triclinic space */
+    /**< The DD cell lower corner, in triclinic space */
+    rvec cell_x0 = { };
+    /**< The DD cell upper corner, in triclinic space */
+    rvec cell_x1 = { };
 
-    rvec old_cell_x0;             /**< The old \p cell_x0, to check cg displacements */
-    rvec old_cell_x1;             /**< The old \p cell_x1, to check cg displacements */
+    /**< The old \p cell_x0, to check cg displacements */
+    rvec old_cell_x0 = { };
+    /**< The old \p cell_x1, to check cg displacements */
+    rvec old_cell_x1 = { };
 
     /** The communication setup and charge group boundaries for the zones */
     gmx_domdec_zones_t zones;
@@ -467,21 +557,23 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
      * cell boundaries of neighboring cells for staggered grids when using
      * dynamic load balancing.
      */
-    gmx_ddzone_t zone_d1[2];          /**< Zone limits for dim 1 with staggered grids */
-    gmx_ddzone_t zone_d2[2][2];       /**< Zone limits for dim 2 with staggered grids */
+    /**< Zone limits for dim 1 with staggered grids */
+    gmx_ddzone_t zone_d1[2];
+    /**< Zone limits for dim 2 with staggered grids */
+    gmx_ddzone_t zone_d2[2][2];
 
     /** The coordinate/force communication setup and indices */
     gmx_domdec_comm_dim_t cd[DIM];
     /** The maximum number of cells to communicate with in one dimension */
-    int                   maxpulse;
+    int                   maxpulse = 0;
 
     /** Which cg distribution is stored on the master node,
      *  stored as DD partitioning call count.
      */
-    int64_t master_cg_ddp_count;
+    int64_t master_cg_ddp_count = 0;
 
     /** The number of cg's received from the direct neighbors */
-    int  zone_ncg1[DD_MAXZONE];
+    int  zone_ncg1[DD_MAXZONE] = {0};
 
     /** The atom ranges in the local state */
     DDAtomRanges atomRanges;
@@ -496,71 +588,103 @@ struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
     DDBuffer<gmx::RVec> rvecBuffer;
 
     /* Temporary storage for thread parallel communication setup */
-    std::vector<dd_comm_setup_work_t> dth; /**< Thread-local work data */
+    /**< Thread-local work data */
+    std::vector<dd_comm_setup_work_t> dth;
 
     /* Communication buffer only used with multiple grid pulses */
-    DDBuffer<gmx::RVec> rvecBuffer2; /**< Another rvec comm. buffer */
+    /**< Another rvec comm. buffer */
+    DDBuffer<gmx::RVec> rvecBuffer2;
 
     /* Communication buffers for local redistribution */
-    std::array<std::vector<int>, DIM*2>       cggl_flag;  /**< Charge group flag comm. buffers */
-    std::array<std::vector<gmx::RVec>, DIM*2> cgcm_state; /**< Charge group center comm. buffers */
+    /**< Charge group flag comm. buffers */
+    std::array<std::vector<int>, DIM*2>       cggl_flag;
+    /**< Charge group center comm. buffers */
+    std::array<std::vector<gmx::RVec>, DIM*2> cgcm_state;
 
     /* Cell sizes for dynamic load balancing */
     std::vector<DDCellsizesWithDlb> cellsizesWithDlb;
 
     /* Stuff for load communication */
-    gmx_bool        bRecordLoad;         /**< Should we record the load */
-    domdec_load_t  *load;                /**< The recorded load data */
-    int             nrank_gpu_shared;    /**< The number of MPI ranks sharing the GPU our rank is using */
+    /**< Should we record the load */
+    gmx_bool        bRecordLoad = false;
+    /**< The recorded load data */
+    domdec_load_t  *load = nullptr;
+    /**< The number of MPI ranks sharing the GPU our rank is using */
+    int             nrank_gpu_shared = 0;
 #if GMX_MPI
-    MPI_Comm       *mpi_comm_load;       /**< The MPI load communicator */
-    MPI_Comm        mpi_comm_gpu_shared; /**< The MPI load communicator for ranks sharing a GPU */
+    /**< The MPI load communicator */
+    MPI_Comm       *mpi_comm_load = nullptr;
+    /**< The MPI load communicator for ranks sharing a GPU */
+    MPI_Comm        mpi_comm_gpu_shared;
 #endif
 
     /* Information for managing the dynamic load balancing */
-    int            dlb_scale_lim;      /**< Maximum DLB scaling per load balancing step in percent */
+    /**< Maximum DLB scaling per load balancing step in percent */
+    int            dlb_scale_lim = 0;
 
-    BalanceRegion *balanceRegion;      /**< Struct for timing the force load balancing region */
+    /**< Struct for timing the force load balancing region */
+    BalanceRegion *balanceRegion = nullptr;
 
     /* Cycle counters over nstlist steps */
-    float  cycl[ddCyclNr];             /**< Total cycles counted */
-    int    cycl_n[ddCyclNr];           /**< The number of cycle recordings */
-    float  cycl_max[ddCyclNr];         /**< The maximum cycle count */
+    /**< Total cycles counted */
+    float  cycl[ddCyclNr] = { };
+    /**< The number of cycle recordings */
+    int    cycl_n[ddCyclNr] = { };
+    /**< The maximum cycle count */
+    float  cycl_max[ddCyclNr] = { };
     /** Flop counter (0=no,1=yes,2=with (eFlop-1)*5% noise */
-    int    eFlop;
-    double flop;                       /**< Total flops counted */
-    int    flop_n;                     /**< The number of flop recordings */
+    int    eFlop = 0;
+    /**< Total flops counted */
+    double flop = 0.0;
+    /**< The number of flop recordings */
+    int    flop_n = 0;
     /** How many times did we have load measurements */
-    int    n_load_have;
+    int    n_load_have = 0;
     /** How many times have we collected the load measurements */
-    int    n_load_collect;
+    int    n_load_collect = 0;
 
     /* Cycle count history for DLB checks */
-    float       cyclesPerStepBeforeDLB;     /**< The averaged cycles per step over the last nstlist step before turning on DLB */
-    float       cyclesPerStepDlbExpAverage; /**< The running average of the cycles per step during DLB */
-    bool        haveTurnedOffDlb;           /**< Have we turned off DLB (after turning DLB on)? */
-    int64_t     dlbSlowerPartitioningCount; /**< The DD step at which we last measured that DLB off was faster than DLB on, 0 if there was no such step */
+    /**< The averaged cycles per step over the last nstlist step before turning on DLB */
+    float       cyclesPerStepBeforeDLB = 0;
+    /**< The running average of the cycles per step during DLB */
+    float       cyclesPerStepDlbExpAverage = 0;
+    /**< Have we turned off DLB (after turning DLB on)? */
+    bool        haveTurnedOffDlb = false;
+    /**< The DD step at which we last measured that DLB off was faster than DLB on, 0 if there was no such step */
+    int64_t     dlbSlowerPartitioningCount = 0;
 
     /* Statistics for atoms */
-    double sum_nat[static_cast<int>(DDAtomRanges::Type::Number)]; /**< The atoms per range, summed over the steps */
+    /**< The atoms per range, summed over the steps */
+    double sum_nat[static_cast<int>(DDAtomRanges::Type::Number)] = { };
 
     /* Statistics for calls and times */
-    int    ndecomp;                    /**< The number of partioning calls */
-    int    nload;                      /**< The number of load recordings */
-    double load_step;                  /**< Total MD step time */
-    double load_sum;                   /**< Total PP force time */
-    double load_max;                   /**< Max \p load_sum over the ranks */
-    ivec   load_lim;                   /**< Was load balancing limited, per DD dim */
-    double load_mdf;                   /**< Total time on PP done during PME overlap time */
-    double load_pme;                   /**< Total time on our PME-only rank */
+    /**< The number of partioning calls */
+    int    ndecomp = 0;
+    /**< The number of load recordings */
+    int    nload = 0;
+    /**< Total MD step time */
+    double load_step = 0.0;
+    /**< Total PP force time */
+    double load_sum = 0.0;
+    /**< Max \p load_sum over the ranks */
+    double load_max = 0.0;
+    /**< Was load balancing limited, per DD dim */
+    ivec   load_lim = { };
+    /**< Total time on PP done during PME overlap time */
+    double load_mdf = 0.0;
+    /**< Total time on our PME-only rank */
+    double load_pme = 0.0;
 
     /** The last partition step */
-    int64_t partition_step;
+    int64_t partition_step = 0;
 
     /* Debugging */
-    int  nstDDDump;                    /**< Step interval for dumping the local+non-local atoms to pdb */
-    int  nstDDDumpGrid;                /**< Step interval for duming the DD grid to pdb */
-    int  DD_debug;                     /**< DD debug print level: 0, 1, 2 */
+    /**< Step interval for dumping the local+non-local atoms to pdb */
+    int  nstDDDump = 0;
+    /**< Step interval for duming the DD grid to pdb */
+    int  nstDDDumpGrid = 0;
+    /**< DD debug print level: 0, 1, 2 */
+    int  DD_debug = 0;
 };
 
 /*! \brief DD zone permutation
index 766089b5b328c61def2be54cc69aa2077c44ad97..4121a84eb5abf767a2bf185d85d54672cbdee2c3 100644 (file)
@@ -76,37 +76,48 @@ class LocalAtomSetManager;
 }
 
 typedef struct {
-    int  j0;     /* j-zone start               */
-    int  j1;     /* j-zone end                 */
-    int  cg1;    /* i-charge-group end         */
-    int  jcg0;   /* j-charge-group start       */
-    int  jcg1;   /* j-charge-group end         */
-    ivec shift0; /* Minimum shifts to consider */
-    ivec shift1; /* Maximum shifts to consider */
+    /* j-zone start               */
+    int  j0 = 0;
+    /* j-zone end                 */
+    int  j1 = 0;
+    /* i-charge-group end         */
+    int  cg1 = 0;
+    /* j-charge-group start       */
+    int  jcg0 = 0;
+    /* j-charge-group end         */
+    int  jcg1 = 0;
+    /* Minimum shifts to consider */
+    ivec shift0 = { };
+    /* Maximum shifts to consider */
+    ivec shift1 = { };
 } gmx_domdec_ns_ranges_t;
 
 typedef struct {
-    rvec x0;     /* Zone lower corner in triclinic coordinates         */
-    rvec x1;     /* Zone upper corner in triclinic coordinates         */
-    rvec bb_x0;  /* Zone bounding box lower corner in Cartesian coords */
-    rvec bb_x1;  /* Zone bounding box upper corner in Cartesian coords */
+    /* Zone lower corner in triclinic coordinates         */
+    rvec x0 = { };
+    /* Zone upper corner in triclinic coordinates         */
+    rvec x1 = { };
+    /* Zone bounding box lower corner in Cartesian coords */
+    rvec bb_x0 = { };
+    /* Zone bounding box upper corner in Cartesian coords */
+    rvec bb_x1 = { };
 } gmx_domdec_zone_size_t;
 
 struct gmx_domdec_zones_t {
     /* The number of zones including the home zone */
-    int                    n;
+    int                    n = 0;
     /* The shift of the zones with respect to the home zone */
-    ivec                   shift[DD_MAXZONE];
+    ivec                   shift[DD_MAXZONE] = { };
     /* The charge group boundaries for the zones */
-    int                    cg_range[DD_MAXZONE+1];
+    int                    cg_range[DD_MAXZONE+1] = { };
     /* The number of neighbor search zones with i-particles */
-    int                    nizone;
+    int                    nizone = 0;
     /* The neighbor search charge group ranges for each i-zone */
     gmx_domdec_ns_ranges_t izone[DD_MAXIZONE];
     /* Boundaries of the zones */
     gmx_domdec_zone_size_t size[DD_MAXZONE];
     /* The cg density of the home zone */
-    real                   dens_zone0;
+    real                   dens_zone0 = 0;
 };
 
 struct gmx_ddbox_t {
index 1a9cc81276e72070ae6dfd012b621969a8811db8..39197ecaf3d8078b3e830d850f44b125c687a91c 100644 (file)
@@ -176,6 +176,9 @@ bool pme_gpu_supports_hardware(const gmx_hw_info_t &hwinfo,
         {
             errorReasons.emplace_back("non-AMD devices");
         }
+#ifdef __APPLE__
+        errorReasons.emplace_back("Apple OS X operating system");
+#endif
     }
     return addMessageIfNotSupported(errorReasons, error);
 }
index 8b19f474a5f9f9d5dbf82dfee8edde8d1647c472..ba3ed180fe0951e80ff0fe065cb4a029f5591d00 100644 (file)
@@ -122,11 +122,15 @@ inline void reduce_atom_forces(__local float * __restrict__  sm_forces,
         int elementIndex = smemReserved + lineIndex;
         // Store input force contributions
         sm_forceReduction[elementIndex] = (dimIndex == XX) ? fx : (dimIndex == YY) ? fy : fz;
-        /* This barrier was not needed in CUDA. Different OpenCL compilers might have different ideas
+
+#if !defined(_AMD_SOURCE_)
+        /* This barrier was not needed in CUDA, nor is it needed on AMD GPUs.
+         * Different OpenCL compilers might have different ideas
          * about #pragma unroll, though. OpenCL 2 has _attribute__((opencl_unroll_hint)).
          * #2519
          */
         barrier(CLK_LOCAL_MEM_FENCE);
+#endif
 
         // Reduce to fit into smemPerDim (warp size)
 #pragma unroll
index 4a6d33c3c040f75ffe1484d06468fad4736d5508..22760e89ffe93dd90a46a24cf02e0976911769d6 100644 (file)
@@ -580,7 +580,7 @@ bool_t xdr_double(XDR * xdrs, double * dp)
 
         case XDR_ENCODE:
             ip     = reinterpret_cast<int *>(dp);
-            tmp[0] = ip[!bool(LSW)];
+            tmp[0] = ip[bool(LSW == 0)];
             tmp[1] = ip[LSW];
             return static_cast<bool_t>(bool(xdr_putint32(xdrs, tmp)) &&
                                        bool(xdr_putint32(xdrs, tmp+1)));
index 14efea694f678d1fe27bc1bea435b48b0340d3fc..aff5e005cd74caf22ee5e3d30ca28e8e52ebee9d 100644 (file)
@@ -287,7 +287,7 @@ gmx_fprintf_pqr_atomline(FILE *            fp,
     res_seq_number  = res_seq_number % 10000;
 
     int n = fprintf(fp,
-                    "%s %d %s %s %c %d %8.3f %8.3f %8.3f %6.2f %6.2f\n",
+                    "%-6s%5d %-4.4s%4.4s%c%4d %8.3f %8.3f %8.3f %6.2f %6.2f\n",
                     pdbtp[record],
                     atom_seq_number,
                     atom_name,
index a6cdf15f6b00f0402cb1d4f8984789aba7adeb35..0609301b504d41b806dd816e5ed74f14f892eb75 100644 (file)
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -1166,10 +1166,11 @@ int gmx_anaeig(int argc, char *argv[])
         {
             gmx_fatal(FARGS, "Need a second eigenvector file to do this analysis.");
         }
-        read_eigenvectors(Vec2File, &neig2, &bFit2,
+        int natoms2;
+        read_eigenvectors(Vec2File, &natoms2, &bFit2,
                           &xref2, &bDMR2, &xav2, &bDMA2, &nvec2, &eignr2, &eigvec2, &eigval2);
 
-        neig2 = DIM*neig2;
+        neig2 = std::min(nvec2, DIM*natoms2);
         if (neig2 != neig1)
         {
             gmx_fatal(FARGS, "Dimensions in the eigenvector files don't match");
index d50e22915cc44d73c1e5fb325e8a5cb7364bd59b..9fe99f540fc317ce724bda6e52d7b81fc79316d4 100644 (file)
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -750,7 +750,7 @@ static void gromos(int n1, real **mat, real rmsdcut, t_clusters *clust)
 }
 
 static rvec **read_whole_trj(const char *fn, int isize, const int index[], int skip,
-                             int *nframe, real **time,  matrix  **boxes, int **frameindexes, const gmx_output_env_t *oenv, gmx_bool bPBC, gmx_rmpbc_t gpbc)
+                             int *nframe, real **time,  matrix  **boxes, int **frameindices, const gmx_output_env_t *oenv, gmx_bool bPBC, gmx_rmpbc_t gpbc)
 {
     rvec       **xx, *x;
     matrix       box;
@@ -778,7 +778,7 @@ static rvec **read_whole_trj(const char *fn, int isize, const int index[], int s
             srenew(xx, max_nf);
             srenew(*time, max_nf);
             srenew(*boxes, max_nf);
-            srenew(*frameindexes, max_nf);
+            srenew(*frameindices, max_nf);
         }
         if ((i % skip) == 0)
         {
@@ -790,7 +790,7 @@ static rvec **read_whole_trj(const char *fn, int isize, const int index[], int s
             }
             (*time)[clusterIndex] = t;
             copy_mat(box, (*boxes)[clusterIndex]);
-            (*frameindexes)[clusterIndex] = nframes_read(status);
+            (*frameindices)[clusterIndex] = nframes_read(status);
             clusterIndex++;
         }
         i++;
@@ -976,7 +976,7 @@ static void ana_trans(t_clusters *clust, int nf,
 static void analyze_clusters(int nf, t_clusters *clust, real **rmsd,
                              int natom, t_atoms *atoms, rvec *xtps,
                              real *mass, rvec **xx, real *time,
-                             matrix *boxes, int *frameindexes,
+                             matrix *boxes, int *frameindices,
                              int ifsize, int *fitidx,
                              int iosize, int *outidx,
                              const char *trxfn, const char *sizefn,
@@ -1084,7 +1084,7 @@ static void analyze_clusters(int nf, t_clusters *clust, real **rmsd,
             fprintf(size_fp, "@g%d type %s\n", 0, "bar");
         }
     }
-    if (clustndxfn && frameindexes)
+    if (clustndxfn && frameindices)
     {
         ndxfn = gmx_ffopen(clustndxfn, "w");
     }
@@ -1209,7 +1209,7 @@ static void analyze_clusters(int nf, t_clusters *clust, real **rmsd,
             fprintf(log, "%s %6g", buf, time[i1]);
             if (ndxfn)
             {
-                fprintf(ndxfn, " %6d", frameindexes[i1]);
+                fprintf(ndxfn, " %6d", frameindices[i1] + 1);
             }
         }
         fprintf(log, "\n");
@@ -1432,7 +1432,7 @@ int gmx_cluster(int argc, char *argv[])
     real              *eigenvectors;
 
     int                isize = 0, ifsize = 0, iosize = 0;
-    int               *index = nullptr, *fitidx = nullptr, *outidx = nullptr, *frameindexes = nullptr;
+    int               *index = nullptr, *fitidx = nullptr, *outidx = nullptr, *frameindices = nullptr;
     char              *grpname;
     real               rmsd, **d1, **d2, *time = nullptr, time_invfac, *mass = nullptr;
     char               buf[STRLEN], buf1[80];
@@ -1696,7 +1696,7 @@ int gmx_cluster(int argc, char *argv[])
         /* Loop over first coordinate file */
         fn = opt2fn("-f", NFILE, fnm);
 
-        xx = read_whole_trj(fn, isize, index, skip, &nf, &time, &boxes, &frameindexes, oenv, bPBC, gpbc);
+        xx = read_whole_trj(fn, isize, index, skip, &nf, &time, &boxes, &frameindices, oenv, bPBC, gpbc);
         output_env_conv_times(oenv, nf, time);
         if (!bRMSdist || bAnalyze)
         {
@@ -1924,7 +1924,7 @@ int gmx_cluster(int argc, char *argv[])
             copy_rvec(xtps[index[i]], usextps[i]);
         }
         useatoms.nr = isize;
-        analyze_clusters(nf, &clust, rms->mat, isize, &useatoms, usextps, mass, xx, time, boxes, frameindexes,
+        analyze_clusters(nf, &clust, rms->mat, isize, &useatoms, usextps, mass, xx, time, boxes, frameindices,
                          ifsize, fitidx, iosize, outidx,
                          bReadTraj ? trx_out_fn : nullptr,
                          opt2fn_null("-sz", NFILE, fnm),
@@ -1935,7 +1935,7 @@ int gmx_cluster(int argc, char *argv[])
                          bAverage, write_ncl, write_nst, rmsmin, bFit, log,
                          rlo_bot, rhi_bot, oenv);
         sfree(boxes);
-        sfree(frameindexes);
+        sfree(frameindices);
     }
     gmx_ffclose(log);
 
index 066921f858e585f9e792f938ce402acfc1133de0..9fa96015cd0d27c8824df92938669262c9254f6b 100644 (file)
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2008, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -884,9 +884,15 @@ static t_gridcell ***init_grid(gmx_bool bBox, rvec box[], real rcut, ivec ngrid)
     }
     else
     {
-        printf("\nWill do grid-seach on %dx%dx%d grid, rcut=%g\n",
+        printf("\nWill do grid-search on %dx%dx%d grid, rcut=%3.8f\n",
                ngrid[XX], ngrid[YY], ngrid[ZZ], rcut);
     }
+    if (((ngrid[XX]*ngrid[YY]*ngrid[ZZ]) * sizeof(grid)) > INT_MAX)
+    {
+        gmx_fatal(FARGS, "Failed to allocate memory for %d x %d x %d grid points, which is larger than the maximum of %zu. "
+                  "You are likely either using a box that is too large (box dimensions are %3.8f nm x%3.8f nm x%3.8f nm) or a cutoff (%3.8f nm) that is too small.",
+                  ngrid[XX], ngrid[YY], ngrid[ZZ], INT_MAX/sizeof(grid), box[XX][XX], box[YY][YY], box[ZZ][ZZ], rcut);
+    }
     snew(grid, ngrid[ZZ]);
     for (z = 0; z < ngrid[ZZ]; z++)
     {
index 12d9d1d2276e4e6bb68c64725d744c161e52ed07..c8175ac1b880f57dd79b21d1d33b2ecd54c7d87b 100644 (file)
@@ -2050,7 +2050,8 @@ static void read_tpr_header(const char *fn, t_UmbrellaHeader* header, t_Umbrella
         header->pcrd[i].pull_type     = ir->pull->coord[i].eType;
         header->pcrd[i].geometry      = ir->pull->coord[i].eGeom;
         header->pcrd[i].ngroup        = ir->pull->coord[i].ngroup;
-        header->pcrd[i].k             = ir->pull->coord[i].k;
+        /* Angle type coordinates are handled fully in degrees in gmx wham */
+        header->pcrd[i].k             = ir->pull->coord[i].k*pull_conversion_factor_internal2userinput(&ir->pull->coord[i]);
         header->pcrd[i].init_dist     = ir->pull->coord[i].init;
 
         copy_ivec(ir->pull->coord[i].dim, header->pcrd[i].dim);
index caba460cd06edff3f482b0cfc26c487b0b734548..8a881e1640d1266f848cfd9d34a76c02d6452b1d 100644 (file)
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -124,9 +124,9 @@ static void get_params(const char *mpin, const char *mpout, t_psrec *psr)
 
     wi = init_warning(FALSE, 0);
 
-    std::string libmpin = gmx::findLibraryFile(mpin);
-    if (!libmpin.empty())
+    if (mpin != nullptr)
     {
+        std::string        libmpin = gmx::findLibraryFile(mpin);
         gmx::TextInputFile stream(libmpin);
         inp = read_inpfile(&stream, libmpin.c_str(), wi);
     }
@@ -134,6 +134,7 @@ static void get_params(const char *mpin, const char *mpout, t_psrec *psr)
     {
         inp.clear();
     }
+
     psr->bw        = get_eenum(&inp, "black&white",             gmx_bools);
     psr->linewidth = get_ereal(&inp, "linewidth",      1.0, wi);
     setStringEntry(&inp, "titlefont",      psr->titfont,        "Helvetica");
index 089d178113058187f373a42bee990368fa68f604..bbec8d7fc9a21c317e6ecb524a9cbeef2661aa69 100644 (file)
@@ -781,15 +781,15 @@ int gmx_editconf(int argc, char *argv[])
         printf("Incompatible options -mead and -grasp. Turning off -grasp\n");
         bGrasp = FALSE;
     }
-    if ((bGrasp || bCONECT) && (outftp != efPDB))
+    if (bGrasp && (outftp != efPDB))
     {
         gmx_fatal(FARGS, "Output file should be a .pdb file"
-                  " when using the -grasp or -connect options\n");
+                  " when using the -grasp option\n");
     }
-    if ((bMead || bGrasp || bCONECT) && (fn2ftp(infile) != efTPR))
+    if ((bMead || bGrasp) && (fn2ftp(infile) != efTPR))
     {
         gmx_fatal(FARGS, "Input file should be a .tpr file"
-                  " when using the -mead or -connect options\n");
+                  " when using the -mead option\n");
     }
 
     t_symtab  symtab;
index a8e621ca342a16c925ce0a0cf2bc4f39e872b68f..d232b55dbf109e4f0aae681c6b3cc9999a114990 100644 (file)
@@ -453,6 +453,11 @@ process_directive(gmx_cpp_t         *handlep,
         {
             return eCPP_SYNTAX;
         }
+        // An include needs to be followed by either a '"' or a '<' as a first character.
+        if ((dval[0] != '"') && (dval[0] != '<'))
+        {
+            return eCPP_INVALID_INCLUDE_DELIMITER;
+        }
         for (size_t i1 = 0; i1 < dval.size(); i1++)
         {
             if ((dval[i1] == '"') || (dval[i1] == '<') || (dval[i1] == '>'))
@@ -575,7 +580,7 @@ int cpp_read_line(gmx_cpp_t *handlep, int n, char buf[])
             if (!bEOF)
             {
                 /* Something strange happened, fgets returned NULL,
-                 * but we are not at EOF.
+                 * but we are not at EOF. Maybe wrong line endings?
                  */
                 return eCPP_UNKNOWN;
             }
@@ -734,8 +739,8 @@ char *cpp_error(gmx_cpp_t *handlep, int status)
     char        buf[256];
     const char *ecpp[] = {
         "OK", "File not found", "End of file", "Syntax error", "Interrupted",
-        "Invalid file handle",
-        "File not open", "Unknown error", "Error status out of range"
+        "Invalid file handle", "Invalid delimiter for filename in #include statement",
+        "File not open", "Unknown error, perhaps your text file uses wrong line endings?", "Error status out of range"
     };
     gmx_cpp_t   handle = *handlep;
 
index f85fe06b54287f0df10410a03ee3c3ea1869fa38..988dd1a23a59b9b7497808d018707683fff4ee18 100644 (file)
@@ -45,7 +45,7 @@ typedef struct gmx_cpp *gmx_cpp_t;
 /* The possible return codes for these functions */
 enum {
     eCPP_OK, eCPP_FILE_NOT_FOUND, eCPP_EOF, eCPP_SYNTAX, eCPP_INTERRUPT,
-    eCPP_INVALID_HANDLE,
+    eCPP_INVALID_HANDLE, eCPP_INVALID_INCLUDE_DELIMITER,
     eCPP_FILE_NOT_OPEN, eCPP_UNKNOWN, eCPP_NR
 };
 
index 553cd9f232ac7024db4f3d1cccb82f87424c6e47..43a317cc211793837f3f4531f07e310b3c5bc733 100644 (file)
@@ -2356,8 +2356,14 @@ int gmx_grompp(int argc, char *argv[])
 
     if (ir->bDoAwh)
     {
+        tensor compressibility = { { 0 } };
+        if (ir->epc != epcNO)
+        {
+            copy_mat(ir->compress, compressibility);
+        }
         setStateDependentAwhParams(ir->awhParams, ir->pull, pull,
-                                   state.box, ir->ePBC, &ir->opts, wi);
+                                   state.box, ir->ePBC, compressibility,
+                                   &ir->opts, wi);
     }
 
     if (ir->bPull)
index d49e1dfa8fb877e906a4236794456b894fbb48b6..bd598b7b46d0d42c1bb4565c2a384df6adfffa59 100644 (file)
@@ -547,7 +547,14 @@ void push_at (t_symtab *symtab, PreprocessingAtomTypes *at, PreprocessingBondAto
 
     if ((nr = at->atomTypeFromName(type)) != NOTSET)
     {
-        auto message = gmx::formatString("Overriding atomtype %s", type);
+        auto message = gmx::formatString
+                ("Atomtype %s was defined previously (e.g. in the forcefield files), "
+                "and has now been defined again. This could happen e.g. if you would "
+                "use a self-contained molecule .itp file that duplicates or replaces "
+                "the contents of the standard force-field files. You should check "
+                "the contents of your files and remove such repetition. If you know "
+                "you should override the previous definition, then you could choose "
+                "to suppress this warning with -maxwarn.", type);
         warning(wi, message);
         if ((nr = at->setType(nr, symtab, *atom, type, interactionType, batype_nr,
                               atomnr)) == NOTSET)
@@ -655,7 +662,13 @@ static void push_bondtype(InteractionsOfType              *bt,
                 else if (!haveWarned)
                 {
                     auto message = gmx::formatString
-                            ("Overriding %s parameters.%s",
+                            ("Bondtype %s was defined previously (e.g. in the forcefield files), "
+                            "and has now been defined again. This could happen e.g. if you would "
+                            "use a self-contained molecule .itp file that duplicates or replaces "
+                            "the contents of the standard force-field files. You should check "
+                            "the contents of your files and remove such repetition. If you know "
+                            "you should override the previous definition, then you could choose "
+                            "to suppress this warning with -maxwarn.%s",
                             interaction_function[ftype].longname,
                             (ftype == F_PDIHS) ?
                             "\nUse dihedraltype 9 to allow several multiplicity terms. Only consecutive "
@@ -1128,7 +1141,14 @@ void push_nbt(Directive d, t_nbparam **nbt, PreprocessingAtomTypes *atypes,
         }
         if (!bId)
         {
-            auto message = gmx::formatString("Overriding non-bonded parameters,");
+            auto message = gmx::formatString
+                    ("Non-bonded parameters were defined previously (e.g. in the forcefield files), "
+                    "and have now been defined again. This could happen e.g. if you would "
+                    "use a self-contained molecule .itp file that duplicates or replaces "
+                    "the contents of the standard force-field files. You should check "
+                    "the contents of your files and remove such repetition. If you know "
+                    "you should override the previous definitions, then you could choose "
+                    "to suppress this warning with -maxwarn.");
             warning(wi, message);
             fprintf(stderr, "  old:");
             for (i = 0; i < nrfp; i++)
index 8f8afc9936443fe5a698ec746fbc8ad10a465e3c..61afde793559aa24612ca948c58980197879e310 100644 (file)
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2015,2016,2018, by the GROMACS development team, led by
+ * Copyright (c) 2015,2016,2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -49,6 +49,8 @@
 
 #include <gtest/gtest.h>
 
+#include "gromacs/utility/stringutil.h"
+
 namespace
 {
 
@@ -97,31 +99,39 @@ TEST(HardwareTopologyTest, ProcessorSelfconsistency)
 
     if (hwTop.supportLevel() >= gmx::HardwareTopology::SupportLevel::Basic)
     {
-        int socketsInMachine = hwTop.machine().sockets.size();
-        int coresPerSocket   = hwTop.machine().sockets[0].cores.size();
-        int hwThreadsPerCore = hwTop.machine().sockets[0].cores[0].hwThreads.size();
+        SCOPED_TRACE(gmx::formatString("Logical Processor count %d", hwTop.machine().logicalProcessorCount));
 
-        // Check that logical processor information is reasonable
-        for (auto &l : hwTop.machine().logicalProcessors)
-        {
-            EXPECT_TRUE(l.socketRankInMachine >= 0 && l.socketRankInMachine < socketsInMachine);
-            EXPECT_TRUE(l.coreRankInSocket >= 0 && l.coreRankInSocket < coresPerSocket);
-            EXPECT_TRUE(l.hwThreadRankInCore >= 0 && l.hwThreadRankInCore < hwThreadsPerCore);
-        }
+        int  socketsInMachine = hwTop.machine().sockets.size();
+        int  coresPerSocket   = hwTop.machine().sockets[0].cores.size();
+        int  hwThreadsPerCore = hwTop.machine().sockets[0].cores[0].hwThreads.size();
 
-        // Double-check that the tree is self-consistent with logical processor info
-        for (int s = 0; s < socketsInMachine; s++)
+        auto logicalProcessors = hwTop.machine().logicalProcessors;
+        for (auto logicalProcessorIt = logicalProcessors.begin();
+             logicalProcessorIt != logicalProcessors.end();
+             ++logicalProcessorIt)
         {
-            for (int c = 0; c < coresPerSocket; c++)
+            // Check that logical processor information contains
+            // reasonable values.
+            SCOPED_TRACE(gmx::formatString("Socket rank in machine: %d", logicalProcessorIt->socketRankInMachine));
+            SCOPED_TRACE(gmx::formatString("Core rank in socket:    %d", logicalProcessorIt->coreRankInSocket));
+            SCOPED_TRACE(gmx::formatString("Hw thread rank in core: %d", logicalProcessorIt->hwThreadRankInCore));
+            EXPECT_TRUE(logicalProcessorIt->socketRankInMachine >= 0 && logicalProcessorIt->socketRankInMachine < socketsInMachine);
+            EXPECT_TRUE(logicalProcessorIt->coreRankInSocket >= 0 && logicalProcessorIt->coreRankInSocket < coresPerSocket);
+            EXPECT_TRUE(logicalProcessorIt->hwThreadRankInCore >= 0 && logicalProcessorIt->hwThreadRankInCore < hwThreadsPerCore);
+            // Check that logical processor information is distinct
+            // for each logical processor.
+
+            for (auto remainingLogicalProcessorIt = logicalProcessorIt + 1;
+                 remainingLogicalProcessorIt != logicalProcessors.end();
+                 ++remainingLogicalProcessorIt)
             {
-                for (int t = 0; t < hwThreadsPerCore; t++)
-                {
-                    int idx = hwTop.machine().sockets[s].cores[c].hwThreads[t].logicalProcessorId;
-                    EXPECT_LT(idx, hwTop.machine().logicalProcessorCount);
-                    EXPECT_EQ(hwTop.machine().logicalProcessors[idx].socketRankInMachine, s);
-                    EXPECT_EQ(hwTop.machine().logicalProcessors[idx].coreRankInSocket, c) << "logical:" << idx;
-                    EXPECT_EQ(hwTop.machine().logicalProcessors[idx].hwThreadRankInCore, t);
-                }
+                SCOPED_TRACE(gmx::formatString("Other socket rank in machine: %d", remainingLogicalProcessorIt->socketRankInMachine));
+                SCOPED_TRACE(gmx::formatString("Other core rank in socket:    %d", remainingLogicalProcessorIt->coreRankInSocket));
+                SCOPED_TRACE(gmx::formatString("Other hw thread rank in core: %d", remainingLogicalProcessorIt->hwThreadRankInCore));
+                EXPECT_TRUE((logicalProcessorIt->socketRankInMachine != remainingLogicalProcessorIt->socketRankInMachine) ||
+                            (logicalProcessorIt->coreRankInSocket    != remainingLogicalProcessorIt->coreRankInSocket) ||
+                            (logicalProcessorIt->hwThreadRankInCore  != remainingLogicalProcessorIt->hwThreadRankInCore)) <<
+                "This pair of logical processors have the same descriptive information, which is an error";
             }
         }
     }
index 3ad73d5dc1aa76a8be397c46ae93d4d58f287a24..d10dd20f60c52cb3d27ff479539c9bb058fc14df 100644 (file)
@@ -340,7 +340,7 @@ calc_bonded_reduction_mask(int                       natoms,
     if (bondedThreading.nthreads > BITMASK_SIZE)
     {
 #pragma omp master
-        gmx_fatal(FARGS, "You are using %d OpenMP threads, which is larger than GMX_OPENMP_MAX_THREADS (%d). Decrease the number of OpenMP threads or rebuild GROMACS with a larger value for GMX_OPENMP_MAX_THREADS.",
+        gmx_fatal(FARGS, "You are using %d OpenMP threads, which is larger than GMX_OPENMP_MAX_THREADS (%d). Decrease the number of OpenMP threads or rebuild GROMACS with a larger value for GMX_OPENMP_MAX_THREADS passed to CMake.",
                   bondedThreading.nthreads, GMX_OPENMP_MAX_THREADS);
 #pragma omp barrier
     }
index eebf97f2f475b3b25f4119fea5597e222b73f7b4..6eb325684e794c5ad1aaa75e125f169bafc5a5b6 100644 (file)
@@ -836,7 +836,10 @@ void bcast_ir_mtop(const t_commrec *cr, t_inputrec *inputrec, gmx_mtop_t *mtop)
     block_bc(cr, mtop->bIntermolecularInteractions);
     if (mtop->bIntermolecularInteractions)
     {
-        mtop->intermolecular_ilist = std::make_unique<InteractionLists>();
+        if (!MASTER(cr))
+        {
+            mtop->intermolecular_ilist = std::make_unique<InteractionLists>();
+        }
         bc_ilists(cr, mtop->intermolecular_ilist.get());
     }
 
index c3ccdafec3ec2dcc715884efdad874f00dc85b09..94747562493edcf54cb9b6cb096241ebc254f2b3 100644 (file)
@@ -530,6 +530,7 @@ static int gen_rm_list(rm_t *rm_p, t_block *ins_at, t_block *rest_at, t_pbc *pbc
     r_min_rad = probe_rad*probe_rad;
     gmx::RangePartitioning molecules = gmx_mtop_molecules(*mtop);
     snew(rm_p->block, molecules.numBlocks());
+    snew(rm_p->mol, molecules.numBlocks());
     nrm    = nupper = 0;
     nlower = low_up_rm;
     for (i = 0; i < ins_at->nr; i++)
@@ -689,7 +690,7 @@ static void rm_group(SimulationGroups *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_s
     for (int i = 0; i < rm_p->nr; i++)
     {
         mol_id = rm_p->mol[i];
-        at     = molecules.block(mol_id).size();
+        at     = molecules.block(mol_id).begin();
         block  = rm_p->block[i];
         mtop->molblock[block].nmol--;
         for (j = 0; j < mtop->moltype[mtop->molblock[block].type].atoms.nr; j++)
@@ -700,23 +701,27 @@ static void rm_group(SimulationGroups *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_s
     }
 
     mtop->natoms    -= n;
-    state_change_natoms(state, state->natoms - n);
-    snew(x_tmp, state->natoms);
-    snew(v_tmp, state->natoms);
+    /* We cannot change the size of the state datastructures here
+     * because we still access the coordinate arrays for all positions
+     * before removing the molecules we want to remove.
+     */
+    const int newStateAtomNumber = state->natoms - n;
+    snew(x_tmp, newStateAtomNumber);
+    snew(v_tmp, newStateAtomNumber);
 
     for (auto group : keysOf(groups->groupNumbers))
     {
         if (!groups->groupNumbers[group].empty())
         {
-            groups->groupNumbers[group].resize(state->natoms);
-            new_egrp[group].resize(state->natoms);
+            groups->groupNumbers[group].resize(newStateAtomNumber);
+            new_egrp[group].resize(newStateAtomNumber);
         }
     }
 
     auto x = makeArrayRef(state->x);
     auto v = makeArrayRef(state->v);
     rm = 0;
-    for (int i = 0; i < state->natoms+n; i++)
+    for (int i = 0; i < state->natoms; i++)
     {
         bRM = FALSE;
         for (j = 0; j < n; j++)
@@ -759,6 +764,7 @@ static void rm_group(SimulationGroups *groups, gmx_mtop_t *mtop, rm_t *rm_p, t_s
             }
         }
     }
+    state_change_natoms(state, newStateAtomNumber);
     for (int i = 0; i < state->natoms; i++)
     {
         copy_rvec(x_tmp[i], x[i]);
index fa49e80294d24b1634c6206984af6840b346e40c..7d8a05465e1b5f1d7921bc60a30f122276608e1a 100644 (file)
@@ -237,7 +237,7 @@ make_shake_sblock_serial(shakedata *shaked,
         fprintf(debug, "Going to sort constraints\n");
     }
 
-    qsort(sb, ncons, sizeof(*sb), pcomp);
+    std::qsort(sb, ncons, sizeof(*sb), pcomp);
 
     if (debug)
     {
index 5ff0d7a32afc33b10df3525cf0285ba1e405754e..26b4fb93c8212d084cf09aaf0a2e44433c31f990 100644 (file)
@@ -235,7 +235,7 @@ class LegacyMdrunOptions
             { "-append",  FALSE, etBOOL, {&appendOption},
               "Append to previous output files when continuing from checkpoint instead of adding the simulation part number to all file names" },
             { "-nsteps",  FALSE, etINT64, {&mdrunOptions.numStepsCommandline},
-              "Run this number of steps, overrides .mdp file option (-1 means infinite, -2 means use mdp option, smaller is invalid)" },
+              "Run this number of steps (-1 means infinite, -2 means use mdp option, smaller is invalid)" },
             { "-maxh",   FALSE, etREAL, {&mdrunOptions.maximumHoursToRun},
               "Terminate after 0.99 times this time (hours)" },
             { "-replex",  FALSE, etINT, {&replExParams.exchangeInterval},
index 3cfd3f77c65304a901359f0358b22d3c0c72eac6..76addf8c6552c6cdcd0d8e01b39dcb23038b1119 100644 (file)
@@ -2125,15 +2125,15 @@ LegacySimulator::do_lbfgs()
                 {
                     /* Replace c endpoint with b */
                     c   = b;
-                    /* swap states b and c */
-                    swap_em_state(&sb, &sc);
+                    /* copy state b to c */
+                    *sc = *sb;
                 }
                 else
                 {
                     /* Replace a endpoint with b */
                     a   = b;
-                    /* swap states a and b */
-                    swap_em_state(&sa, &sb);
+                    /* copy state b to a */
+                    *sa = *sb;
                 }
 
                 /*
@@ -2339,7 +2339,7 @@ LegacySimulator::do_lbfgs()
         }
 
         // Reset stepsize in we are doing more iterations
-        stepsize = 1.0/ems.fnorm;
+        stepsize = 1.0;
 
         /* Stop when the maximum force lies below tolerance.
          * If we have reached machine precision, converged is already set to true.
@@ -2597,8 +2597,15 @@ LegacySimulator::do_steep()
             }
         }
 
-        /* Determine new step  */
-        stepsize = ustep/s_min->fmax;
+        // If the force is very small after finishing minimization,
+        // we risk dividing by zero when calculating the step size.
+        // So we check first if the minimization has stopped before
+        // trying to obtain a new step size.
+        if (!bDone)
+        {
+            /* Determine new step  */
+            stepsize = ustep/s_min->fmax;
+        }
 
         /* Check if stepsize is too small, with 1 nm as a characteristic length */
 #if GMX_DOUBLE
index 0c354711d517219d46c83c8484cf60453547ab45..79ebf9fd33fc50042233446d27e68de5c780b7f6 100644 (file)
@@ -919,6 +919,18 @@ int Mdrunner::mdrunner()
 
     if (startingBehavior != StartingBehavior::NewSimulation)
     {
+        /* Check if checkpoint file exists before doing continuation.
+         * This way we can use identical input options for the first and subsequent runs...
+         */
+        if (mdrunOptions.numStepsCommandline > -2)
+        {
+            /* Temporarily set the number of steps to unmlimited to avoid
+             * triggering the nsteps check in load_checkpoint().
+             * This hack will go away soon when the -nsteps option is removed.
+             */
+            inputrec->nsteps = -1;
+        }
+
         load_checkpoint(opt2fn_master("-cpi", filenames.size(), filenames.data(), cr),
                         logFileHandle,
                         cr, domdecOptions.numCells,
index ee890ba7a4367913cf86180efbde8d697f35d90c..144f194cb5d9630337142beda35b25415e007392 100644 (file)
@@ -3763,6 +3763,31 @@ static void copySelectedListRange(const nbnxn_ci_t * gmx_restrict srcCi,
     }
 }
 
+#if defined(__GNUC__) && !defined(__clang__) && !defined(__ICC) && __GNUC__ == 7
+/* Avoid gcc 7 avx512 loop vectorization bug (actually only needed with -mavx512f) */
+#pragma GCC push_options
+#pragma GCC optimize ("no-tree-vectorize")
+#endif
+
+/* Returns the number of cluster pairs that are in use summed over all lists */
+static int countClusterpairs(gmx::ArrayRef<const NbnxnPairlistCpu> pairlists)
+{
+    /* gcc 7 with -mavx512f can miss the contributions of 16 consecutive
+     * elements to the sum calculated in this loop. Above we have disabled
+     * loop vectorization to avoid this bug.
+     */
+    int ncjTotal = 0;
+    for (const auto &pairlist : pairlists)
+    {
+        ncjTotal += pairlist.ncjInUse;
+    }
+    return ncjTotal;
+}
+
+#if defined(__GNUC__) && !defined(__clang__) && !defined(__ICC) && __GNUC__ == 7
+#pragma GCC pop_options
+#endif
+
 /* This routine re-balances the pairlists such that all are nearly equally
  * sized. Only whole i-entries are moved between lists. These are moved
  * between the ends of the lists, such that the buffer reduction cost should
@@ -3775,11 +3800,7 @@ static void rebalanceSimpleLists(gmx::ArrayRef<const NbnxnPairlistCpu> srcSet,
                                  gmx::ArrayRef<NbnxnPairlistCpu>       destSet,
                                  gmx::ArrayRef<PairsearchWork>         searchWork)
 {
-    int ncjTotal = 0;
-    for (auto &src : srcSet)
-    {
-        ncjTotal += src.ncjInUse;
-    }
+    const int ncjTotal  = countClusterpairs(srcSet);
     const int numLists  = srcSet.ssize();
     const int ncjTarget = (ncjTotal + numLists - 1)/numLists;
 
@@ -3848,11 +3869,7 @@ static void rebalanceSimpleLists(gmx::ArrayRef<const NbnxnPairlistCpu> srcSet,
     }
 
 #ifndef NDEBUG
-    int ncjTotalNew = 0;
-    for (auto &dest : destSet)
-    {
-        ncjTotalNew += dest.ncjInUse;
-    }
+    const int ncjTotalNew = countClusterpairs(destSet);
     GMX_RELEASE_ASSERT(ncjTotalNew == ncjTotal, "The total size of the lists before and after rebalancing should match");
 #endif
 }
index 35a917ef8b29b31d5ba720bcb5620ba7746f2733..1b0663bb1d919289b68bdc6da2c5d70251e60f35 100644 (file)
@@ -1597,7 +1597,7 @@ void dd_make_local_pull_groups(const t_commrec *cr, struct pull_t *pull)
 
     for (pull_group_work_t &group : pull->group)
     {
-        if (group.epgrppbc == epgrppbcCOS || !group.globalWeights.empty())
+        if (!group.globalWeights.empty())
         {
             group.localWeights.resize(group.atomSet.numAtomsLocal());
             for (size_t i = 0; i < group.atomSet.numAtomsLocal(); ++i)
index af6b2ed5d0e9ea388b89f9fd2db415bc1bd8c3c1..95ed32a279753bb062838b1aa18913aa630eac79 100644 (file)
@@ -579,6 +579,15 @@ void pull_calc_coms(const t_commrec *cr,
     {
         pull_group_work_t *pgrp      = &pull->group[g];
 
+        /* Cosine-weighted COMs behave different from all other weighted COMs
+         * in the sense that the weights depend on instantaneous coordinates,
+         * not on pre-set weights. Thus we resize the local weight buffer here.
+         */
+        if (pgrp->epgrppbc == epgrppbcCOS)
+        {
+            pgrp->localWeights.resize(pgrp->atomSet.localIndex().size());
+        }
+
         auto               comBuffer =
             gmx::arrayRefFromArray(comm->comBuffer.data() + g*c_comBufferStride, c_comBufferStride);
 
index b514513fd9e479f138423c30aa58838b83818df7..23911799c27671beb6e14a88f41ba7a1ff64ff32 100644 (file)
@@ -885,21 +885,23 @@ gmx_ana_index_make_block(t_blocka *t, const gmx_mtop_t *top, gmx_ana_index_t *g,
                     {
                         int            molnr, atnr_mol;
                         mtopGetMolblockIndex(top, ai, &molb, &molnr, &atnr_mol);
-                        const t_atoms &mol_atoms = top->moltype[top->molblock[molb].type].atoms;
-                        int            last_atom = atnr_mol + 1;
+                        const t_atoms &mol_atoms    = top->moltype[top->molblock[molb].type].atoms;
+                        int            last_atom    = atnr_mol + 1;
+                        const int      currentResid = mol_atoms.atom[atnr_mol].resind;
                         while (last_atom < mol_atoms.nr
-                               && mol_atoms.atom[last_atom].resind == id)
+                               && mol_atoms.atom[last_atom].resind == currentResid)
                         {
                             ++last_atom;
                         }
                         int first_atom = atnr_mol - 1;
                         while (first_atom >= 0
-                               && mol_atoms.atom[first_atom].resind == id)
+                               && mol_atoms.atom[first_atom].resind == currentResid)
                         {
                             --first_atom;
                         }
-                        int first_mol_atom = top->moleculeBlockIndices[molb].globalAtomStart;
-                        first_mol_atom += molnr*top->moleculeBlockIndices[molb].numAtomsPerMolecule;
+                        const MoleculeBlockIndices &molBlock = top->moleculeBlockIndices[molb];
+                        int first_mol_atom                   = molBlock.globalAtomStart;
+                        first_mol_atom += molnr*molBlock.numAtomsPerMolecule;
                         first_atom      = first_mol_atom + first_atom + 1;
                         last_atom       = first_mol_atom + last_atom - 1;
                         for (int j = first_atom; j <= last_atom; ++j)
@@ -910,11 +912,8 @@ gmx_ana_index_make_block(t_blocka *t, const gmx_mtop_t *top, gmx_ana_index_t *g,
                     }
                     case INDEX_MOL:
                     {
-                        size_t molb = 0;
-                        while (molb + 1 < top->molblock.size() && id >= top->moleculeBlockIndices[molb].moleculeIndexStart)
-                        {
-                            ++molb;
-                        }
+                        int                         molnr, atnr_mol;
+                        mtopGetMolblockIndex(top, ai, &molb, &molnr, &atnr_mol);
                         const MoleculeBlockIndices &blockIndices  = top->moleculeBlockIndices[molb];
                         const int                   atomStart     = blockIndices.globalAtomStart + (id - blockIndices.moleculeIndexStart)*blockIndices.numAtomsPerMolecule;
                         for (int j = 0; j < blockIndices.numAtomsPerMolecule; ++j)
index 6e3440d9c3bcc62033b8795e71dc13b673f852dc..70a1a51b916cd40fe0c5ff51da7fba4a54732a37 100644 (file)
@@ -220,7 +220,7 @@ static const char *const help_resindex[] = {
     "[TT]resnr[tt] selects atoms using the residue numbering in the input",
     "file. [TT]resid[tt] is synonym for this keyword for VMD compatibility.",
     "",
-    "[TT]resindex N[tt] selects the [TT]N[tt]th residue starting from the",
+    "[TT]resindex N[tt] selects the [TT]N[tt] th residue starting from the",
     "beginning of the input file. This is useful for uniquely identifying",
     "residues if there are duplicate numbers in the input file (e.g., in",
     "multiple chains).",
index f0f6e4e735f541801c7f966ad62b7397c5a60dc1..87e9dd7ac7d1d0a3462ea61288e042802e64906d 100644 (file)
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2014,2015,2017,2018, by the GROMACS development team, led by
+ * Copyright (c) 2014,2015,2017,2018,2019, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -232,8 +232,13 @@ namespace gmx
 //! \brief The SIMD4 type is always four units wide, but this makes code more explicit
 #define GMX_SIMD4_WIDTH                                          4
 
-//! \brief Required alignment in bytes for aligned load/store (always defined, even without SIMD)
-#define GMX_SIMD_ALIGNMENT                                       8 // 8 (1*double)
+/*! \brief Maximum required alignment in bytes for aligned load/store of multiple
+ * values (maximum required for either float or double). */
+#if GMX_SIMD_DOUBLE_WIDTH >= 2*GMX_SIMD_FLOAT_WIDTH
+#define GMX_SIMD_ALIGNMENT                                       (GMX_SIMD_DOUBLE_WIDTH*8)
+#else
+#define GMX_SIMD_ALIGNMENT                                       (GMX_SIMD_FLOAT_WIDTH*4)
+#endif
 
 //! \brief Accuracy of SIMD 1/sqrt(x) lookup. Used to determine number of iterations.
 #define GMX_SIMD_RSQRT_BITS                                      23
index 9b9d517aa68068b3e9a627ae66544aa1f4bfd2b3..6255e0195fa68a9d668715d4fd0d769a8506860a 100644 (file)
@@ -280,8 +280,8 @@ Other options:
            Append to previous output files when continuing from checkpoint
            instead of adding the simulation part number to all file names
  -nsteps &lt;int&gt;              (-2)
-           Run this number of steps, overrides .mdp file option (-1 means
-           infinite, -2 means use mdp option, smaller is invalid)
+           Run this number of steps (-1 means infinite, -2 means use mdp
+           option, smaller is invalid)
  -maxh   &lt;real&gt;             (-1)
            Terminate after 0.99 times this time (hours)
  -replex &lt;int&gt;              (0)
index e63e2321c77622b5be6ca0c8b98315fc65653e88..9db995ac38389dee76a5fd1d7beaeee2aa0d1123 100644 (file)
@@ -6,7 +6,7 @@
       <Energy Name="Potential">
         <Real Name="Time 0.000000 Step 0 in frame 0">2195.7786482024485</Real>
         <Real Name="Time 0.000000 Step 0 in frame 1">1848.1873657020258</Real>
-        <Real Name="Time 4.000000 Step 4 in frame 2">1847.227343785434</Real>
+        <Real Name="Time 4.000000 Step 4 in frame 2">561.1160975330763</Real>
       </Energy>
     </Minimizer>
   </Simulation>
diff --git a/src/testutils/simulationdatabase/spc_and_methane.gro b/src/testutils/simulationdatabase/spc_and_methane.gro
new file mode 100644 (file)
index 0000000..ddf626f
--- /dev/null
@@ -0,0 +1,11 @@
+water and methane
+  8
+    1SOL     OW    1   0.005   0.600   0.244  0.1823 -0.4158  0.4875
+    1SOL    HW1    2  -0.017   0.690   0.270 -1.7457 -0.5883 -0.4604
+    1SOL    HW2    3   0.051   0.610   0.161  2.5085 -0.1501  1.7627
+    1CH4      C  649  -0.024  -0.222  -0.640  0.0000  0.0000  0.0000
+    1CH4     H1  650  -0.083  -0.303  -0.646  0.0000  0.0000  0.0000
+    1CH4     H2  651  -0.080  -0.140  -0.642  0.0000  0.0000  0.0000
+    1CH4     H3  652   0.040  -0.221  -0.716  0.0000  0.0000  0.0000
+    1CH4     H4  653   0.027  -0.225  -0.553  0.0000  0.0000  0.0000
+   3.72412   3.72412   3.72412
\ No newline at end of file
diff --git a/src/testutils/simulationdatabase/spc_and_methane.ndx b/src/testutils/simulationdatabase/spc_and_methane.ndx
new file mode 100644 (file)
index 0000000..8c7e887
--- /dev/null
@@ -0,0 +1,6 @@
+[ System ]
+   1    2    3    4    5    6    7    8
+[ SOL ]
+   1    2    3
+[ CH4 ]
+   4    5    6    7    8
diff --git a/src/testutils/simulationdatabase/spc_and_methane.top b/src/testutils/simulationdatabase/spc_and_methane.top
new file mode 100644 (file)
index 0000000..d3ab903
--- /dev/null
@@ -0,0 +1,33 @@
+#include "oplsaa.ff/forcefield.itp"
+
+[ atomtypes ]
+;name  at.num      mass        charge   ptype       sigma        epsilon
+  CH4    0          0.0        0.000       A        0.371        1.26
+
+#include "oplsaa.ff/spce.itp"
+
+[ moleculetype ]
+; Name            nrexcl
+ methane          3
+
+[ atoms ]
+;   nr       type  resnr residue  atom   cgnr     charge       mass
+     1   opls_138      1    METH     C      1      -0.24     12.011
+     2   opls_140      1    METH    H1      1       0.06      1.008
+     3   opls_140      1    METH    H2      1       0.06      1.008
+     4   opls_140      1    METH    H3      1       0.06      1.008
+     5   opls_140      1    METH    H4      1       0.06      1.008
+
+[ bonds ]
+;  ai    aj funct            c0            c1            c2            c3
+    1     2     1
+    1     3     1
+    1     4     1
+    1     5     1
+
+[ system ]
+Water and methane
+
+[ molecules ]
+SOL       1
+methane   1