Finally move OpenMM to contrib
authorMark Abraham <mark.j.abraham@gmail.com>
Mon, 14 Jan 2013 00:47:02 +0000 (01:47 +0100)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Thu, 17 Jan 2013 00:33:02 +0000 (01:33 +0100)
I've paid lip service to attempting to leave the build functional.
There's a few places GMX_OPENMM still has to appear in the main code,
because of how CMake works and how badly GROMACS libraries work. I'm
prepared to fix obvious major issues anyone spots. I'm proposing not
to do any testing of whether the build works, since the 4.5 build is
still available. If later on someone wants both OpenMM and some
feature in 4.6 they can ask then for help fixing the damage here.

Fixes #1079, for certain values of "fixes" :-)

Change-Id: Ic27b8c3bdc1255cc3edb01bcfa11f47d6a03b44d

26 files changed:
CMakeLists.txt
COPYING
INSTALL-GPU [deleted file]
README-GPU [deleted file]
admin/installguide/installguide.tex
include/gpu_utils.h
include/mdrun.h
src/contrib/BuildMdrunOpenMM.cmake [new file with mode: 0644]
src/contrib/CMakeLists.txt
src/contrib/FindOpenMM.cmake [moved from cmake/FindOpenMM.cmake with 100% similarity]
src/contrib/md_openmm.c [moved from src/kernel/md_openmm.c with 100% similarity]
src/contrib/md_openmm.h [moved from src/kernel/md_openmm.h with 100% similarity]
src/contrib/mdrun_openmm.c [new file with mode: 0644]
src/contrib/openmm_gpu_utils.cu [new file with mode: 0644]
src/contrib/openmm_gpu_utils.h [new file with mode: 0644]
src/contrib/openmm_wrapper.cpp [moved from src/kernel/openmm_wrapper.cpp with 100% similarity]
src/contrib/openmm_wrapper.h [moved from src/kernel/openmm_wrapper.h with 76% similarity]
src/contrib/runner_openmm.c [new file with mode: 0644]
src/gmxlib/gpu_utils/gpu_utils.cu
src/kernel/CMakeLists.txt
src/kernel/mdrun.c
src/kernel/membed.c
src/kernel/runner.c
src/mdlib/mdebin.c
src/tools/CMakeLists.txt
src/tools/gmx_membed.c

index 68c3306d4eae4c33802133ce9e2d281f004468ab..0a39a467718a2abd4efab81605a40d873fc8bfcf 100644 (file)
@@ -142,17 +142,13 @@ endif()
 # decide on GPU settings based on user-settings and GPU/CUDA detection
 include(gmxManageGPU)
 
-# TODO: move OpenMM to contrib
-option(GMX_OPENMM "Accelerated execution on GPUs through the OpenMM library (not fully supported)" OFF)
-mark_as_advanced(GMX_OPENMM)
-
 option(GMX_FORCE_CXX "Enable C++ compilation even if not necessary" OFF)
 mark_as_advanced(GMX_FORCE_CXX)
 
 option(GMX_COOL_QUOTES "Enable Gromacs cool quotes" ON)
 mark_as_advanced(GMX_COOL_QUOTES)
 
-if(GMX_GPU OR GMX_OPENMM OR GMX_FORCE_CXX)
+if(GMX_GPU OR GMX_FORCE_CXX)
     enable_language(CXX)
 endif()
 set(CMAKE_PREFIX_PATH "" CACHE STRING "Extra locations to search for external libraries and tools (give directory without lib, bin, or include)")
@@ -317,10 +313,6 @@ if (GMX_DEFAULT_SUFFIX)
     set (GMX_BINARY_SUFFIX "${GMX_BINARY_SUFFIX}_d")
     set (GMX_LIBS_SUFFIX "${GMX_LIBS_SUFFIX}_d")
   endif(GMX_DOUBLE)
-  if (GMX_OPENMM)
-    set (GMX_BINARY_SUFFIX "-openmm")
-    set (GMX_LIBS_SUFFIX "_openmm")
-  endif(GMX_OPENMM)
   mark_as_advanced(FORCE GMX_BINARY_SUFFIX GMX_LIBS_SUFFIX)
   if (NOT SUFFIX_QUIETLY)
     message(STATUS "Using default binary suffix: \"${GMX_BINARY_SUFFIX}\"")
@@ -343,53 +335,6 @@ if(GMX_SOFTWARE_INVSQRT)
   set(PKG_CFLAGS "${PKG_CFLAGS} -DGMX_SOFTWARE_INVSQRT")
 endif(GMX_SOFTWARE_INVSQRT)
 
-#######################################################################
-# Check for options incompatible with OpenMM build                    #
-#######################################################################
-if(GMX_OPENMM)
-    # we'll use the built-in fft to avoid unnecessary dependencies
-    string(TOUPPER ${GMX_FFT_LIBRARY} GMX_FFT_LIBRARY)
-    if(NOT ${GMX_FFT_LIBRARY} STREQUAL "FFTPACK")
-        message(STATUS "No external FFT libraries needed for the OpenMM build, switching to fftpack!")
-        set(GMX_FFT_LIBRARY "fftpack" CACHE STRING 
-               "No external FFT libraries needed for the OpenMM build, switching to  fftpack!" FORCE)
-    endif()
-    if(GMX_MPI)
-        message(FATAL_ERROR "The OpenMM build is not compatible with MPI!")
-    endif(GMX_MPI)
-    if(GMX_THREAD_MPI)
-        message(STATUS "Thread-MPI not compatible with OpenMM, disabled!")
-        set(GMX_THREAD_MPI OFF CACHE BOOL
-               "Thread-MPI not compatible with OpenMM build, disabled!" FORCE)
-    endif(GMX_THREAD_MPI)
-    if(GMX_OPENMP)
-        message(STATUS "OpenMP multithreading not compatible with OpenMM, disabled")
-        set(GMX_OPENMP OFF CACHE BOOL
-            "OpenMP multithreading not compatible with OpenMM, disabled!" FORCE)
-    endif()
-    if(GMX_SOFTWARE_INVSQRT)
-        set(GMX_SOFTWARE_INVSQRT OFF CACHE STRING 
-                "The OpenMM build does not need GROMACS software 1/sqrt!" FORCE)
-    endif(GMX_SOFTWARE_INVSQRT)
-    string(TOUPPER ${GMX_CPU_ACCELERATION} GMX_CPU_ACCELERATION)
-    if(NOT GMX_CPU_ACCELERATION STREQUAL "NONE")
-        message(STATUS "Switching off CPU-based acceleration, the OpenMM build does not support/need any!")
-        set(GMX_CPU_ACCELERATION "None" CACHE STRING
-            "Switching off CPU-based acceleration, the OpenMM build does not support/need any!" FORCE)
-    endif()
-    if(GMX_FAHCORE)
-        message(FATAL_ERROR "The OpenMM build does not support FAH build!")
-    endif(GMX_FAHCORE)
-    if(GMX_DOUBLE)
-        message(FATAL_ERROR  "The OpenMM-build does not support double precision calculations!")
-    endif()
-    # mark as advanced the unused variables
-    mark_as_advanced(FORCE GMX_CPU_ACCELERATION GMX_MPI GMX_FFT_LIBRARY 
-        GMX_QMMM_PROGRAM GMX_THREAD_MPI GMX_DOUBLE)
-else(GMX_OPENMM)
-     mark_as_advanced(CLEAR GMX_CPU_ACCELERATION GMX_MPI GMX_FFT_LIBRARY 
-        GMX_QMMM_PROGRAM GMX_THREAD_MPI GMX_DOUBLE)
-endif(GMX_OPENMM)
 
 
 ########################################################################
@@ -609,16 +554,6 @@ else(GMX_THREAD_MPI)
     tmpi_get_source_list(THREAD_MPI_SRC NOMPI)
 endif(GMX_THREAD_MPI)
 
-if(GMX_OPENMM)
-    set(CUDA_BUILD_EMULATION OFF)
-    find_package(CUDA 3.1 REQUIRED)
-    add_definitions(-DGMX_OPENMM)
-    if(CMAKE_BUILD_TYPE STREQUAL "DEBUG")    
-        set(CUDA_VERBOSE_BUILD ON)
-    endif()
-    find_package(OpenMM) 
-endif(GMX_OPENMM)
-
 if(GMX_GPU)
     # now that we have detected the dependencies, do the second configure pass
     gmx_gpu_setup()
@@ -940,7 +875,7 @@ else(${GMX_QMMM_PROGRAM} STREQUAL "GAUSSIAN")
     MESSAGE(FATAL_ERROR "Invalid QM/MM program option: ${GMX_QMMM_PROGRAM}. Choose one of: Gaussian, Mopac, Gamess, Orca, None")
 endif(${GMX_QMMM_PROGRAM} STREQUAL "GAUSSIAN")
 
-# Process FFT library settings - if not OpenMM build 
+# Process FFT library settings
 string(TOUPPER ${GMX_FFT_LIBRARY} ${GMX_FFT_LIBRARY})
 set(PKG_FFT "")
 set(PKG_FFT_LIBS "")
diff --git a/COPYING b/COPYING
index fb5e6371e5fdf7a98aa8d06d2f3b8306a8790ec6..43c15989e4edd30cb685f8b5500ae209203c93e2 100644 (file)
--- a/COPYING
+++ b/COPYING
@@ -15,7 +15,6 @@ This file contains the licenses for the following bodies of code:
 5. thread_mpi
 6. Blas
 7. Lapack
-8. OpenMM
 
 Our chosen method for packaging distributions (CPack) only permits a
 package to have a single license file, so we are unfortunately forced
@@ -1071,9 +1070,3 @@ in Gromacs (primarily full & sparse matrix diagonalization), so in most cases it
 better idea to use the full reference implementation.
 
 Erik Lindahl, 2008-10-07.
-
-8. OpenMM (binary distributions only)
-=====================================
-
-There are several licenses which cover different parts of OpenMM as described
-in the file openmm/licenses/Licenses.txt accompanying the binary distribution.
diff --git a/INSTALL-GPU b/INSTALL-GPU
deleted file mode 100644 (file)
index 3b9326a..0000000
+++ /dev/null
@@ -1,73 +0,0 @@
-
-                    Welcome to GROMACS-GPU!
-
-Note: Detailed, step-by-step installation instructions as 
-well as additional information are available on the website 
-http://www.gromacs.org/gpu
-
-
-* INSTALLING FROM BINARY DISTRIBUTION:
-
-0. Prerequisites: 
-    - OpenMM (included in the binary release)
-    - NVIDIA CUDA libraries (version >=3.1);
-    - NVIDIA driver (for details on compatiblity consult
-         http://www.nvidia.com/Download/index5.aspx);
-    - NVIDIA CUDA-enabled GPU (for compatiblity list see 
-         http://www.gromacs.org/gpu).
-
-
-1. Download and unpack the binary package for the respective 
-OS and architecture. Copy the content of the package to your 
-normal GROMACS installation directory (or to a custom location). 
-Note that as the distributed Gromacs-GPU packages do not contain 
-the entire set of tools and utilities included in a full GROMACS 
-installation. Therefore, it is recommended to have a ≥v4.5 
-standard Gromacs installation along the GPU accelerated one.
-e.g. on unix: 
-
-tar -xvf gromacs-4.5-GPU.tar.gz
-cp -R gromacs-4.5-GPU/* PATH_TO_GROMACS_INSTALLATION
-
-
-2. Add the openmm/lib directory to your library path, e.g. in bash:
-
-export LD_LIBRARY_PATH=PATH_TO_GROMACS/openmm/lib:$LD_LIBRARY_PATH
-
-If there are other OpenMM versions installed, make sure that the 
-supplied libraries have preference when running mdrun-gpu. 
-Also, make sure that the installed CUDA libraries match the version 
-of CUDA with which GROMACS-GPU has been compiled.
-
-
-3. Set the OPENMM_PLUGIN_DIR environment variable to contain the path 
-to the openmm/lib/plugins directory, e.g. in bash:
-
-export OPENMM_PLUGIN_DIR=PATH_TO_GROMACS/openmm/lib/plugins
-
-
-4. At this point, running the command: 
-
-PATH_TO_GROMACS/bin/mdrun-gpu -h 
-
-should display the standard mdrun help which means that all the 
-necessary libraries are accessible.
-
-
-
-* INSTALLING FROM SOURCE DISTRIBUTION:
-
-GROMACS-GPU uses a cmake build-generator and makefiles on unix. 
-All you have to do is run:
-
-cmake PATH_TO_SOURCE_DIRECTORY -DGMX_OPENMM=ON -DGMX_THREADS=OFF 
-make mdrun
-make install-mdrun
-
-
-
-* ARE YOU STILL HAVING PROBLEMS?
-
-Post it to the GROMACS mailing lists. We read these on a regular basis,
-and in many cases another user might already be familiar with the
-task you're trying to perform!
diff --git a/README-GPU b/README-GPU
deleted file mode 100644 (file)
index 057f721..0000000
+++ /dev/null
@@ -1,21 +0,0 @@
-
-               Welcome to the official version of GROMACS-GPU!
-
-
-GROMACS-GPU provides GPU-accelerated MD simulations.
-The current release uses the OpenMM library. 
-
-This release is targeted at developers and advanced users,
-and care should be taken before production use. 
-
-Installation instructions are available in the INSTALL-GPU file. 
-More details on the GPU accelerated release are available on the 
-GROMACS-GPU website: http://www.gromacs.org/gpu.
-
-See also the README file for general GROMACS information.
-
-                               * * * * *
-
-GROMACS is free software, distributed under the GNU General Public License.
-The GROMACS-GPU release however includes / relies on code covered by several 
-different licences. For details see the COPYING-GPU file.
index 34320c41ec5118241beff5ddf70d71bdf2e0a108..16d09c7e4acbee4b7e9a5c725c95043c3ea92c15 100644 (file)
@@ -94,10 +94,14 @@ version is encouraged. \nvidia{} GPUs with at least \nvidia{} compute
 capability 2.0 are required, e.g. Fermi or Kepler cards.
 
 The GPU support from \gromacs{} version 4.5 using \openmm{}
-\url{https://simtk.org/home/openmm} is still available, also requires
-\cuda{}, and remains the only hardware-based acceleration available
-for implicit solvent simulations in \gromacs{}. This parallelization
-path may not be maintained in the future.
+\url{https://simtk.org/home/openmm} is still contained in the code,
+but in the ``user contributions'' section (\texttt{src/contrib}). You
+will need to edit \texttt{src/contrib/CMakeLists.txt} to enable it. It
+also requires \cuda{}, and remains the only hardware-based
+acceleration available for implicit solvent simulations in
+\gromacs{}. There are no plans to maintain this implementation for
+this parallelization path in the future, and it may not even build
+now.
 
 If you wish to run in parallel on multiple machines across a network,
 you will need to have
index 91cac1c8a074763823967a4459f8d10e47865ea8..3024d156e66dbaebae322c8d5fcff75a716e5f81 100644 (file)
@@ -65,9 +65,6 @@ int do_full_memtest(int dev_id) FUNC_TERM_INT
 FUNC_QUALIFIER
 int do_timed_memtest(int dev_id, int time_limit) FUNC_TERM_INT
 
-FUNC_QUALIFIER
-gmx_bool is_gmx_openmm_supported_gpu(int dev_id, char *gpu_name) FUNC_TERM_INT
-
 FUNC_QUALIFIER
 int detect_cuda_gpus(gmx_gpu_info_t *gpu_info, char *err_str) FUNC_TERM_INT
 
index 3908b574df71b50de426632a52541d1b48f7e4b0..1d0adc5143aff3f5b0121bc618c12617265aa470 100644 (file)
@@ -144,9 +144,6 @@ typedef double gmx_integrator_t(FILE *log,t_commrec *cr,
 
 gmx_integrator_t do_md;
 
-gmx_integrator_t do_md_openmm;
-
-
 
 /* ROUTINES from minimize.c */
 
diff --git a/src/contrib/BuildMdrunOpenMM.cmake b/src/contrib/BuildMdrunOpenMM.cmake
new file mode 100644 (file)
index 0000000..52aaec6
--- /dev/null
@@ -0,0 +1,63 @@
+#
+# This file is part of the GROMACS molecular simulation package.
+#
+# Copyright (c) 2012, by the GROMACS development team, led by
+# David van der Spoel, Berk Hess, Erik Lindahl, and including many
+# others, as listed in the AUTHORS file in the top-level source
+# directory and at http://www.gromacs.org.
+#
+# GROMACS is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public License
+# as published by the Free Software Foundation; either version 2.1
+# of the License, or (at your option) any later version.
+#
+# GROMACS is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with GROMACS; if not, see
+# http://www.gnu.org/licenses, or write to the Free Software Foundation,
+# Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+#
+# If you want to redistribute modifications to GROMACS, please
+# consider that scientific software is very special. Version
+# control is crucial - bugs must be traceable. We will be happy to
+# consider code for inclusion in the official distribution, but
+# derived work must not be called official GROMACS. Details are found
+# in the README & COPYING files - if they are missing, get the
+# official version at http://www.gromacs.org.
+#
+# To help us fund GROMACS development, we humbly ask that you cite
+# the research papers on the package. Check out http://www.gromacs.org.
+#
+
+include_directories(${OpenMM_INCLUDE_DIR})
+link_directories(${OpenMM_LIBRARY_DIR}) 
+# with this define no evn.var. is needed with OPENMM_PLUGIN_DIR
+# if the same OpenMM installation is used for running and building 
+add_definitions( -DOPENMM_PLUGIN_DIR="${OpenMM_PLUGIN_DIR}" ) 
+file(TO_CMAKE_PATH ${OpenMM_PLUGIN_DIR} _path)
+add_library(openmm_api_wrapper STATIC openmm_wrapper.cpp)
+target_link_libraries(openmm_api_wrapper ${OpenMM_LIBRARIES})
+list(APPEND GMX_EXTRA_LIBRARIES openmm_api_wrapper ${OpenMM_LIBRARIES})   
+
+list(REMOVE_ITEM MDRUN_SOURCES mdrun.c)
+list(APPEND MDRUN_SOURCES
+    ${CMAKE_CURRENT_SOURCE_DIR}/src/contrib/md_openmm.c
+    ${CMAKE_CURRENT_SOURCE_DIR}/src/contrib/mdrun_openmm.c)
+
+# this is to circumvent the following MSVC error: 
+# warning LNK4098: defaultlib 'LIBCMT' conflicts with use of other libs
+# fatal error LNK1169: one or more multiply defined symbols found
+if(GMX_OPENMM AND MSVC)
+    set_target_properties(mdrun PROPERTIES LINK_FLAGS "/NODEFAULTLIB:LIBCMT")
+endif()
+
+include_directories(${CMAKE_SOURCE_DIR}/src/gmxlib/gpu_utils)
+
+set_source_files_properties(main.c PROPERTIES LANGUAGE CXX)
+if ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang")
+    set_source_files_properties(main.c PROPERTIES COMPILE_FLAGS "-x c++")
+endif()
index ed3084afa828c92b2754386dbc2d80e00b83318d..6e3d99a73d6a25a95db37eb02566ec2bd18a1082 100644 (file)
 #
 set(CONTRIB_PROGRAMS 
      #add here any programs you want to compile
+
 )
 
+# Uncomment the next line to build OpenMM:
+#option(GMX_OPENMM "Accelerated execution on GPUs through the OpenMM library (no longer supported" ON)
+
+# At run time, you may need to set the environment variable
+# OPENMM_PLUGIN_DIR=PATH_TO_GROMACS/openmm/lib/plugins
+# to make things work
+
+if(GMX_OPENMM)
+    if(GMX_GPU)
+        message(FATAL_ERROR "The OpenMM build is not compatible with the native GPU build")
+    endif()
+
+    enable_language(CXX)
+    set (GMX_BINARY_SUFFIX "-openmm")
+    set (GMX_LIBS_SUFFIX "_openmm")
+
+#######################################################################
+# Check for options incompatible with OpenMM build                    #
+#######################################################################
+    # we'll use the built-in fft to avoid unnecessary dependencies
+    string(TOUPPER ${GMX_FFT_LIBRARY} GMX_FFT_LIBRARY)
+    if(NOT ${GMX_FFT_LIBRARY} STREQUAL "FFTPACK")
+        message(STATUS "No external FFT libraries needed for the OpenMM build, switching to fftpack!")
+        set(GMX_FFT_LIBRARY "fftpack" CACHE STRING 
+               "No external FFT libraries needed for the OpenMM build, switching to  fftpack!" FORCE)
+    endif()
+    if(GMX_MPI)
+        message(FATAL_ERROR "The OpenMM build is not compatible with MPI!")
+    endif(GMX_MPI)
+    if(GMX_THREAD_MPI)
+        message(STATUS "Thread-MPI not compatible with OpenMM, disabled!")
+        set(GMX_THREAD_MPI OFF CACHE BOOL
+               "Thread-MPI not compatible with OpenMM build, disabled!" FORCE)
+    endif(GMX_THREAD_MPI)
+    if(GMX_OPENMP)
+        message(STATUS "OpenMP multithreading not compatible with OpenMM, disabled")
+        set(GMX_OPENMP OFF CACHE BOOL
+            "OpenMP multithreading not compatible with OpenMM, disabled!" FORCE)
+    endif()
+    if(GMX_SOFTWARE_INVSQRT)
+        set(GMX_SOFTWARE_INVSQRT OFF CACHE STRING 
+                "The OpenMM build does not need GROMACS software 1/sqrt!" FORCE)
+    endif(GMX_SOFTWARE_INVSQRT)
+    string(TOUPPER ${GMX_CPU_ACCELERATION} GMX_CPU_ACCELERATION)
+    if(NOT GMX_CPU_ACCELERATION STREQUAL "NONE")
+        message(STATUS "Switching off CPU-based acceleration, the OpenMM build does not support/need any!")
+        set(GMX_CPU_ACCELERATION "None" CACHE STRING
+            "Switching off CPU-based acceleration, the OpenMM build does not support/need any!" FORCE)
+    endif()
+    if(GMX_FAHCORE)
+        message(FATAL_ERROR "The OpenMM build does not support FAH build!")
+    endif(GMX_FAHCORE)
+    if(GMX_DOUBLE)
+        message(FATAL_ERROR  "The OpenMM-build does not support double precision calculations!")
+    endif()
+
+    set(CUDA_BUILD_EMULATION OFF)
+    find_package(CUDA 3.1 REQUIRED)
+    add_definitions(-DGMX_OPENMM)
+    if(CMAKE_BUILD_TYPE STREQUAL "DEBUG")    
+        set(CUDA_VERBOSE_BUILD ON)
+    endif()
+    list(APPEND CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR}/src/contrib)
+    find_package(OpenMM) 
+
+    # mark as advanced the unused variables
+    mark_as_advanced(FORCE GMX_CPU_ACCELERATION GMX_MPI GMX_FFT_LIBRARY 
+        GMX_QMMM_PROGRAM GMX_THREAD_MPI GMX_DOUBLE)
+
+else(GMX_OPENMM)
+
+     mark_as_advanced(CLEAR GMX_CPU_ACCELERATION GMX_MPI GMX_FFT_LIBRARY 
+        GMX_QMMM_PROGRAM GMX_THREAD_MPI GMX_DOUBLE)
+
+endif(GMX_OPENMM)
+
 foreach(PROG ${CONTRIB_PROGRAMS})
         add_executable(${PROG} ${PROG}.c ${NGMX_COMMON_SOURCE})
        set_target_properties(${PROG} PROPERTIES OUTPUT_NAME "${PROG}${GMX_BINARY_SUFFIX}")
diff --git a/src/contrib/mdrun_openmm.c b/src/contrib/mdrun_openmm.c
new file mode 100644 (file)
index 0000000..53bba75
--- /dev/null
@@ -0,0 +1,509 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
+ * Copyright (c) 2001-2004, The GROMACS development team,
+ * check out http://www.gromacs.org for more information.
+ * Copyright (c) 2012, by the GROMACS development team, led by
+ * David van der Spoel, Berk Hess, Erik Lindahl, and including many
+ * others, as listed in the AUTHORS file in the top-level source
+ * directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include "typedefs.h"
+#include "macros.h"
+#include "copyrite.h"
+#include "main.h"
+#include "statutil.h"
+#include "smalloc.h"
+#include "futil.h"
+#include "smalloc.h"
+#include "edsam.h"
+#include "mdrun.h"
+#include "xmdrun.h"
+#include "checkpoint.h"
+#ifdef GMX_THREAD_MPI
+#include "thread_mpi.h"
+#endif
+
+/* afm stuf */
+#include "pull.h"
+
+int cmain(int argc,char *argv[])
+{
+  const char *desc[] = {
+    "This is an experimental release of GROMACS for accelerated",
+       "Molecular Dynamics simulations on GPU processors. Support is provided",
+       "by the OpenMM library (https://simtk.org/home/openmm).[PAR]",
+       "*Warning*[BR]",
+       "This release is targeted at developers and advanced users and",
+       "care should be taken before production use. The following should be",
+       "noted before using the program:[PAR]",
+       " * The current release runs only on modern nVidia GPU hardware with CUDA support.",
+       "Make sure that the necessary CUDA drivers and libraries for your operating system",
+       "are already installed. The CUDA SDK also should be installed in order to compile",
+       "the program from source (http://www.nvidia.com/object/cuda_home.html).[PAR]",
+       " * Multiple GPU cards are not supported.[PAR]",
+       " * Only a small subset of the GROMACS features and options are supported on the GPUs.",
+       "See below for a detailed list.[PAR]",
+       " * Consumer level GPU cards are known to often have problems with faulty memory.",
+       "It is recommended that a full memory check of the cards is done at least once",
+       "(for example, using the memtest=full option).",
+       "A partial memory check (for example, memtest=15) before and",
+       "after the simulation run would help spot",
+       "problems resulting from processor overheating.[PAR]",
+       " * The maximum size of the simulated systems depends on the available",
+       "GPU memory,for example, a GTX280 with 1GB memory has been tested with systems",
+       "of up to about 100,000 atoms.[PAR]",
+       " * In order to take a full advantage of the GPU platform features, many algorithms",
+       "have been implemented in a very different way than they are on the CPUs.",
+       "Therefore numercal correspondence between properties of the state of",
+       "simulated systems should not be expected. Moreover, the values will likely vary",
+       "when simulations are done on different GPU hardware.[PAR]",
+       " * Frequent retrieval of system state information such as",
+       "trajectory coordinates and energies can greatly influence the performance",
+       "of the program due to slow CPU<->GPU memory transfer speed.[PAR]",
+       " * MD algorithms are complex, and although the Gromacs code is highly tuned for them,",
+       "they often do not translate very well onto the streaming architetures.",
+       "Realistic expectations about the achievable speed-up from test with GTX280:",
+       "For small protein systems in implicit solvent using all-vs-all kernels the acceleration",
+       "can be as high as 20 times, but in most other setups involving cutoffs and PME the",
+       "acceleration is usually only ~4 times relative to a 3GHz CPU.[PAR]",
+       "Supported features:[PAR]",
+       " * Integrators: md/md-vv/md-vv-avek, sd/sd1 and bd.\n",
+       " * Long-range interactions (option coulombtype): Reaction-Field, Ewald, PME, and cut-off (for Implicit Solvent only)\n",
+       " * Temperature control: Supported only with the md/md-vv/md-vv-avek, sd/sd1 and bd integrators.\n",
+       " * Pressure control: Supported.\n",
+       " * Implicit solvent: Supported.\n",
+       "A detailed description can be found on the GROMACS website:\n",
+       "http://www.gromacs.org/gpu[PAR]",
+/* From the original mdrun documentaion */
+    "The [TT]mdrun[tt] program reads the run input file ([TT]-s[tt])",
+    "and distributes the topology over nodes if needed.",
+    "[TT]mdrun[tt] produces at least four output files.",
+    "A single log file ([TT]-g[tt]) is written, unless the option",
+    "[TT]-seppot[tt] is used, in which case each node writes a log file.",
+    "The trajectory file ([TT]-o[tt]), contains coordinates, velocities and",
+    "optionally forces.",
+    "The structure file ([TT]-c[tt]) contains the coordinates and",
+    "velocities of the last step.",
+    "The energy file ([TT]-e[tt]) contains energies, the temperature,",
+    "pressure, etc, a lot of these things are also printed in the log file.",
+    "Optionally coordinates can be written to a compressed trajectory file",
+    "([TT]-x[tt]).[PAR]",
+/* openmm specific information */
+       "Usage with OpenMM:[BR]",
+       "[TT]mdrun -device \"OpenMM:platform=Cuda,memtest=15,deviceid=0,force-device=no\"[tt][PAR]",
+       "Options:[PAR]",
+       "      [TT]platform[tt] = Cuda\t\t:\tThe only available value. OpenCL support will be available in future.\n",
+       "      [TT]memtest[tt] = 15\t\t:\tRun a partial, random GPU memory test for the given amount of seconds. A full test",
+       "(recommended!) can be run with \"memtest=full\". Memory testing can be disabled with \"memtest=off\".\n",
+       "      [TT]deviceid[tt] = 0\t\t:\tSpecify the target device when multiple cards are present.",
+       "Only one card can be used at any given time though.\n",
+       "      [TT]force-device[tt] = no\t\t:\tIf set to \"yes\" [TT]mdrun[tt]  will be forced to execute on",
+       "hardware that is not officially supported. GPU acceleration can also be achieved on older",
+       "but Cuda capable cards, although the simulation might be too slow, and the memory limits too strict.",
+  };
+  t_commrec    *cr;
+  t_filenm fnm[] = {
+    { efTPX, NULL,      NULL,       ffREAD },
+    { efTRN, "-o",      NULL,       ffWRITE },
+    { efXTC, "-x",      NULL,       ffOPTWR },
+    { efCPT, "-cpi",    NULL,       ffOPTRD },
+    { efCPT, "-cpo",    NULL,       ffOPTWR },
+    { efSTO, "-c",      "confout",  ffWRITE },
+    { efEDR, "-e",      "ener",     ffWRITE },
+    { efLOG, "-g",      "md",       ffWRITE },
+    { efXVG, "-dhdl",   "dhdl",     ffOPTWR },
+    { efXVG, "-field",  "field",    ffOPTWR },
+    { efXVG, "-table",  "table",    ffOPTRD },
+    { efXVG, "-tabletf", "tabletf",    ffOPTRD },
+    { efXVG, "-tablep", "tablep",   ffOPTRD },
+    { efXVG, "-tableb", "table",    ffOPTRD },
+    { efTRX, "-rerun",  "rerun",    ffOPTRD },
+    { efXVG, "-tpi",    "tpi",      ffOPTWR },
+    { efXVG, "-tpid",   "tpidist",  ffOPTWR },
+    { efEDI, "-ei",     "sam",      ffOPTRD },
+    { efEDO, "-eo",     "sam",      ffOPTWR },
+    { efGCT, "-j",      "wham",     ffOPTRD },
+    { efGCT, "-jo",     "bam",      ffOPTWR },
+    { efXVG, "-ffout",  "gct",      ffOPTWR },
+    { efXVG, "-devout", "deviatie", ffOPTWR },
+    { efXVG, "-runav",  "runaver",  ffOPTWR },
+    { efXVG, "-px",     "pullx",    ffOPTWR },
+    { efXVG, "-pf",     "pullf",    ffOPTWR },
+    { efXVG, "-ro",     "rotation", ffOPTWR },
+    { efLOG, "-ra",     "rotangles",ffOPTWR },
+    { efLOG, "-rs",     "rotslabs", ffOPTWR },
+    { efLOG, "-rt",     "rottorque",ffOPTWR },
+    { efMTX, "-mtx",    "nm",       ffOPTWR },
+    { efNDX, "-dn",     "dipole",   ffOPTWR },
+    { efRND, "-multidir",NULL,      ffOPTRDMULT},
+    { efDAT, "-membed", "membed",   ffOPTRD },
+    { efTOP, "-mp",     "membed",   ffOPTRD },
+    { efNDX, "-mn",     "membed",   ffOPTRD }
+  };
+#define NFILE asize(fnm)
+
+  /* Command line options ! */
+  gmx_bool bCart        = FALSE;
+  gmx_bool bPPPME       = FALSE;
+  gmx_bool bPartDec     = FALSE;
+  gmx_bool bDDBondCheck = TRUE;
+  gmx_bool bDDBondComm  = TRUE;
+  gmx_bool bTunePME     = TRUE;
+  gmx_bool bTestVerlet  = FALSE;
+  gmx_bool bVerbose     = FALSE;
+  gmx_bool bCompact     = TRUE;
+  gmx_bool bSepPot      = FALSE;
+  gmx_bool bRerunVSite  = FALSE;
+  gmx_bool bIonize      = FALSE;
+  gmx_bool bConfout     = TRUE;
+  gmx_bool bReproducible = FALSE;
+    
+  int  npme=-1;
+  int  nmultisim=0;
+  int  nstglobalcomm=-1;
+  int  repl_ex_nst=0;
+  int  repl_ex_seed=-1;
+  int  repl_ex_nex=0;
+  int  nstepout=100;
+  int  resetstep=-1;
+  int  nsteps=-2; /* the value -2 means that the mdp option will be used */
+  
+  rvec realddxyz={0,0,0};
+  const char *ddno_opt[ddnoNR+1] =
+    { NULL, "interleave", "pp_pme", "cartesian", NULL };
+  const char *dddlb_opt[] =
+    { NULL, "auto", "no", "yes", NULL };
+  const char *nbpu_opt[] =
+    { NULL, "auto", "cpu", "gpu", "gpu_cpu", NULL };
+  real rdd=0.0,rconstr=0.0,dlb_scale=0.8,pforce=-1;
+  char *ddcsx=NULL,*ddcsy=NULL,*ddcsz=NULL;
+  real cpt_period=15.0,max_hours=-1;
+  gmx_bool bAppendFiles=TRUE;
+  gmx_bool bKeepAndNumCPT=FALSE;
+  gmx_bool bResetCountersHalfWay=FALSE;
+  output_env_t oenv=NULL;
+  const char *deviceOptions = "";
+
+  gmx_hw_opt_t hw_opt={0,0,0,0,TRUE,FALSE,0,NULL};
+
+  t_pargs pa[] = {
+
+    { "-pd",      FALSE, etBOOL,{&bPartDec},
+      "Use particle decompostion" },
+    { "-dd",      FALSE, etRVEC,{&realddxyz},
+      "Domain decomposition grid, 0 is optimize" },
+    { "-ddorder", FALSE, etENUM, {ddno_opt},
+      "DD node order" },
+    { "-npme",    FALSE, etINT, {&npme},
+      "Number of separate nodes to be used for PME, -1 is guess" },
+    { "-nt",      FALSE, etINT, {&hw_opt.nthreads_tot},
+      "Total number of threads to start (0 is guess)" },
+    { "-ntmpi",   FALSE, etINT, {&hw_opt.nthreads_tmpi},
+      "Number of thread-MPI threads to start (0 is guess)" },
+    { "-ntomp",   FALSE, etINT, {&hw_opt.nthreads_omp},
+      "Number of OpenMP threads per MPI process/thread to start (0 is guess)" },
+    { "-ntomp_pme", FALSE, etINT, {&hw_opt.nthreads_omp_pme},
+      "Number of OpenMP threads per MPI process/thread to start (0 is -ntomp)" },
+    { "-pin",     FALSE, etBOOL, {&hw_opt.bThreadPinning},
+      "Pin OpenMP threads to cores" },
+    { "-pinht",   FALSE, etBOOL, {&hw_opt.bPinHyperthreading},
+      "Always pin threads to Hyper-Threading cores" },
+    { "-pinoffset", FALSE, etINT, {&hw_opt.core_pinning_offset},
+      "Core offset for pinning (for running multiple mdrun processes on a single physical node)" },
+    { "-gpu_id",  FALSE, etSTR, {&hw_opt.gpu_id},
+      "List of GPU id's to use" },
+    { "-ddcheck", FALSE, etBOOL, {&bDDBondCheck},
+      "Check for all bonded interactions with DD" },
+    { "-ddbondcomm", FALSE, etBOOL, {&bDDBondComm},
+      "HIDDENUse special bonded atom communication when [TT]-rdd[tt] > cut-off" },
+    { "-rdd",     FALSE, etREAL, {&rdd},
+      "The maximum distance for bonded interactions with DD (nm), 0 is determine from initial coordinates" },
+    { "-rcon",    FALSE, etREAL, {&rconstr},
+      "Maximum distance for P-LINCS (nm), 0 is estimate" },
+    { "-dlb",     FALSE, etENUM, {dddlb_opt},
+      "Dynamic load balancing (with DD)" },
+    { "-dds",     FALSE, etREAL, {&dlb_scale},
+      "Minimum allowed dlb scaling of the DD cell size" },
+    { "-ddcsx",   FALSE, etSTR, {&ddcsx},
+      "HIDDENThe DD cell sizes in x" },
+    { "-ddcsy",   FALSE, etSTR, {&ddcsy},
+      "HIDDENThe DD cell sizes in y" },
+    { "-ddcsz",   FALSE, etSTR, {&ddcsz},
+      "HIDDENThe DD cell sizes in z" },
+    { "-gcom",    FALSE, etINT,{&nstglobalcomm},
+      "Global communication frequency" },
+    { "-nb",      FALSE, etENUM, {&nbpu_opt},
+      "Calculate non-bonded interactions on" },
+    { "-tunepme", FALSE, etBOOL, {&bTunePME},  
+      "Optimize PME load between PP/PME nodes or GPU/CPU" },
+    { "-testverlet", FALSE, etBOOL, {&bTestVerlet},
+      "Test the Verlet non-bonded scheme" },
+    { "-v",       FALSE, etBOOL,{&bVerbose},  
+      "Be loud and noisy" },
+    { "-compact", FALSE, etBOOL,{&bCompact},  
+      "Write a compact log file" },
+    { "-seppot",  FALSE, etBOOL, {&bSepPot},
+      "Write separate V and dVdl terms for each interaction type and node to the log file(s)" },
+    { "-pforce",  FALSE, etREAL, {&pforce},
+      "Print all forces larger than this (kJ/mol nm)" },
+    { "-reprod",  FALSE, etBOOL,{&bReproducible},  
+      "Try to avoid optimizations that affect binary reproducibility" },
+    { "-cpt",     FALSE, etREAL, {&cpt_period},
+      "Checkpoint interval (minutes)" },
+    { "-cpnum",   FALSE, etBOOL, {&bKeepAndNumCPT},
+      "Keep and number checkpoint files" },
+    { "-append",  FALSE, etBOOL, {&bAppendFiles},
+      "Append to previous output files when continuing from checkpoint instead of adding the simulation part number to all file names" },
+    { "-nsteps",  FALSE, etINT, {&nsteps},
+      "Run this number of steps, overrides .mdp file option" },
+    { "-maxh",   FALSE, etREAL, {&max_hours},
+      "Terminate after 0.99 times this time (hours)" },
+    { "-multi",   FALSE, etINT,{&nmultisim}, 
+      "Do multiple simulations in parallel" },
+    { "-replex",  FALSE, etINT, {&repl_ex_nst}, 
+      "Attempt replica exchange periodically with this period (steps)" },
+    { "-nex",  FALSE, etINT, {&repl_ex_nex},
+      "Number of random exchanges to carry out each exchange interval (N^3 is one suggestion).  -nex zero or not specified gives neighbor replica exchange." },
+    { "-reseed",  FALSE, etINT, {&repl_ex_seed}, 
+      "Seed for replica exchange, -1 is generate a seed" },
+    { "-rerunvsite", FALSE, etBOOL, {&bRerunVSite},
+      "HIDDENRecalculate virtual site coordinates with [TT]-rerun[tt]" },
+    { "-ionize",  FALSE, etBOOL,{&bIonize},
+      "Do a simulation including the effect of an X-Ray bombardment on your system" },
+    { "-confout", FALSE, etBOOL, {&bConfout},
+      "HIDDENWrite the last configuration with [TT]-c[tt] and force checkpointing at the last step" },
+    { "-stepout", FALSE, etINT, {&nstepout},
+      "HIDDENFrequency of writing the remaining runtime" },
+    { "-resetstep", FALSE, etINT, {&resetstep},
+      "HIDDENReset cycle counters after these many time steps" },
+    { "-resethway", FALSE, etBOOL, {&bResetCountersHalfWay},
+      "HIDDENReset the cycle counters after half the number of steps or halfway [TT]-maxh[tt]" },
+    { "-device",  FALSE, etSTR, {&deviceOptions},
+      "Device option string" }
+  };
+  gmx_edsam_t  ed;
+  unsigned long Flags, PCA_Flags;
+  ivec     ddxyz;
+  int      dd_node_order;
+  gmx_bool     bAddPart;
+  FILE     *fplog,*fptest;
+  int      sim_part,sim_part_fn;
+  const char *part_suffix=".part";
+  char     suffix[STRLEN];
+  int      rc;
+  char **multidir=NULL;
+
+
+  cr = init_par(&argc,&argv);
+
+  if (MASTER(cr))
+    CopyRight(stderr, argv[0]);
+
+  PCA_Flags = (PCA_CAN_SET_DEFFNM | (MASTER(cr) ? 0 : PCA_QUIET));
+  
+  /* Comment this in to do fexist calls only on master
+   * works not with rerun or tables at the moment
+   * also comment out the version of init_forcerec in md.c 
+   * with NULL instead of opt2fn
+   */
+  /*
+     if (!MASTER(cr))
+     {
+     PCA_Flags |= PCA_NOT_READ_NODE;
+     }
+     */
+
+  parse_common_args(&argc,argv,PCA_Flags, NFILE,fnm,asize(pa),pa,
+                    asize(desc),desc,0,NULL, &oenv);
+
+
+
+  /* we set these early because they might be used in init_multisystem() 
+     Note that there is the potential for npme>nnodes until the number of
+     threads is set later on, if there's thread parallelization. That shouldn't
+     lead to problems. */ 
+  dd_node_order = nenum(ddno_opt);
+  cr->npmenodes = npme;
+
+  /* now check the -multi and -multidir option */
+  if (opt2bSet("-multidir", NFILE, fnm))
+  {
+      int i;
+      if (nmultisim > 0)
+      {
+          gmx_fatal(FARGS, "mdrun -multi and -multidir options are mutually exclusive.");
+      }
+      nmultisim = opt2fns(&multidir, "-multidir", NFILE, fnm);
+  }
+
+
+  if (repl_ex_nst != 0 && nmultisim < 2)
+      gmx_fatal(FARGS,"Need at least two replicas for replica exchange (option -multi)");
+
+  if (repl_ex_nex < 0)
+      gmx_fatal(FARGS,"Replica exchange number of exchanges needs to be positive");
+
+  if (nmultisim > 1) {
+#ifndef GMX_THREAD_MPI
+    gmx_bool bParFn = (multidir == NULL);
+    init_multisystem(cr, nmultisim, multidir, NFILE, fnm, bParFn);
+#else
+    gmx_fatal(FARGS,"mdrun -multi is not supported with the thread library.Please compile GROMACS with MPI support");
+#endif
+  }
+
+  bAddPart = !bAppendFiles;
+
+  /* Check if there is ANY checkpoint file available */        
+  sim_part    = 1;
+  sim_part_fn = sim_part;
+  if (opt2bSet("-cpi",NFILE,fnm))
+  {
+      if (bSepPot && bAppendFiles)
+      {
+          gmx_fatal(FARGS,"Output file appending is not supported with -seppot");
+      }
+
+      bAppendFiles =
+                read_checkpoint_simulation_part(opt2fn_master("-cpi", NFILE,
+                                                              fnm,cr),
+                                                &sim_part_fn,NULL,cr,
+                                                bAppendFiles,NFILE,fnm,
+                                                part_suffix,&bAddPart);
+      if (sim_part_fn==0 && MASTER(cr))
+      {
+          fprintf(stdout,"No previous checkpoint file present, assuming this is a new run.\n");
+      }
+      else
+      {
+          sim_part = sim_part_fn + 1;
+      }
+
+      if (MULTISIM(cr) && MASTER(cr))
+      {
+          check_multi_int(stdout,cr->ms,sim_part,"simulation part");
+      }
+  } 
+  else
+  {
+      bAppendFiles = FALSE;
+  }
+
+  if (!bAppendFiles)
+  {
+      sim_part_fn = sim_part;
+  }
+
+  if (bAddPart)
+  {
+      /* Rename all output files (except checkpoint files) */
+      /* create new part name first (zero-filled) */
+      sprintf(suffix,"%s%04d",part_suffix,sim_part_fn);
+
+      add_suffix_to_output_names(fnm,NFILE,suffix);
+      if (MASTER(cr))
+      {
+          fprintf(stdout,"Checkpoint file is from part %d, new output files will be suffixed '%s'.\n",sim_part-1,suffix);
+      }
+  }
+
+  Flags = opt2bSet("-rerun",NFILE,fnm) ? MD_RERUN : 0;
+  Flags = Flags | (bSepPot       ? MD_SEPPOT       : 0);
+  Flags = Flags | (bIonize       ? MD_IONIZE       : 0);
+  Flags = Flags | (bPartDec      ? MD_PARTDEC      : 0);
+  Flags = Flags | (bDDBondCheck  ? MD_DDBONDCHECK  : 0);
+  Flags = Flags | (bDDBondComm   ? MD_DDBONDCOMM   : 0);
+  Flags = Flags | (bTunePME      ? MD_TUNEPME      : 0);
+  Flags = Flags | (bTestVerlet   ? MD_TESTVERLET   : 0);
+  Flags = Flags | (bConfout      ? MD_CONFOUT      : 0);
+  Flags = Flags | (bRerunVSite   ? MD_RERUN_VSITE  : 0);
+  Flags = Flags | (bReproducible ? MD_REPRODUCIBLE : 0);
+  Flags = Flags | (bAppendFiles  ? MD_APPENDFILES  : 0); 
+  Flags = Flags | (opt2parg_bSet("-append", asize(pa),pa) ? MD_APPENDFILESSET : 0); 
+  Flags = Flags | (bKeepAndNumCPT ? MD_KEEPANDNUMCPT : 0); 
+  Flags = Flags | (sim_part>1    ? MD_STARTFROMCPT : 0); 
+  Flags = Flags | (bResetCountersHalfWay ? MD_RESETCOUNTERSHALFWAY : 0);
+
+
+  /* We postpone opening the log file if we are appending, so we can 
+     first truncate the old log file and append to the correct position 
+     there instead.  */
+  if ((MASTER(cr) || bSepPot) && !bAppendFiles) 
+  {
+      gmx_log_open(ftp2fn(efLOG,NFILE,fnm),cr,
+                   !bSepPot,Flags & MD_APPENDFILES,&fplog);
+      CopyRight(fplog,argv[0]);
+      please_cite(fplog,"Hess2008b");
+      please_cite(fplog,"Spoel2005a");
+      please_cite(fplog,"Lindahl2001a");
+      please_cite(fplog,"Berendsen95a");
+  }
+  else if (!MASTER(cr) && bSepPot)
+  {
+      gmx_log_open(ftp2fn(efLOG,NFILE,fnm),cr,!bSepPot,Flags,&fplog);
+  }
+  else
+  {
+      fplog = NULL;
+  }
+
+  ddxyz[XX] = (int)(realddxyz[XX] + 0.5);
+  ddxyz[YY] = (int)(realddxyz[YY] + 0.5);
+  ddxyz[ZZ] = (int)(realddxyz[ZZ] + 0.5);
+
+  rc = mdrunner(&hw_opt, fplog,cr,NFILE,fnm,oenv,bVerbose,bCompact,
+                nstglobalcomm, ddxyz,dd_node_order,rdd,rconstr,
+                dddlb_opt[0],dlb_scale,ddcsx,ddcsy,ddcsz,
+                nbpu_opt[0],
+                nsteps,nstepout,resetstep,
+                nmultisim,repl_ex_nst,repl_ex_nex,repl_ex_seed,
+                pforce, cpt_period,max_hours,deviceOptions,Flags);
+
+  gmx_finalize_par();
+
+  if (MULTIMASTER(cr)) {
+      thanx(stderr);
+  }
+
+  /* Log file has to be closed in mdrunner if we are appending to it 
+     (fplog not set here) */
+  if (MASTER(cr) && !bAppendFiles) 
+  {
+      gmx_log_close(fplog);
+  }
+
+  return rc;
+}
+
diff --git a/src/contrib/openmm_gpu_utils.cu b/src/contrib/openmm_gpu_utils.cu
new file mode 100644 (file)
index 0000000..cf10dfb
--- /dev/null
@@ -0,0 +1,145 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
+ * Copyright (c) 2001-2010,2012 The GROMACS development team,
+ * check out http://www.gromacs.org for more information.
+ * Copyright (c) 2012, by the GROMACS development team, led by
+ * David van der Spoel, Berk Hess, Erik Lindahl, and including many
+ * others, as listed in the AUTHORS file in the top-level source
+ * directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+
+#include <stdio.h>
+#include <stdlib.h>
+
+#include "smalloc.h"
+#include "string2.h"
+#include "types/hw_info.h"
+
+#include "openmm_gpu_utils.h"
+#include "../../src/gmxlib/cuda_tools/cudautils.cuh"
+
+// TODO put this list into an external file and include it so that the list is easily accessible
+/*! List of supported GPUs. */
+static const char * const SupportedGPUs[] = {
+    /* GT400 */
+    "Geforce GTX 480",
+    "Geforce GTX 470",
+    "Geforce GTX 465",
+    "Geforce GTX 460",
+
+    "Tesla C2070",
+    "Tesla C2050",
+    "Tesla S2070",
+    "Tesla S2050",
+    "Tesla M2070",
+    "Tesla M2050",
+
+    "Quadro 5000",
+    "Quadro 6000",
+
+    /* GT200 */
+    "Geforce GTX 295",
+    "Geforce GTX 285",
+    "Geforce GTX 280",
+    "Geforce GTX 275",
+    "Geforce GTX 260",
+    "GeForce GTS 250",
+    "GeForce GTS 150",
+
+    "GeForce GTX 285M",
+    "GeForce GTX 280M",
+
+    "Tesla S1070",
+    "Tesla C1060",
+    "Tesla M1060",
+
+    "Quadro FX 5800",
+    "Quadro FX 4800",
+    "Quadro CX",
+    "Quadro Plex 2200 D2",
+    "Quadro Plex 2200 S4",
+
+    /* G90 */
+    "GeForce 9800 G", /* GX2, GTX, GTX+, GT */
+    "GeForce 9800M GTX",
+
+    "Quadro FX 4700",
+    "Quadro Plex 2100 D4"
+};
+
+/*! Number of supported GPUs */
+#define NB_GPUS (sizeof(SupportedGPUs)/sizeof(SupportedGPUs[0]))
+
+FUNC_QUALIFIER
+gmx_bool is_gmx_openmm_supported_gpu(int dev_id, char *gpu_name) FUNC_TERM_INT
+
+/*! 
+ * \brief Checks whether the GPU with the given name is supported in Gromacs-OpenMM.
+ * 
+ * \param[in] gpu_name  the name of the CUDA device
+ * \returns             TRUE if the device is supported, otherwise FALSE
+ */
+static bool is_gmx_openmm_supported_gpu_name(char *gpuName)
+{
+    size_t i;
+    for (i = 0; i < NB_GPUS; i++)
+    {
+        trim(gpuName);
+        if (gmx_strncasecmp(gpuName, SupportedGPUs[i], strlen(SupportedGPUs[i])) == 0)
+            return 1;
+    }
+    return 0;
+}
+
+/*! \brief Checks whether the GPU with the given device id is supported in Gromacs-OpenMM.
+ *
+ * \param[in] dev_id    the device id of the GPU or -1 if the device has already been selected
+ * \param[out] gpu_name Set to contain the name of the CUDA device, if NULL passed, no device name is set. 
+ * \returns             TRUE if the device is supported, otherwise FALSE
+ * 
+ */
+gmx_bool is_gmx_openmm_supported_gpu(int dev_id, char *gpu_name)
+{
+    cudaDeviceProp dev_prop;
+
+    if (debug) fprintf(debug, "Checking compatibility with device #%d, %s\n", dev_id, gpu_name);
+
+    if (do_sanity_checks(dev_id, &dev_prop) != 0)
+        return -1;
+
+    if (gpu_name != NULL)
+    { 
+        strcpy(gpu_name, dev_prop.name);
+    }
+    return is_gmx_openmm_supported_gpu_name(dev_prop.name);
+}
+
+
diff --git a/src/contrib/openmm_gpu_utils.h b/src/contrib/openmm_gpu_utils.h
new file mode 100644 (file)
index 0000000..c5874d6
--- /dev/null
@@ -0,0 +1,64 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
+ * Copyright (c) 2001-2010, The GROMACS development team,
+ * check out http://www.gromacs.org for more information.
+ * Copyright (c) 2012, by the GROMACS development team, led by
+ * David van der Spoel, Berk Hess, Erik Lindahl, and including many
+ * others, as listed in the AUTHORS file in the top-level source
+ * directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+
+#ifndef _OPENMM_GPU_UTILS_H_
+#define _OPENMM_GPU_UTILS_H_
+
+#include "types/simple.h"
+#include "types/hw_info.h"
+
+#define FUNC_TERM_INT {return -1;}
+#define FUNC_TERM_VOID {}
+#define FUNC_QUALIFIER static
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+FUNC_QUALIFIER
+gmx_bool is_gmx_openmm_supported_gpu(int dev_id, char *gpu_name) FUNC_TERM_INT
+
+#ifdef __cplusplus
+}
+#endif
+
+#undef FUNC_TERM_INT
+#undef FUNC_TERM_VOID
+#undef FUNC_QUALIFIER
+
+#endif /* _OPENMM_GPU_UTILS_H_ */
similarity index 76%
rename from src/kernel/openmm_wrapper.h
rename to src/contrib/openmm_wrapper.h
index ed30f71579669da64c7cb71b179dbb7d09736f16..7013426aa40dff399c5045e343588a620c07d1db 100644 (file)
@@ -44,7 +44,6 @@ extern "C"
 {
 #endif
 
-#ifdef GMX_OPENMM
 void* openmm_init(FILE *fplog, const char *platformOptStr,
                     t_inputrec *ir,
                     gmx_mtop_t *top_global, gmx_localtop_t *top,
@@ -60,27 +59,6 @@ void openmm_copy_state(void *data,
                         gmx_bool includePos, gmx_bool includeVel, gmx_bool includeForce, gmx_bool includeEnergy);
 
 void openmm_cleanup(FILE *fplog, void* data);
-#else 
-/* dummy versions of the wrapper functions to enable compilation of 
-   do_md_openmm even when OpenMM is not used */ 
-void* openmm_init(FILE *fplog, const char *platformOptStr,
-                    t_inputrec *ir,
-                    gmx_mtop_t *top_global, gmx_localtop_t *top,
-                    t_mdatoms *mdatoms, t_forcerec *fr, t_state *state){return NULL;}
-
-void openmm_take_one_step(void* data){}
-
-void openmm_take_steps(void* data, int nsteps){}
-
-void openmm_copy_state(void *data,
-                        t_state *state, double *time,
-                        rvec f[], gmx_enerdata_t *enerd,
-                        gmx_bool includePos, gmx_bool includeVel, gmx_bool includeForce, gmx_bool includeEnergy){}
-
-void openmm_cleanup(FILE *fplog, void* data){}
-
-#endif /*GMX_OPENMM*/
-
 
 #ifdef __cplusplus
 } // extern "C"
diff --git a/src/contrib/runner_openmm.c b/src/contrib/runner_openmm.c
new file mode 100644 (file)
index 0000000..a04b618
--- /dev/null
@@ -0,0 +1,1998 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
+ * Copyright (c) 2001-2004, The GROMACS development team,
+ * check out http://www.gromacs.org for more information.
+ * Copyright (c) 2012, by the GROMACS development team, led by
+ * David van der Spoel, Berk Hess, Erik Lindahl, and including many
+ * others, as listed in the AUTHORS file in the top-level source
+ * directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#if defined(HAVE_SCHED_H) && defined(HAVE_SCHED_GETAFFINITY)
+#define _GNU_SOURCE
+#include <sched.h>
+#include <sys/syscall.h>
+#endif
+#include <signal.h>
+#include <stdlib.h>
+#include <string.h>
+#include <assert.h>
+
+#include "typedefs.h"
+#include "smalloc.h"
+#include "sysstuff.h"
+#include "statutil.h"
+#include "mdrun.h"
+#include "md_logging.h"
+#include "md_support.h"
+#include "network.h"
+#include "pull.h"
+#include "names.h"
+#include "disre.h"
+#include "orires.h"
+#include "pme.h"
+#include "mdatoms.h"
+#include "repl_ex.h"
+#include "qmmm.h"
+#include "mpelogging.h"
+#include "domdec.h"
+#include "partdec.h"
+#include "coulomb.h"
+#include "constr.h"
+#include "mvdata.h"
+#include "checkpoint.h"
+#include "mtop_util.h"
+#include "sighandler.h"
+#include "tpxio.h"
+#include "txtdump.h"
+#include "gmx_detect_hardware.h"
+#include "gmx_omp_nthreads.h"
+#include "pull_rotation.h"
+#include "calc_verletbuf.h"
+#include "../mdlib/nbnxn_search.h"
+#include "../mdlib/nbnxn_consts.h"
+#include "gmx_fatal_collective.h"
+#include "membed.h"
+#include "md_openmm.h"
+#include "gmx_omp.h"
+
+#include "thread_mpi/threads.h"
+
+#ifdef GMX_LIB_MPI
+#include <mpi.h>
+#endif
+#ifdef GMX_THREAD_MPI
+#include "tmpi.h"
+#endif
+
+#ifdef GMX_FAHCORE
+#include "corewrap.h"
+#endif
+
+#include "gpu_utils.h"
+#include "nbnxn_cuda_data_mgmt.h"
+
+typedef struct { 
+    gmx_integrator_t *func;
+} gmx_intp_t;
+
+/* The array should match the eI array in include/types/enums.h */
+const gmx_intp_t integrator[eiNR] = { {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm},{do_md_openmm}};
+
+gmx_large_int_t     deform_init_init_step_tpx;
+matrix              deform_init_box_tpx;
+#ifdef GMX_THREAD_MPI
+tMPI_Thread_mutex_t deform_init_box_mutex=TMPI_THREAD_MUTEX_INITIALIZER;
+#endif
+
+
+#ifdef GMX_THREAD_MPI
+struct mdrunner_arglist
+{
+    gmx_hw_opt_t *hw_opt;
+    FILE *fplog;
+    t_commrec *cr;
+    int nfile;
+    const t_filenm *fnm;
+    output_env_t oenv;
+    gmx_bool bVerbose;
+    gmx_bool bCompact;
+    int nstglobalcomm;
+    ivec ddxyz;
+    int dd_node_order;
+    real rdd;
+    real rconstr;
+    const char *dddlb_opt;
+    real dlb_scale;
+    const char *ddcsx;
+    const char *ddcsy;
+    const char *ddcsz;
+    const char *nbpu_opt;
+    int nsteps_cmdline;
+    int nstepout;
+    int resetstep;
+    int nmultisim;
+    int repl_ex_nst;
+    int repl_ex_nex;
+    int repl_ex_seed;
+    real pforce;
+    real cpt_period;
+    real max_hours;
+    const char *deviceOptions;
+    unsigned long Flags;
+    int ret; /* return value */
+};
+
+
+/* The function used for spawning threads. Extracts the mdrunner() 
+   arguments from its one argument and calls mdrunner(), after making
+   a commrec. */
+static void mdrunner_start_fn(void *arg)
+{
+    struct mdrunner_arglist *mda=(struct mdrunner_arglist*)arg;
+    struct mdrunner_arglist mc=*mda; /* copy the arg list to make sure 
+                                        that it's thread-local. This doesn't
+                                        copy pointed-to items, of course,
+                                        but those are all const. */
+    t_commrec *cr;  /* we need a local version of this */
+    FILE *fplog=NULL;
+    t_filenm *fnm;
+
+    fnm = dup_tfn(mc.nfile, mc.fnm);
+
+    cr = init_par_threads(mc.cr);
+
+    if (MASTER(cr))
+    {
+        fplog=mc.fplog;
+    }
+
+    mda->ret=mdrunner(mc.hw_opt, fplog, cr, mc.nfile, fnm, mc.oenv, 
+                      mc.bVerbose, mc.bCompact, mc.nstglobalcomm, 
+                      mc.ddxyz, mc.dd_node_order, mc.rdd,
+                      mc.rconstr, mc.dddlb_opt, mc.dlb_scale, 
+                      mc.ddcsx, mc.ddcsy, mc.ddcsz,
+                      mc.nbpu_opt,
+                      mc.nsteps_cmdline, mc.nstepout, mc.resetstep,
+                      mc.nmultisim, mc.repl_ex_nst, mc.repl_ex_nex, mc.repl_ex_seed, mc.pforce, 
+                      mc.cpt_period, mc.max_hours, mc.deviceOptions, mc.Flags);
+}
+
+/* called by mdrunner() to start a specific number of threads (including 
+   the main thread) for thread-parallel runs. This in turn calls mdrunner()
+   for each thread. 
+   All options besides nthreads are the same as for mdrunner(). */
+static t_commrec *mdrunner_start_threads(gmx_hw_opt_t *hw_opt, 
+              FILE *fplog,t_commrec *cr,int nfile, 
+              const t_filenm fnm[], const output_env_t oenv, gmx_bool bVerbose,
+              gmx_bool bCompact, int nstglobalcomm,
+              ivec ddxyz,int dd_node_order,real rdd,real rconstr,
+              const char *dddlb_opt,real dlb_scale,
+              const char *ddcsx,const char *ddcsy,const char *ddcsz,
+              const char *nbpu_opt,
+              int nsteps_cmdline, int nstepout,int resetstep,
+              int nmultisim,int repl_ex_nst,int repl_ex_nex, int repl_ex_seed,
+              real pforce,real cpt_period, real max_hours, 
+              const char *deviceOptions, unsigned long Flags)
+{
+    int ret;
+    struct mdrunner_arglist *mda;
+    t_commrec *crn; /* the new commrec */
+    t_filenm *fnmn;
+
+    /* first check whether we even need to start tMPI */
+    if (hw_opt->nthreads_tmpi < 2)
+    {
+        return cr;
+    }
+
+    /* a few small, one-time, almost unavoidable memory leaks: */
+    snew(mda,1);
+    fnmn=dup_tfn(nfile, fnm);
+
+    /* fill the data structure to pass as void pointer to thread start fn */
+    mda->hw_opt=hw_opt;
+    mda->fplog=fplog;
+    mda->cr=cr;
+    mda->nfile=nfile;
+    mda->fnm=fnmn;
+    mda->oenv=oenv;
+    mda->bVerbose=bVerbose;
+    mda->bCompact=bCompact;
+    mda->nstglobalcomm=nstglobalcomm;
+    mda->ddxyz[XX]=ddxyz[XX];
+    mda->ddxyz[YY]=ddxyz[YY];
+    mda->ddxyz[ZZ]=ddxyz[ZZ];
+    mda->dd_node_order=dd_node_order;
+    mda->rdd=rdd;
+    mda->rconstr=rconstr;
+    mda->dddlb_opt=dddlb_opt;
+    mda->dlb_scale=dlb_scale;
+    mda->ddcsx=ddcsx;
+    mda->ddcsy=ddcsy;
+    mda->ddcsz=ddcsz;
+    mda->nbpu_opt=nbpu_opt;
+    mda->nsteps_cmdline=nsteps_cmdline;
+    mda->nstepout=nstepout;
+    mda->resetstep=resetstep;
+    mda->nmultisim=nmultisim;
+    mda->repl_ex_nst=repl_ex_nst;
+    mda->repl_ex_nex=repl_ex_nex;
+    mda->repl_ex_seed=repl_ex_seed;
+    mda->pforce=pforce;
+    mda->cpt_period=cpt_period;
+    mda->max_hours=max_hours;
+    mda->deviceOptions=deviceOptions;
+    mda->Flags=Flags;
+
+    /* now spawn new threads that start mdrunner_start_fn(), while 
+       the main thread returns */
+    ret=tMPI_Init_fn(TRUE, hw_opt->nthreads_tmpi,
+                     (hw_opt->bThreadPinning ? TMPI_AFFINITY_ALL_CORES : TMPI_AFFINITY_NONE),
+                     mdrunner_start_fn, (void*)(mda) );
+    if (ret!=TMPI_SUCCESS)
+        return NULL;
+
+    /* make a new comm_rec to reflect the new situation */
+    crn=init_par_threads(cr);
+    return crn;
+}
+
+
+static int get_tmpi_omp_thread_division(const gmx_hw_info_t *hwinfo,
+                                        const gmx_hw_opt_t *hw_opt,
+                                        int nthreads_tot,
+                                        int ngpu)
+{
+    int nthreads_tmpi;
+
+    /* There are no separate PME nodes here, as we ensured in
+     * check_and_update_hw_opt that nthreads_tmpi>0 with PME nodes
+     * and a conditional ensures we would not have ended up here.
+     * Note that separate PME nodes might be switched on later.
+     */
+    if (ngpu > 0)
+    {
+        nthreads_tmpi = ngpu;
+        if (nthreads_tot > 0 && nthreads_tot < nthreads_tmpi)
+        {
+            nthreads_tmpi = nthreads_tot;
+        }
+    }
+    else if (hw_opt->nthreads_omp > 0)
+    {
+        /* Here we could oversubscribe, when we do, we issue a warning later */
+        nthreads_tmpi = max(1,nthreads_tot/hw_opt->nthreads_omp);
+    }
+    else
+    {
+        /* TODO choose nthreads_omp based on hardware topology
+           when we have a hardware topology detection library */
+        /* In general, when running up to 4 threads, OpenMP should be faster.
+         * Note: on AMD Bulldozer we should avoid running OpenMP over two dies.
+         * On Intel>=Nehalem running OpenMP on a single CPU is always faster,
+         * even on two CPUs it's usually faster (but with many OpenMP threads
+         * it could be faster not to use HT, currently we always use HT).
+         * On Nehalem/Westmere we want to avoid running 16 threads over
+         * two CPUs with HT, so we need a limit<16; thus we use 12.
+         * A reasonable limit for Intel Sandy and Ivy bridge,
+         * not knowing the topology, is 16 threads.
+         */
+        const int nthreads_omp_always_faster             =  4;
+        const int nthreads_omp_always_faster_Nehalem     = 12;
+        const int nthreads_omp_always_faster_SandyBridge = 16;
+        const int first_model_Nehalem     = 0x1A;
+        const int first_model_SandyBridge = 0x2A;
+        gmx_bool bIntel_Family6;
+
+        bIntel_Family6 =
+            (gmx_cpuid_vendor(hwinfo->cpuid_info) == GMX_CPUID_VENDOR_INTEL &&
+             gmx_cpuid_family(hwinfo->cpuid_info) == 6);
+
+        if (nthreads_tot <= nthreads_omp_always_faster ||
+            (bIntel_Family6 &&
+             ((gmx_cpuid_model(hwinfo->cpuid_info) >= nthreads_omp_always_faster_Nehalem && nthreads_tot <= nthreads_omp_always_faster_Nehalem) ||
+              (gmx_cpuid_model(hwinfo->cpuid_info) >= nthreads_omp_always_faster_SandyBridge && nthreads_tot <= nthreads_omp_always_faster_SandyBridge))))
+        {
+            /* Use pure OpenMP parallelization */
+            nthreads_tmpi = 1;
+        }
+        else
+        {
+            /* Don't use OpenMP parallelization */
+            nthreads_tmpi = nthreads_tot;
+        }
+    }
+
+    return nthreads_tmpi;
+}
+
+
+/* Get the number of threads to use for thread-MPI based on how many
+ * were requested, which algorithms we're using,
+ * and how many particles there are.
+ * At the point we have already called check_and_update_hw_opt.
+ * Thus all options should be internally consistent and consistent
+ * with the hardware, except that ntmpi could be larger than #GPU.
+ */
+static int get_nthreads_mpi(gmx_hw_info_t *hwinfo,
+                            gmx_hw_opt_t *hw_opt,
+                            t_inputrec *inputrec, gmx_mtop_t *mtop,
+                            const t_commrec *cr,
+                            FILE *fplog)
+{
+    int nthreads_hw,nthreads_tot_max,nthreads_tmpi,nthreads_new,ngpu;
+    int min_atoms_per_mpi_thread;
+    char *env;
+    char sbuf[STRLEN];
+    gmx_bool bCanUseGPU;
+
+    if (hw_opt->nthreads_tmpi > 0)
+    {
+        /* Trivial, return right away */
+        return hw_opt->nthreads_tmpi;
+    }
+
+    nthreads_hw = hwinfo->nthreads_hw_avail;
+
+    /* How many total (#tMPI*#OpenMP) threads can we start? */ 
+    if (hw_opt->nthreads_tot > 0)
+    {
+        nthreads_tot_max = hw_opt->nthreads_tot;
+    }
+    else
+    {
+        nthreads_tot_max = nthreads_hw;
+    }
+
+    bCanUseGPU = (inputrec->cutoff_scheme == ecutsVERLET && hwinfo->bCanUseGPU);
+    if (bCanUseGPU)
+    {
+        ngpu = hwinfo->gpu_info.ncuda_dev_use;
+    }
+    else
+    {
+        ngpu = 0;
+    }
+
+    nthreads_tmpi =
+        get_tmpi_omp_thread_division(hwinfo,hw_opt,nthreads_tot_max,ngpu);
+
+    if (inputrec->eI == eiNM || EI_TPI(inputrec->eI))
+    {
+        /* Steps are divided over the nodes iso splitting the atoms */
+        min_atoms_per_mpi_thread = 0;
+    }
+    else
+    {
+        if (bCanUseGPU)
+        {
+            min_atoms_per_mpi_thread = MIN_ATOMS_PER_GPU;
+        }
+        else
+        {
+            min_atoms_per_mpi_thread = MIN_ATOMS_PER_MPI_THREAD;
+        }
+    }
+
+    /* Check if an algorithm does not support parallel simulation.  */
+    if (nthreads_tmpi != 1 &&
+        ( inputrec->eI == eiLBFGS ||
+          inputrec->coulombtype == eelEWALD ) )
+    {
+        nthreads_tmpi = 1;
+
+        md_print_warn(cr,fplog,"The integration or electrostatics algorithm doesn't support parallel runs. Using a single thread-MPI thread.\n");
+        if (hw_opt->nthreads_tmpi > nthreads_tmpi)
+        {
+            gmx_fatal(FARGS,"You asked for more than 1 thread-MPI thread, but an algorithm doesn't support that");
+        }
+    }
+    else if (mtop->natoms/nthreads_tmpi < min_atoms_per_mpi_thread)
+    {
+        /* the thread number was chosen automatically, but there are too many
+           threads (too few atoms per thread) */
+        nthreads_new = max(1,mtop->natoms/min_atoms_per_mpi_thread);
+
+        /* Avoid partial use of Hyper-Threading */
+        if (gmx_cpuid_x86_smt(hwinfo->cpuid_info) == GMX_CPUID_X86_SMT_ENABLED &&
+            nthreads_new > nthreads_hw/2 && nthreads_new < nthreads_hw)
+        {
+            nthreads_new = nthreads_hw/2;
+        }
+
+        /* Avoid large prime numbers in the thread count */
+        if (nthreads_new >= 6)
+        {
+            /* Use only 6,8,10 with additional factors of 2 */
+            int fac;
+
+            fac = 2;
+            while (3*fac*2 <= nthreads_new)
+            {
+                fac *= 2;
+            }
+
+            nthreads_new = (nthreads_new/fac)*fac;
+        }
+        else
+        {
+            /* Avoid 5 */
+            if (nthreads_new == 5)
+            {
+                nthreads_new = 4;
+            }
+        }
+
+        nthreads_tmpi = nthreads_new;
+
+        fprintf(stderr,"\n");
+        fprintf(stderr,"NOTE: Parallelization is limited by the small number of atoms,\n");
+        fprintf(stderr,"      only starting %d thread-MPI threads.\n",nthreads_tmpi);
+        fprintf(stderr,"      You can use the -nt and/or -ntmpi option to optimize the number of threads.\n\n");
+    }
+
+    return nthreads_tmpi;
+}
+#endif /* GMX_THREAD_MPI */
+
+
+/* Environment variable for setting nstlist */
+static const char*  NSTLIST_ENVVAR          =  "GMX_NSTLIST";
+/* Try to increase nstlist when using a GPU with nstlist less than this */
+static const int    NSTLIST_GPU_ENOUGH      = 20;
+/* Increase nstlist until the non-bonded cost increases more than this factor */
+static const float  NBNXN_GPU_LIST_OK_FAC   = 1.25;
+/* Don't increase nstlist beyond a non-bonded cost increases of this factor */
+static const float  NBNXN_GPU_LIST_MAX_FAC  = 1.40;
+
+/* Try to increase nstlist when running on a GPU */
+static void increase_nstlist(FILE *fp,t_commrec *cr,
+                             t_inputrec *ir,const gmx_mtop_t *mtop,matrix box)
+{
+    char *env;
+    int  nstlist_orig,nstlist_prev;
+    verletbuf_list_setup_t ls;
+    real rlist_inc,rlist_ok,rlist_max,rlist_new,rlist_prev;
+    int  i;
+    t_state state_tmp;
+    gmx_bool bBox,bDD,bCont;
+    const char *nstl_fmt="\nFor optimal performance with a GPU nstlist (now %d) should be larger.\nThe optimum depends on your CPU and GPU resources.\nYou might want to try several nstlist values.\n";
+    const char *vbd_err="Can not increase nstlist for GPU run because verlet-buffer-drift is not set or used";
+    const char *box_err="Can not increase nstlist for GPU run because the box is too small";
+    const char *dd_err ="Can not increase nstlist for GPU run because of domain decomposition limitations";
+    char buf[STRLEN];
+
+    /* Number of + nstlist alternative values to try when switching  */
+    const int nstl[]={ 20, 25, 40, 50 };
+#define NNSTL  sizeof(nstl)/sizeof(nstl[0])
+
+    env = getenv(NSTLIST_ENVVAR);
+    if (env == NULL)
+    {
+        if (fp != NULL)
+        {
+            fprintf(fp,nstl_fmt,ir->nstlist);
+        }
+    }
+
+    if (ir->verletbuf_drift == 0)
+    {
+        gmx_fatal(FARGS,"You are using an old tpr file with a GPU, please generate a new tpr file with an up to date version of grompp");
+    }
+
+    if (ir->verletbuf_drift < 0)
+    {
+        if (MASTER(cr))
+        {
+            fprintf(stderr,"%s\n",vbd_err);
+        }
+        if (fp != NULL)
+        {
+            fprintf(fp,"%s\n",vbd_err);
+        }
+
+        return;
+    }
+
+    nstlist_orig = ir->nstlist;
+    if (env != NULL)
+    {
+        sprintf(buf,"Getting nstlist from environment variable GMX_NSTLIST=%s",env);
+        if (MASTER(cr))
+        {
+            fprintf(stderr,"%s\n",buf);
+        }
+        if (fp != NULL)
+        {
+            fprintf(fp,"%s\n",buf);
+        }
+        sscanf(env,"%d",&ir->nstlist);
+    }
+
+    verletbuf_get_list_setup(TRUE,&ls);
+
+    /* Allow rlist to make the list double the size of the cut-off sphere */
+    rlist_inc = nbnxn_get_rlist_effective_inc(NBNXN_GPU_CLUSTER_SIZE,mtop->natoms/det(box));
+    rlist_ok  = (max(ir->rvdw,ir->rcoulomb) + rlist_inc)*pow(NBNXN_GPU_LIST_OK_FAC,1.0/3.0) - rlist_inc;
+    rlist_max = (max(ir->rvdw,ir->rcoulomb) + rlist_inc)*pow(NBNXN_GPU_LIST_MAX_FAC,1.0/3.0) - rlist_inc;
+    if (debug)
+    {
+        fprintf(debug,"GPU nstlist tuning: rlist_inc %.3f rlist_max %.3f\n",
+                rlist_inc,rlist_max);
+    }
+
+    i = 0;
+    nstlist_prev = nstlist_orig;
+    rlist_prev   = ir->rlist;
+    do
+    {
+        if (env == NULL)
+        {
+            ir->nstlist = nstl[i];
+        }
+
+        /* Set the pair-list buffer size in ir */
+        calc_verlet_buffer_size(mtop,det(box),ir,ir->verletbuf_drift,&ls,
+                                NULL,&rlist_new);
+
+        /* Does rlist fit in the box? */
+        bBox = (sqr(rlist_new) < max_cutoff2(ir->ePBC,box));
+        bDD  = TRUE;
+        if (bBox && DOMAINDECOMP(cr))
+        {
+            /* Check if rlist fits in the domain decomposition */
+            if (inputrec2nboundeddim(ir) < DIM)
+            {
+                gmx_incons("Changing nstlist with domain decomposition and unbounded dimensions is not implemented yet");
+            }
+            copy_mat(box,state_tmp.box);
+            bDD = change_dd_cutoff(cr,&state_tmp,ir,rlist_new);
+        }
+
+        bCont = FALSE;
+
+        if (env == NULL)
+        {
+            if (bBox && bDD && rlist_new <= rlist_max)
+            {
+                /* Increase nstlist */
+                nstlist_prev = ir->nstlist;
+                rlist_prev   = rlist_new;
+                bCont = (i+1 < NNSTL && rlist_new < rlist_ok);
+            }
+            else
+            {
+                /* Stick with the previous nstlist */
+                ir->nstlist = nstlist_prev;
+                rlist_new   = rlist_prev;
+                bBox = TRUE;
+                bDD  = TRUE;
+            }
+        }
+
+        i++;
+    }
+    while (bCont);
+
+    if (!bBox || !bDD)
+    {
+        gmx_warning(!bBox ? box_err : dd_err);
+        if (fp != NULL)
+        {
+            fprintf(fp,"\n%s\n",bBox ? box_err : dd_err);
+        }
+        ir->nstlist = nstlist_orig;
+    }
+    else if (ir->nstlist != nstlist_orig || rlist_new != ir->rlist)
+    {
+        sprintf(buf,"Changing nstlist from %d to %d, rlist from %g to %g",
+                nstlist_orig,ir->nstlist,
+                ir->rlist,rlist_new);
+        if (MASTER(cr))
+        {
+            fprintf(stderr,"%s\n\n",buf);
+        }
+        if (fp != NULL)
+        {
+            fprintf(fp,"%s\n\n",buf);
+        }
+        ir->rlist     = rlist_new;
+        ir->rlistlong = rlist_new;
+    }
+}
+
+static void prepare_verlet_scheme(FILE *fplog,
+                                  gmx_hw_info_t *hwinfo,
+                                  t_commrec *cr,
+                                  gmx_hw_opt_t *hw_opt,
+                                  const char *nbpu_opt,
+                                  t_inputrec *ir,
+                                  const gmx_mtop_t *mtop,
+                                  matrix box,
+                                  gmx_bool *bUseGPU)
+{
+    /* Here we only check for GPU usage on the MPI master process,
+     * as here we don't know how many GPUs we will use yet.
+     * We check for a GPU on all processes later.
+     */
+    *bUseGPU = hwinfo->bCanUseGPU || (getenv("GMX_EMULATE_GPU") != NULL);
+
+    if (ir->verletbuf_drift > 0)
+    {
+        /* Update the Verlet buffer size for the current run setup */
+        verletbuf_list_setup_t ls;
+        real rlist_new;
+
+        /* Here we assume CPU acceleration is on. But as currently
+         * calc_verlet_buffer_size gives the same results for 4x8 and 4x4
+         * and 4x2 gives a larger buffer than 4x4, this is ok.
+         */
+        verletbuf_get_list_setup(*bUseGPU,&ls);
+
+        calc_verlet_buffer_size(mtop,det(box),ir,
+                                ir->verletbuf_drift,&ls,
+                                NULL,&rlist_new);
+        if (rlist_new != ir->rlist)
+        {
+            if (fplog != NULL)
+            {
+                fprintf(fplog,"\nChanging rlist from %g to %g for non-bonded %dx%d atom kernels\n\n",
+                        ir->rlist,rlist_new,
+                        ls.cluster_size_i,ls.cluster_size_j);
+            }
+            ir->rlist     = rlist_new;
+            ir->rlistlong = rlist_new;
+        }
+    }
+
+    /* With GPU or emulation we should check nstlist for performance */
+    if ((EI_DYNAMICS(ir->eI) &&
+         *bUseGPU &&
+         ir->nstlist < NSTLIST_GPU_ENOUGH) ||
+        getenv(NSTLIST_ENVVAR) != NULL)
+    {
+        /* Choose a better nstlist */
+        increase_nstlist(fplog,cr,ir,mtop,box);
+    }
+}
+
+static void convert_to_verlet_scheme(FILE *fplog,
+                                     t_inputrec *ir,
+                                     gmx_mtop_t *mtop,real box_vol)
+{
+    char *conv_mesg="Converting input file with group cut-off scheme to the Verlet cut-off scheme";
+
+    md_print_warn(NULL,fplog,"%s\n",conv_mesg);
+
+    ir->cutoff_scheme   = ecutsVERLET;
+    ir->verletbuf_drift = 0.005;
+
+    if (ir->rcoulomb != ir->rvdw)
+    {
+        gmx_fatal(FARGS,"The VdW and Coulomb cut-offs are different, whereas the Verlet scheme only supports equal cut-offs");
+    }
+
+    if (ir->vdwtype == evdwUSER || EEL_USER(ir->coulombtype))
+    {
+        gmx_fatal(FARGS,"User non-bonded potentials are not (yet) supported with the Verlet scheme");
+    }
+    else if (EVDW_SWITCHED(ir->vdwtype) || EEL_SWITCHED(ir->coulombtype))
+    {
+        md_print_warn(NULL,fplog,"Converting switched or shifted interactions to a shifted potential (without force shift), this will lead to slightly different interaction potentials");
+
+        if (EVDW_SWITCHED(ir->vdwtype))
+        {
+            ir->vdwtype = evdwCUT;
+        }
+        if (EEL_SWITCHED(ir->coulombtype))
+        {
+            if (EEL_FULL(ir->coulombtype))
+            {
+                /* With full electrostatic only PME can be switched */
+                ir->coulombtype = eelPME;
+            }
+            else
+            {
+                md_print_warn(NULL,fplog,"NOTE: Replacing %s electrostatics with reaction-field with epsilon-rf=inf\n",eel_names[ir->coulombtype]);
+                ir->coulombtype = eelRF;
+                ir->epsilon_rf  = 0.0;
+            }
+        }
+
+        /* We set the target energy drift to a small number.
+         * Note that this is only for testing. For production the user
+         * should think about this and set the mdp options.
+         */
+        ir->verletbuf_drift = 1e-4;
+    }
+
+    if (inputrec2nboundeddim(ir) != 3)
+    {
+        gmx_fatal(FARGS,"Can only convert old tpr files to the Verlet cut-off scheme with 3D pbc");
+    }
+
+    if (ir->efep != efepNO || ir->implicit_solvent != eisNO)
+    {
+        gmx_fatal(FARGS,"Will not convert old tpr files to the Verlet cut-off scheme with free-energy calculations or implicit solvent");
+    }
+
+    if (EI_DYNAMICS(ir->eI) && !(EI_MD(ir->eI) && ir->etc == etcNO))
+    {
+        verletbuf_list_setup_t ls;
+
+        verletbuf_get_list_setup(FALSE,&ls);
+        calc_verlet_buffer_size(mtop,box_vol,ir,ir->verletbuf_drift,&ls,
+                                NULL,&ir->rlist);
+    }
+    else
+    {
+        ir->verletbuf_drift = -1;
+        ir->rlist           = 1.05*max(ir->rvdw,ir->rcoulomb);
+    }
+
+    gmx_mtop_remove_chargegroups(mtop);
+}
+
+/* Check the process affinity mask. If it is non-zero, something
+ * else has set the affinity, and mdrun should honor that and
+ * not attempt to do its own thread pinning.
+ *
+ * This function should be called twice. Once before the OpenMP
+ * library gets initialized with bAfterOpenMPInit=FALSE (which will
+ * detect affinity set by external tools like taskset), and again
+ * later, after the OpenMP initialization, with bAfterOpenMPInit=TRUE
+ * (which will detect affinity changes made by the OpenMP library).
+ *
+ * Note that this will only work on Linux, because we use a GNU
+ * feature. */
+static void check_cpu_affinity_set(FILE *fplog, const t_commrec *cr,
+                                   gmx_hw_opt_t *hw_opt, int ncpus,
+                                   gmx_bool bAfterOpenmpInit)
+{
+#ifdef HAVE_SCHED_GETAFFINITY
+    cpu_set_t mask_current;
+    int       i, ret, cpu_count, cpu_set;
+    gmx_bool  bAllSet;
+
+    assert(hw_opt);
+    if (!hw_opt->bThreadPinning)
+    {
+        /* internal affinity setting is off, don't bother checking process affinity */
+        return;
+    }
+
+    CPU_ZERO(&mask_current);
+    if ((ret = sched_getaffinity(0, sizeof(cpu_set_t), &mask_current)) != 0)
+    {
+        /* failed to query affinity mask, will just return */
+        if (debug)
+        {
+            fprintf(debug, "Failed to query affinity mask (error %d)", ret);
+        }
+        return;
+    }
+
+    /* Before proceeding with the actual check, make sure that the number of
+     * detected CPUs is >= the CPUs in the current set.
+     * We need to check for CPU_COUNT as it was added only in glibc 2.6. */
+#ifdef CPU_COUNT
+    if (ncpus < CPU_COUNT(&mask_current))
+    {
+        if (debug)
+        {
+            fprintf(debug, "%d CPUs detected, but %d was returned by CPU_COUNT",
+                    ncpus, CPU_COUNT(&mask_current));
+        }
+        return;
+    }
+#endif /* CPU_COUNT */
+
+    bAllSet = TRUE;
+    for (i = 0; (i < ncpus && i < CPU_SETSIZE); i++)
+    {
+        bAllSet = bAllSet && (CPU_ISSET(i, &mask_current) != 0);
+    }
+
+    if (!bAllSet)
+    {
+        if (!bAfterOpenmpInit)
+        {
+            md_print_warn(cr, fplog,
+                          "%s detected a non-default process affinity, "
+                          "so it will not attempt to pin its threads", ShortProgram());
+        }
+        else
+        {
+            md_print_warn(cr, fplog,
+                          "%s detected a non-default process affinity, "
+                          "probably set by the OpenMP library, "
+                          "so it will not attempt to pin its threads", ShortProgram());
+        }
+        hw_opt->bThreadPinning = FALSE;
+
+        if (debug)
+        {
+            fprintf(debug, "Non-default affinity mask found, mdrun will not pin threads\n");
+        }
+    }
+    else
+    {
+        if (debug)
+        {
+            fprintf(debug, "Default affinity mask found\n");
+        }
+    }
+#endif /* HAVE_SCHED_GETAFFINITY */
+}
+
+/* Set CPU affinity. Can be important for performance.
+   On some systems (e.g. Cray) CPU Affinity is set by default.
+   But default assigning doesn't work (well) with only some ranks
+   having threads. This causes very low performance.
+   External tools have cumbersome syntax for setting affinity
+   in the case that only some ranks have threads.
+   Thus it is important that GROMACS sets the affinity internally
+   if only PME is using threads.
+*/
+static void set_cpu_affinity(FILE *fplog,
+                             const t_commrec *cr,
+                             gmx_hw_opt_t *hw_opt,
+                             int nthreads_pme,
+                             const gmx_hw_info_t *hwinfo,
+                             const t_inputrec *inputrec)
+{
+#if defined GMX_THREAD_MPI
+    /* With the number of TMPI threads equal to the number of cores
+     * we already pinned in thread-MPI, so don't pin again here.
+     */
+    if (hw_opt->nthreads_tmpi == tMPI_Thread_get_hw_number())
+    {
+        return;
+    }
+#endif
+
+#ifndef __APPLE__
+    /* If the tMPI thread affinity setting is not supported encourage the user
+     * to report it as it's either a bug or an exotic platform which we might
+     * want to support. */
+    if (tMPI_Thread_setaffinity_support() != TMPI_SETAFFINITY_SUPPORT_YES)
+    {
+        md_print_warn(NULL, fplog,
+                      "Can not set thread affinities on the current plarform. On NUMA systems this\n"
+                      "can cause performance degradation. If you think your platform should support\n"
+                      "setting affinities, contact the GROMACS developers.");
+        return;
+    }
+#endif /* __APPLE__ */
+
+    if (hw_opt->bThreadPinning)
+    {
+        int nth_affinity_set, thread_id_node, thread_id,
+            nthread_local, nthread_node, nthread_hw_max, nphyscore;
+        int offset;
+        char *env;
+
+        /* threads on this MPI process or TMPI thread */
+        if (cr->duty & DUTY_PP)
+        {
+            nthread_local = gmx_omp_nthreads_get(emntNonbonded);
+        }
+        else
+        {
+            nthread_local = gmx_omp_nthreads_get(emntPME);
+        }
+
+        /* map the current process to cores */
+        thread_id_node = 0;
+        nthread_node = nthread_local;
+#ifdef GMX_MPI
+        if (PAR(cr) || MULTISIM(cr))
+        {
+            /* We need to determine a scan of the thread counts in this
+             * compute node.
+             */
+            MPI_Comm comm_intra;
+
+            MPI_Comm_split(MPI_COMM_WORLD,gmx_hostname_num(),cr->rank_intranode,
+                           &comm_intra);
+            MPI_Scan(&nthread_local,&thread_id_node,1,MPI_INT,MPI_SUM,comm_intra);
+            /* MPI_Scan is inclusive, but here we need exclusive */
+            thread_id_node -= nthread_local;
+            /* Get the total number of threads on this physical node */
+            MPI_Allreduce(&nthread_local,&nthread_node,1,MPI_INT,MPI_SUM,comm_intra);
+            MPI_Comm_free(&comm_intra);
+        }
+#endif
+
+        offset = 0;
+        if (hw_opt->core_pinning_offset > 0)
+        {
+            offset = hw_opt->core_pinning_offset;
+            if (SIMMASTER(cr))
+            {
+                fprintf(stderr, "Applying core pinning offset %d\n", offset);
+            }
+            if (fplog)
+            {
+                fprintf(fplog, "Applying core pinning offset %d\n", offset);
+            }
+        }
+
+        /* With Intel Hyper-Threading enabled, we want to pin consecutive
+         * threads to physical cores when using more threads than physical
+         * cores or when the user requests so.
+         */
+        nthread_hw_max = hwinfo->nthreads_hw_avail;
+        nphyscore = -1;
+        if (hw_opt->bPinHyperthreading ||
+            (gmx_cpuid_x86_smt(hwinfo->cpuid_info) == GMX_CPUID_X86_SMT_ENABLED &&
+             nthread_node > nthread_hw_max/2 && getenv("GMX_DISABLE_PINHT") == NULL))
+        {
+            if (gmx_cpuid_x86_smt(hwinfo->cpuid_info) != GMX_CPUID_X86_SMT_ENABLED)
+            {
+                /* We print to stderr on all processes, as we might have
+                 * different settings on different physical nodes.
+                 */
+                if (gmx_cpuid_vendor(hwinfo->cpuid_info) != GMX_CPUID_VENDOR_INTEL)
+                {
+                    md_print_warn(NULL, fplog, "Pinning for Hyper-Threading layout requested, "
+                                  "but non-Intel CPU detected (vendor: %s)\n",
+                                  gmx_cpuid_vendor_string[gmx_cpuid_vendor(hwinfo->cpuid_info)]);
+                }
+                else
+                {
+                    md_print_warn(NULL, fplog, "Pinning for Hyper-Threading layout requested, "
+                                  "but the CPU detected does not have Intel Hyper-Threading support "
+                                  "(or it is turned off)\n");
+                }
+            }
+            nphyscore = nthread_hw_max/2;
+
+            if (SIMMASTER(cr))
+            {
+                fprintf(stderr, "Pinning to Hyper-Threading cores with %d physical cores in a compute node\n",
+                        nphyscore);
+            }
+            if (fplog)
+            {
+                fprintf(fplog, "Pinning to Hyper-Threading cores with %d physical cores in a compute node\n",
+                        nphyscore);
+            }
+        }
+
+        /* Set the per-thread affinity. In order to be able to check the success
+         * of affinity settings, we will set nth_affinity_set to 1 on threads
+         * where the affinity setting succeded and to 0 where it failed.
+         * Reducing these 0/1 values over the threads will give the total number
+         * of threads on which we succeeded.
+         */
+         nth_affinity_set = 0;
+#pragma omp parallel firstprivate(thread_id_node) num_threads(nthread_local) \
+                     reduction(+:nth_affinity_set)
+        {
+            int      core;
+            gmx_bool setaffinity_ret;
+
+            thread_id       = gmx_omp_get_thread_num();
+            thread_id_node += thread_id;
+            if (nphyscore <= 0)
+            {
+                core = offset + thread_id_node;
+            }
+            else
+            {
+                /* Lock pairs of threads to the same hyperthreaded core */
+                core = offset + thread_id_node/2 + (thread_id_node % 2)*nphyscore;
+            }
+
+            setaffinity_ret = tMPI_Thread_setaffinity_single(tMPI_Thread_self(), core);
+
+            /* store the per-thread success-values of the setaffinity */
+            nth_affinity_set = (setaffinity_ret == 0);
+
+            if (debug)
+            {
+                fprintf(debug, "On rank %2d, thread %2d, core %2d the affinity setting returned %d\n",
+                        cr->nodeid, gmx_omp_get_thread_num(), core, setaffinity_ret);
+            }
+        }
+
+        if (nth_affinity_set > nthread_local)
+        {
+            char msg[STRLEN];
+
+            sprintf(msg, "Looks like we have set affinity for more threads than "
+                    "we have (%d > %d)!\n", nth_affinity_set, nthread_local);
+            gmx_incons(msg);
+        }
+        else
+        {
+            /* check & warn if some threads failed to set their affinities */
+            if (nth_affinity_set != nthread_local)
+            {
+                char sbuf1[STRLEN], sbuf2[STRLEN];
+
+                /* sbuf1 contains rank info, while sbuf2 OpenMP thread info */
+                sbuf1[0] = sbuf2[0] = '\0';
+#ifdef GMX_MPI
+#ifdef GMX_THREAD_MPI
+                sprintf(sbuf1, "In thread-MPI thread #%d: ", cr->nodeid);
+#else /* GMX_LIB_MPI */
+                sprintf(sbuf1, "In MPI process #%d: ", cr->nodeid);
+#endif
+#endif /* GMX_MPI */
+
+                if (nthread_local > 1)
+                {
+                    sprintf(sbuf2, "of %d/%d thread%s ",
+                            nthread_local - nth_affinity_set, nthread_local,
+                            (nthread_local - nth_affinity_set) > 1 ? "s" : "");
+                }
+
+                md_print_warn(NULL, fplog,
+                              "NOTE: %sAffinity setting %sfailed.\n"
+                              "      This can cause performance degradation!",
+                              sbuf1, sbuf2);
+            }
+        }
+    }
+}
+
+
+static void check_and_update_hw_opt(gmx_hw_opt_t *hw_opt,
+                                    int cutoff_scheme)
+{
+    gmx_omp_nthreads_read_env(&hw_opt->nthreads_omp);
+
+#ifndef GMX_THREAD_MPI
+    if (hw_opt->nthreads_tot > 0)
+    {
+        gmx_fatal(FARGS,"Setting the total number of threads is only supported with thread-MPI and Gromacs was compiled without thread-MPI");
+    }
+    if (hw_opt->nthreads_tmpi > 0)
+    {
+        gmx_fatal(FARGS,"Setting the number of thread-MPI threads is only supported with thread-MPI and Gromacs was compiled without thread-MPI");
+    }
+#endif
+
+    if (hw_opt->nthreads_tot > 0 && hw_opt->nthreads_omp_pme <= 0)
+    {
+        /* We have the same number of OpenMP threads for PP and PME processes,
+         * thus we can perform several consistency checks.
+         */
+        if (hw_opt->nthreads_tmpi > 0 &&
+            hw_opt->nthreads_omp > 0 &&
+            hw_opt->nthreads_tot != hw_opt->nthreads_tmpi*hw_opt->nthreads_omp)
+        {
+            gmx_fatal(FARGS,"The total number of threads requested (%d) does not match the thread-MPI threads (%d) times the OpenMP threads (%d) requested",
+                      hw_opt->nthreads_tot,hw_opt->nthreads_tmpi,hw_opt->nthreads_omp);
+        }
+
+        if (hw_opt->nthreads_tmpi > 0 &&
+            hw_opt->nthreads_tot % hw_opt->nthreads_tmpi != 0)
+        {
+            gmx_fatal(FARGS,"The total number of threads requested (%d) is not divisible by the number of thread-MPI threads requested (%d)",
+                      hw_opt->nthreads_tot,hw_opt->nthreads_tmpi);
+        }
+
+        if (hw_opt->nthreads_omp > 0 &&
+            hw_opt->nthreads_tot % hw_opt->nthreads_omp != 0)
+        {
+            gmx_fatal(FARGS,"The total number of threads requested (%d) is not divisible by the number of OpenMP threads requested (%d)",
+                      hw_opt->nthreads_tot,hw_opt->nthreads_omp);
+        }
+
+        if (hw_opt->nthreads_tmpi > 0 &&
+            hw_opt->nthreads_omp <= 0)
+        {
+            hw_opt->nthreads_omp = hw_opt->nthreads_tot/hw_opt->nthreads_tmpi;
+        }
+    }
+
+#ifndef GMX_OPENMP
+    if (hw_opt->nthreads_omp > 1)
+    {
+        gmx_fatal(FARGS,"OpenMP threads are requested, but Gromacs was compiled without OpenMP support");
+    }
+#endif
+
+    if (cutoff_scheme == ecutsGROUP)
+    {
+        /* We only have OpenMP support for PME only nodes */
+        if (hw_opt->nthreads_omp > 1)
+        {
+            gmx_fatal(FARGS,"OpenMP threads have been requested with cut-off scheme %s, but these are only supported with cut-off scheme %s",
+                      ecutscheme_names[cutoff_scheme],
+                      ecutscheme_names[ecutsVERLET]);
+        }
+        hw_opt->nthreads_omp = 1;
+    }
+
+    if (hw_opt->nthreads_omp_pme > 0 && hw_opt->nthreads_omp <= 0)
+    {
+        gmx_fatal(FARGS,"You need to specify -ntomp in addition to -ntomp_pme");
+    }
+
+    if (hw_opt->nthreads_tot == 1)
+    {
+        hw_opt->nthreads_tmpi = 1;
+
+        if (hw_opt->nthreads_omp > 1)
+        {
+            gmx_fatal(FARGS,"You requested %d OpenMP threads with %d total threads",
+                      hw_opt->nthreads_tmpi,hw_opt->nthreads_tot);
+        }
+        hw_opt->nthreads_omp = 1;
+    }
+
+    if (hw_opt->nthreads_omp_pme <= 0 && hw_opt->nthreads_omp > 0)
+    {
+        hw_opt->nthreads_omp_pme = hw_opt->nthreads_omp;
+    }
+
+    if (debug)
+    {
+        fprintf(debug,"hw_opt: nt %d ntmpi %d ntomp %d ntomp_pme %d gpu_id '%s'\n",
+                hw_opt->nthreads_tot,
+                hw_opt->nthreads_tmpi,
+                hw_opt->nthreads_omp,
+                hw_opt->nthreads_omp_pme,
+                hw_opt->gpu_id!=NULL ? hw_opt->gpu_id : "");
+                
+    }
+}
+
+
+/* Override the value in inputrec with value passed on the command line (if any) */
+static void override_nsteps_cmdline(FILE *fplog,
+                                    int nsteps_cmdline,
+                                    t_inputrec *ir,
+                                    const t_commrec *cr)
+{
+    assert(ir);
+    assert(cr);
+
+    /* override with anything else than the default -2 */
+    if (nsteps_cmdline > -2)
+    {
+        char stmp[STRLEN];
+
+        ir->nsteps = nsteps_cmdline;
+        if (EI_DYNAMICS(ir->eI))
+        {
+            sprintf(stmp, "Overriding nsteps with value passed on the command line: %d steps, %.3f ps",
+                    nsteps_cmdline, nsteps_cmdline*ir->delta_t);
+        }
+        else
+        {
+            sprintf(stmp, "Overriding nsteps with value passed on the command line: %d steps",
+                    nsteps_cmdline);
+        }
+
+        md_print_warn(cr, fplog, "%s\n", stmp);
+    }
+}
+
+/* Data structure set by SIMMASTER which needs to be passed to all nodes
+ * before the other nodes have read the tpx file and called gmx_detect_hardware.
+ */
+typedef struct {
+    int cutoff_scheme; /* The cutoff scheme from inputrec_t */
+    gmx_bool bUseGPU;       /* Use GPU or GPU emulation          */
+} master_inf_t;
+
+int mdrunner(gmx_hw_opt_t *hw_opt,
+             FILE *fplog,t_commrec *cr,int nfile,
+             const t_filenm fnm[], const output_env_t oenv, gmx_bool bVerbose,
+             gmx_bool bCompact, int nstglobalcomm,
+             ivec ddxyz,int dd_node_order,real rdd,real rconstr,
+             const char *dddlb_opt,real dlb_scale,
+             const char *ddcsx,const char *ddcsy,const char *ddcsz,
+             const char *nbpu_opt,
+             int nsteps_cmdline, int nstepout,int resetstep,
+             int nmultisim,int repl_ex_nst,int repl_ex_nex,
+             int repl_ex_seed, real pforce,real cpt_period,real max_hours,
+             const char *deviceOptions, unsigned long Flags)
+{
+    gmx_bool   bForceUseGPU,bTryUseGPU;
+    double     nodetime=0,realtime;
+    t_inputrec *inputrec;
+    t_state    *state=NULL;
+    matrix     box;
+    gmx_ddbox_t ddbox={0};
+    int        npme_major,npme_minor;
+    real       tmpr1,tmpr2;
+    t_nrnb     *nrnb;
+    gmx_mtop_t *mtop=NULL;
+    t_mdatoms  *mdatoms=NULL;
+    t_forcerec *fr=NULL;
+    t_fcdata   *fcd=NULL;
+    real       ewaldcoeff=0;
+    gmx_pme_t  *pmedata=NULL;
+    gmx_vsite_t *vsite=NULL;
+    gmx_constr_t constr;
+    int        i,m,nChargePerturbed=-1,status,nalloc;
+    char       *gro;
+    gmx_wallcycle_t wcycle;
+    gmx_bool       bReadRNG,bReadEkin;
+    int        list;
+    gmx_runtime_t runtime;
+    int        rc;
+    gmx_large_int_t reset_counters;
+    gmx_edsam_t ed=NULL;
+    t_commrec   *cr_old=cr; 
+    int         nthreads_pme=1;
+    int         nthreads_pp=1;
+    gmx_membed_t membed=NULL;
+    gmx_hw_info_t *hwinfo=NULL;
+    master_inf_t minf={-1,FALSE};
+
+    /* CAUTION: threads may be started later on in this function, so
+       cr doesn't reflect the final parallel state right now */
+    snew(inputrec,1);
+    snew(mtop,1);
+    
+    if (Flags & MD_APPENDFILES) 
+    {
+        fplog = NULL;
+    }
+
+    bForceUseGPU = (strncmp(nbpu_opt, "gpu", 3) == 0);
+    bTryUseGPU   = (strncmp(nbpu_opt, "auto", 4) == 0) || bForceUseGPU;
+
+    snew(state,1);
+    if (SIMMASTER(cr)) 
+    {
+        /* Read (nearly) all data required for the simulation */
+        read_tpx_state(ftp2fn(efTPX,nfile,fnm),inputrec,state,NULL,mtop);
+
+        if (inputrec->cutoff_scheme != ecutsVERLET &&
+            ((Flags & MD_TESTVERLET) || getenv("GMX_VERLET_SCHEME") != NULL))
+        {
+            convert_to_verlet_scheme(fplog,inputrec,mtop,det(state->box));
+        }
+
+        /* Detect hardware, gather information. With tMPI only thread 0 does it
+         * and after threads are started broadcasts hwinfo around. */
+        snew(hwinfo, 1);
+        gmx_detect_hardware(fplog, hwinfo, cr,
+                            bForceUseGPU, bTryUseGPU, hw_opt->gpu_id);
+
+        minf.cutoff_scheme = inputrec->cutoff_scheme;
+        minf.bUseGPU       = FALSE;
+
+        if (inputrec->cutoff_scheme == ecutsVERLET)
+        {
+            prepare_verlet_scheme(fplog,hwinfo,cr,hw_opt,nbpu_opt,
+                                  inputrec,mtop,state->box,
+                                  &minf.bUseGPU);
+        }
+        else if (hwinfo->bCanUseGPU)
+        {
+            md_print_warn(cr,fplog,
+                          "NOTE: GPU(s) found, but the current simulation can not use GPUs\n"
+                          "      To use a GPU, set the mdp option: cutoff-scheme = Verlet\n"
+                          "      (for quick performance testing you can use the -testverlet option)\n");
+
+            if (bForceUseGPU)
+            {
+                gmx_fatal(FARGS,"GPU requested, but can't be used without cutoff-scheme=Verlet");
+            }
+        }
+    }
+#ifndef GMX_THREAD_MPI
+    if (PAR(cr))
+    {
+        gmx_bcast_sim(sizeof(minf),&minf,cr);
+    }
+#endif
+    if (minf.bUseGPU && cr->npmenodes == -1)
+    {
+        /* Don't automatically use PME-only nodes with GPUs */
+        cr->npmenodes = 0;
+    }
+
+    /* Check for externally set OpenMP affinity and turn off internal
+     * pinning if any is found. We need to do this check early to tell
+     * thread-MPI whether it should do pinning when spawning threads.
+     */
+    gmx_omp_check_thread_affinity(fplog, cr, hw_opt);
+
+#ifdef GMX_THREAD_MPI
+    /* With thread-MPI inputrec is only set here on the master thread */
+    if (SIMMASTER(cr))
+#endif
+    {
+        check_and_update_hw_opt(hw_opt,minf.cutoff_scheme);
+
+#ifdef GMX_THREAD_MPI
+        /* Early check for externally set process affinity. Can't do over all
+         * MPI processes because hwinfo is not available everywhere, but with
+         * thread-MPI it's needed as pinning might get turned off which needs
+         * to be known before starting thread-MPI. */
+        check_cpu_affinity_set(fplog,
+                               NULL,
+                               hw_opt, hwinfo->nthreads_hw_avail, FALSE);
+#endif
+
+#ifdef GMX_THREAD_MPI
+        if (cr->npmenodes > 0 && hw_opt->nthreads_tmpi <= 0)
+        {
+            gmx_fatal(FARGS,"You need to explicitly specify the number of MPI threads (-ntmpi) when using separate PME nodes");
+        }
+#endif
+
+        if (hw_opt->nthreads_omp_pme != hw_opt->nthreads_omp &&
+            cr->npmenodes <= 0)
+        {
+            gmx_fatal(FARGS,"You need to explicitly specify the number of PME nodes (-npme) when using different number of OpenMP threads for PP and PME nodes");
+        }
+    }
+
+#ifdef GMX_THREAD_MPI
+    if (SIMMASTER(cr))
+    {
+        /* NOW the threads will be started: */
+        hw_opt->nthreads_tmpi = get_nthreads_mpi(hwinfo,
+                                                 hw_opt,
+                                                 inputrec, mtop,
+                                                 cr, fplog);
+        if (hw_opt->nthreads_tot > 0 && hw_opt->nthreads_omp <= 0)
+        {
+            hw_opt->nthreads_omp = hw_opt->nthreads_tot/hw_opt->nthreads_tmpi;
+        }
+
+        if (hw_opt->nthreads_tmpi > 1)
+        {
+            /* now start the threads. */
+            cr=mdrunner_start_threads(hw_opt, fplog, cr_old, nfile, fnm, 
+                                      oenv, bVerbose, bCompact, nstglobalcomm, 
+                                      ddxyz, dd_node_order, rdd, rconstr, 
+                                      dddlb_opt, dlb_scale, ddcsx, ddcsy, ddcsz,
+                                      nbpu_opt,
+                                      nsteps_cmdline, nstepout, resetstep, nmultisim, 
+                                      repl_ex_nst, repl_ex_nex, repl_ex_seed, pforce,
+                                      cpt_period, max_hours, deviceOptions, 
+                                      Flags);
+            /* the main thread continues here with a new cr. We don't deallocate
+               the old cr because other threads may still be reading it. */
+            if (cr == NULL)
+            {
+                gmx_comm("Failed to spawn threads");
+            }
+        }
+    }
+#endif
+    /* END OF CAUTION: cr is now reliable */
+
+    /* g_membed initialisation *
+     * Because we change the mtop, init_membed is called before the init_parallel *
+     * (in case we ever want to make it run in parallel) */
+    if (opt2bSet("-membed",nfile,fnm))
+    {
+        if (MASTER(cr))
+        {
+            fprintf(stderr,"Initializing membed");
+        }
+        membed = init_membed(fplog,nfile,fnm,mtop,inputrec,state,cr,&cpt_period);
+    }
+
+    if (PAR(cr))
+    {
+        /* now broadcast everything to the non-master nodes/threads: */
+        init_parallel(fplog, cr, inputrec, mtop);
+
+        /* This check needs to happen after get_nthreads_mpi() */
+        if (inputrec->cutoff_scheme == ecutsVERLET && (Flags & MD_PARTDEC))
+        {
+            gmx_fatal_collective(FARGS,cr,NULL,
+                                 "The Verlet cut-off scheme is not supported with particle decomposition.\n"
+                                 "You can achieve the same effect as particle decomposition by running in parallel using only OpenMP threads.");
+        }
+    }
+    if (fplog != NULL)
+    {
+        pr_inputrec(fplog,0,"Input Parameters",inputrec,FALSE);
+    }
+
+#if defined GMX_THREAD_MPI
+    /* With tMPI we detected on thread 0 and we'll just pass the hwinfo pointer
+     * to the other threads  -- slightly uncool, but works fine, just need to
+     * make sure that the data doesn't get freed twice. */
+    if (cr->nnodes > 1)
+    {
+        if (!SIMMASTER(cr))
+        {
+            snew(hwinfo, 1);
+        }
+        gmx_bcast(sizeof(&hwinfo), &hwinfo, cr);
+    }
+#else
+    if (PAR(cr) && !SIMMASTER(cr))
+    {
+        /* now we have inputrec on all nodes, can run the detection */
+        /* TODO: perhaps it's better to propagate within a node instead? */
+        snew(hwinfo, 1);
+        gmx_detect_hardware(fplog, hwinfo, cr,
+                                 bForceUseGPU, bTryUseGPU, hw_opt->gpu_id);
+    }
+
+    /* Now do the affinity check with MPI/no-MPI (done earlier with thread-MPI). */
+    check_cpu_affinity_set(fplog, cr,
+                           hw_opt, hwinfo->nthreads_hw_avail, FALSE);
+#endif
+
+    /* now make sure the state is initialized and propagated */
+    set_state_entries(state,inputrec,cr->nnodes);
+
+    /* remove when vv and rerun works correctly! */
+    if (PAR(cr) && EI_VV(inputrec->eI) && ((Flags & MD_RERUN) || (Flags & MD_RERUN_VSITE)))
+    {
+        gmx_fatal(FARGS,
+                  "Currently can't do velocity verlet with rerun in parallel.");
+    }
+
+    /* A parallel command line option consistency check that we can
+       only do after any threads have started. */
+    if (!PAR(cr) &&
+        (ddxyz[XX] > 1 || ddxyz[YY] > 1 || ddxyz[ZZ] > 1 || cr->npmenodes > 0))
+    {
+        gmx_fatal(FARGS,
+                  "The -dd or -npme option request a parallel simulation, "
+#ifndef GMX_MPI
+                  "but %s was compiled without threads or MPI enabled"
+#else
+#ifdef GMX_THREAD_MPI
+                  "but the number of threads (option -nt) is 1"
+#else
+                  "but %s was not started through mpirun/mpiexec or only one process was requested through mpirun/mpiexec"
+#endif
+#endif
+                  , ShortProgram()
+            );
+    }
+
+    if ((Flags & MD_RERUN) &&
+        (EI_ENERGY_MINIMIZATION(inputrec->eI) || eiNM == inputrec->eI))
+    {
+        gmx_fatal(FARGS, "The .mdp file specified an energy mininization or normal mode algorithm, and these are not compatible with mdrun -rerun");
+    }
+
+    if (can_use_allvsall(inputrec,mtop,TRUE,cr,fplog) && PAR(cr))
+    {
+        /* All-vs-all loops do not work with domain decomposition */
+        Flags |= MD_PARTDEC;
+    }
+
+    if (!EEL_PME(inputrec->coulombtype) || (Flags & MD_PARTDEC))
+    {
+        if (cr->npmenodes > 0)
+        {
+            if (!EEL_PME(inputrec->coulombtype))
+            {
+                gmx_fatal_collective(FARGS,cr,NULL,
+                                     "PME nodes are requested, but the system does not use PME electrostatics");
+            }
+            if (Flags & MD_PARTDEC)
+            {
+                gmx_fatal_collective(FARGS,cr,NULL,
+                                     "PME nodes are requested, but particle decomposition does not support separate PME nodes");
+            }
+        }
+
+        cr->npmenodes = 0;
+    }
+
+#ifdef GMX_FAHCORE
+    fcRegisterSteps(inputrec->nsteps,inputrec->init_step);
+#endif
+
+    /* NMR restraints must be initialized before load_checkpoint,
+     * since with time averaging the history is added to t_state.
+     * For proper consistency check we therefore need to extend
+     * t_state here.
+     * So the PME-only nodes (if present) will also initialize
+     * the distance restraints.
+     */
+    snew(fcd,1);
+
+    /* This needs to be called before read_checkpoint to extend the state */
+    init_disres(fplog,mtop,inputrec,cr,Flags & MD_PARTDEC,fcd,state);
+
+    if (gmx_mtop_ftype_count(mtop,F_ORIRES) > 0)
+    {
+        if (PAR(cr) && !(Flags & MD_PARTDEC))
+        {
+            gmx_fatal(FARGS,"Orientation restraints do not work (yet) with domain decomposition, use particle decomposition (mdrun option -pd)");
+        }
+        /* Orientation restraints */
+        if (MASTER(cr))
+        {
+            init_orires(fplog,mtop,state->x,inputrec,cr->ms,&(fcd->orires),
+                        state);
+        }
+    }
+
+    if (DEFORM(*inputrec))
+    {
+        /* Store the deform reference box before reading the checkpoint */
+        if (SIMMASTER(cr))
+        {
+            copy_mat(state->box,box);
+        }
+        if (PAR(cr))
+        {
+            gmx_bcast(sizeof(box),box,cr);
+        }
+        /* Because we do not have the update struct available yet
+         * in which the reference values should be stored,
+         * we store them temporarily in static variables.
+         * This should be thread safe, since they are only written once
+         * and with identical values.
+         */
+#ifdef GMX_THREAD_MPI
+        tMPI_Thread_mutex_lock(&deform_init_box_mutex);
+#endif
+        deform_init_init_step_tpx = inputrec->init_step;
+        copy_mat(box,deform_init_box_tpx);
+#ifdef GMX_THREAD_MPI
+        tMPI_Thread_mutex_unlock(&deform_init_box_mutex);
+#endif
+    }
+
+    if (opt2bSet("-cpi",nfile,fnm)) 
+    {
+        /* Check if checkpoint file exists before doing continuation.
+         * This way we can use identical input options for the first and subsequent runs...
+         */
+        if( gmx_fexist_master(opt2fn_master("-cpi",nfile,fnm,cr),cr) )
+        {
+            load_checkpoint(opt2fn_master("-cpi",nfile,fnm,cr),&fplog,
+                            cr,Flags & MD_PARTDEC,ddxyz,
+                            inputrec,state,&bReadRNG,&bReadEkin,
+                            (Flags & MD_APPENDFILES),
+                            (Flags & MD_APPENDFILESSET));
+            
+            if (bReadRNG)
+            {
+                Flags |= MD_READ_RNG;
+            }
+            if (bReadEkin)
+            {
+                Flags |= MD_READ_EKIN;
+            }
+        }
+    }
+
+    if (((MASTER(cr) || (Flags & MD_SEPPOT)) && (Flags & MD_APPENDFILES))
+#ifdef GMX_THREAD_MPI
+        /* With thread MPI only the master node/thread exists in mdrun.c,
+         * therefore non-master nodes need to open the "seppot" log file here.
+         */
+        || (!MASTER(cr) && (Flags & MD_SEPPOT))
+#endif
+        )
+    {
+        gmx_log_open(ftp2fn(efLOG,nfile,fnm),cr,!(Flags & MD_SEPPOT),
+                             Flags,&fplog);
+    }
+
+    /* override nsteps with value from cmdline */
+    override_nsteps_cmdline(fplog, nsteps_cmdline, inputrec, cr);
+
+    if (SIMMASTER(cr)) 
+    {
+        copy_mat(state->box,box);
+    }
+
+    if (PAR(cr)) 
+    {
+        gmx_bcast(sizeof(box),box,cr);
+    }
+
+    /* Essential dynamics */
+    if (opt2bSet("-ei",nfile,fnm))
+    {
+        /* Open input and output files, allocate space for ED data structure */
+        ed = ed_open(nfile,fnm,Flags,cr);
+    }
+
+    if (PAR(cr) && !((Flags & MD_PARTDEC) ||
+                     EI_TPI(inputrec->eI) ||
+                     inputrec->eI == eiNM))
+    {
+        cr->dd = init_domain_decomposition(fplog,cr,Flags,ddxyz,rdd,rconstr,
+                                           dddlb_opt,dlb_scale,
+                                           ddcsx,ddcsy,ddcsz,
+                                           mtop,inputrec,
+                                           box,state->x,
+                                           &ddbox,&npme_major,&npme_minor);
+
+        make_dd_communicators(fplog,cr,dd_node_order);
+
+        /* Set overallocation to avoid frequent reallocation of arrays */
+        set_over_alloc_dd(TRUE);
+    }
+    else
+    {
+        /* PME, if used, is done on all nodes with 1D decomposition */
+        cr->npmenodes = 0;
+        cr->duty = (DUTY_PP | DUTY_PME);
+        npme_major = 1;
+        npme_minor = 1;
+        if (!EI_TPI(inputrec->eI))
+        {
+            npme_major = cr->nnodes;
+        }
+        
+        if (inputrec->ePBC == epbcSCREW)
+        {
+            gmx_fatal(FARGS,
+                      "pbc=%s is only implemented with domain decomposition",
+                      epbc_names[inputrec->ePBC]);
+        }
+    }
+
+    if (PAR(cr))
+    {
+        /* After possible communicator splitting in make_dd_communicators.
+         * we can set up the intra/inter node communication.
+         */
+        gmx_setup_nodecomm(fplog,cr);
+    }
+
+    /* Initialize per-physical-node MPI process/thread ID and counters. */
+    gmx_init_intranode_counters(cr);
+
+#ifdef GMX_MPI
+    md_print_info(cr,fplog,"Using %d MPI %s\n",
+                  cr->nnodes,
+#ifdef GMX_THREAD_MPI
+                  cr->nnodes==1 ? "thread" : "threads"
+#else
+                  cr->nnodes==1 ? "process" : "processes"
+#endif
+                  );
+    fflush(stderr);
+#endif
+
+    gmx_omp_nthreads_init(fplog, cr,
+                          hwinfo->nthreads_hw_avail,
+                          hw_opt->nthreads_omp,
+                          hw_opt->nthreads_omp_pme,
+                          (cr->duty & DUTY_PP) == 0,
+                          inputrec->cutoff_scheme == ecutsVERLET);
+
+    gmx_check_hw_runconf_consistency(fplog, hwinfo, cr, hw_opt->nthreads_tmpi, minf.bUseGPU);
+
+    /* getting number of PP/PME threads
+       PME: env variable should be read only on one node to make sure it is 
+       identical everywhere;
+     */
+    /* TODO nthreads_pp is only used for pinning threads.
+     * This is a temporary solution until we have a hw topology library.
+     */
+    nthreads_pp  = gmx_omp_nthreads_get(emntNonbonded);
+    nthreads_pme = gmx_omp_nthreads_get(emntPME);
+
+    wcycle = wallcycle_init(fplog,resetstep,cr,nthreads_pp,nthreads_pme);
+
+    if (PAR(cr))
+    {
+        /* Master synchronizes its value of reset_counters with all nodes 
+         * including PME only nodes */
+        reset_counters = wcycle_get_reset_counters(wcycle);
+        gmx_bcast_sim(sizeof(reset_counters),&reset_counters,cr);
+        wcycle_set_reset_counters(wcycle, reset_counters);
+    }
+
+    snew(nrnb,1);
+    if (cr->duty & DUTY_PP)
+    {
+        /* For domain decomposition we allocate dynamically
+         * in dd_partition_system.
+         */
+        if (DOMAINDECOMP(cr))
+        {
+            bcast_state_setup(cr,state);
+        }
+        else
+        {
+            if (PAR(cr))
+            {
+                bcast_state(cr,state,TRUE);
+            }
+        }
+
+        /* Initiate forcerecord */
+        fr = mk_forcerec();
+        fr->hwinfo = hwinfo;
+        init_forcerec(fplog,oenv,fr,fcd,inputrec,mtop,cr,box,FALSE,
+                      opt2fn("-table",nfile,fnm),
+                      opt2fn("-tabletf",nfile,fnm),
+                      opt2fn("-tablep",nfile,fnm),
+                      opt2fn("-tableb",nfile,fnm),
+                      nbpu_opt,
+                      FALSE,pforce);
+
+        /* version for PCA_NOT_READ_NODE (see md.c) */
+        /*init_forcerec(fplog,fr,fcd,inputrec,mtop,cr,box,FALSE,
+          "nofile","nofile","nofile","nofile",FALSE,pforce);
+          */        
+        fr->bSepDVDL = ((Flags & MD_SEPPOT) == MD_SEPPOT);
+
+        /* Initialize QM-MM */
+        if(fr->bQMMM)
+        {
+            init_QMMMrec(cr,box,mtop,inputrec,fr);
+        }
+
+        /* Initialize the mdatoms structure.
+         * mdatoms is not filled with atom data,
+         * as this can not be done now with domain decomposition.
+         */
+        mdatoms = init_mdatoms(fplog,mtop,inputrec->efep!=efepNO);
+
+        /* Initialize the virtual site communication */
+        vsite = init_vsite(mtop,cr,FALSE);
+
+        calc_shifts(box,fr->shift_vec);
+
+        /* With periodic molecules the charge groups should be whole at start up
+         * and the virtual sites should not be far from their proper positions.
+         */
+        if (!inputrec->bContinuation && MASTER(cr) &&
+            !(inputrec->ePBC != epbcNONE && inputrec->bPeriodicMols))
+        {
+            /* Make molecules whole at start of run */
+            if (fr->ePBC != epbcNONE)
+            {
+                do_pbc_first_mtop(fplog,inputrec->ePBC,box,mtop,state->x);
+            }
+            if (vsite)
+            {
+                /* Correct initial vsite positions are required
+                 * for the initial distribution in the domain decomposition
+                 * and for the initial shell prediction.
+                 */
+                construct_vsites_mtop(fplog,vsite,mtop,state->x);
+            }
+        }
+
+        if (EEL_PME(fr->eeltype))
+        {
+            ewaldcoeff = fr->ewaldcoeff;
+            pmedata = &fr->pmedata;
+        }
+        else
+        {
+            pmedata = NULL;
+        }
+    }
+    else
+    {
+        /* This is a PME only node */
+
+        /* We don't need the state */
+        done_state(state);
+
+        ewaldcoeff = calc_ewaldcoeff(inputrec->rcoulomb, inputrec->ewald_rtol);
+        snew(pmedata,1);
+    }
+
+    /* Before setting affinity, check whether the affinity has changed
+     * - which indicates that probably the OpenMP library has changed it since
+     * we first checked). */
+    check_cpu_affinity_set(fplog, cr, hw_opt, hwinfo->nthreads_hw_avail, TRUE);
+
+    /* Set the CPU affinity */
+    set_cpu_affinity(fplog,cr,hw_opt,nthreads_pme,hwinfo,inputrec);
+
+    /* Initiate PME if necessary,
+     * either on all nodes or on dedicated PME nodes only. */
+    if (EEL_PME(inputrec->coulombtype))
+    {
+        if (mdatoms)
+        {
+            nChargePerturbed = mdatoms->nChargePerturbed;
+        }
+        if (cr->npmenodes > 0)
+        {
+            /* The PME only nodes need to know nChargePerturbed */
+            gmx_bcast_sim(sizeof(nChargePerturbed),&nChargePerturbed,cr);
+        }
+
+        if (cr->duty & DUTY_PME)
+        {
+            status = gmx_pme_init(pmedata,cr,npme_major,npme_minor,inputrec,
+                                  mtop ? mtop->natoms : 0,nChargePerturbed,
+                                  (Flags & MD_REPRODUCIBLE),nthreads_pme);
+            if (status != 0) 
+            {
+                gmx_fatal(FARGS,"Error %d initializing PME",status);
+            }
+        }
+    }
+
+
+    if (integrator[inputrec->eI].func == do_md
+        ||
+        integrator[inputrec->eI].func == do_md_openmm
+        )
+    {
+        /* Turn on signal handling on all nodes */
+        /*
+         * (A user signal from the PME nodes (if any)
+         * is communicated to the PP nodes.
+         */
+        signal_handler_install();
+    }
+
+    if (cr->duty & DUTY_PP)
+    {
+        if (inputrec->ePull != epullNO)
+        {
+            /* Initialize pull code */
+            init_pull(fplog,inputrec,nfile,fnm,mtop,cr,oenv, inputrec->fepvals->init_lambda,
+                      EI_DYNAMICS(inputrec->eI) && MASTER(cr),Flags);
+        }
+        
+        if (inputrec->bRot)
+        {
+           /* Initialize enforced rotation code */
+           init_rot(fplog,inputrec,nfile,fnm,cr,state->x,box,mtop,oenv,
+                    bVerbose,Flags);
+        }
+
+        constr = init_constraints(fplog,mtop,inputrec,ed,state,cr);
+
+        if (DOMAINDECOMP(cr))
+        {
+            dd_init_bondeds(fplog,cr->dd,mtop,vsite,constr,inputrec,
+                            Flags & MD_DDBONDCHECK,fr->cginfo_mb);
+
+            set_dd_parameters(fplog,cr->dd,dlb_scale,inputrec,fr,&ddbox);
+
+            setup_dd_grid(fplog,cr->dd);
+        }
+
+        /* Now do whatever the user wants us to do (how flexible...) */
+        integrator[inputrec->eI].func(fplog,cr,nfile,fnm,
+                                      oenv,bVerbose,bCompact,
+                                      nstglobalcomm,
+                                      vsite,constr,
+                                      nstepout,inputrec,mtop,
+                                      fcd,state,
+                                      mdatoms,nrnb,wcycle,ed,fr,
+                                      repl_ex_nst,repl_ex_nex,repl_ex_seed,
+                                      membed,
+                                      cpt_period,max_hours,
+                                      deviceOptions,
+                                      Flags,
+                                      &runtime);
+
+        if (inputrec->ePull != epullNO)
+        {
+            finish_pull(fplog,inputrec->pull);
+        }
+        
+        if (inputrec->bRot)
+        {
+            finish_rot(fplog,inputrec->rot);
+        }
+
+    } 
+    else 
+    {
+        /* do PME only */
+        gmx_pmeonly(*pmedata,cr,nrnb,wcycle,ewaldcoeff,FALSE,inputrec);
+    }
+
+    if (EI_DYNAMICS(inputrec->eI) || EI_TPI(inputrec->eI))
+    {
+        /* Some timing stats */  
+        if (SIMMASTER(cr))
+        {
+            if (runtime.proc == 0)
+            {
+                runtime.proc = runtime.real;
+            }
+        }
+        else
+        {
+            runtime.real = 0;
+        }
+    }
+
+    wallcycle_stop(wcycle,ewcRUN);
+
+    /* Finish up, write some stuff
+     * if rerunMD, don't write last frame again 
+     */
+    finish_run(fplog,cr,ftp2fn(efSTO,nfile,fnm),
+               inputrec,nrnb,wcycle,&runtime,
+               fr != NULL && fr->nbv != NULL && fr->nbv->bUseGPU ?
+                 nbnxn_cuda_get_timings(fr->nbv->cu_nbv) : NULL,
+               nthreads_pp, 
+               EI_DYNAMICS(inputrec->eI) && !MULTISIM(cr));
+
+    if ((cr->duty & DUTY_PP) && fr->nbv != NULL && fr->nbv->bUseGPU)
+    {
+        char gpu_err_str[STRLEN];
+
+        /* free GPU memory and uninitialize GPU (by destroying the context) */
+        nbnxn_cuda_free(fplog, fr->nbv->cu_nbv);
+
+        if (!free_gpu(gpu_err_str))
+        {
+            gmx_warning("On node %d failed to free GPU #%d: %s",
+                        cr->nodeid, get_current_gpu_device_id(), gpu_err_str);
+        }
+    }
+
+    if (opt2bSet("-membed",nfile,fnm))
+    {
+        sfree(membed);
+    }
+
+#ifdef GMX_THREAD_MPI
+    if (PAR(cr) && SIMMASTER(cr))
+#endif
+    {
+        gmx_hardware_info_free(hwinfo);
+    }
+
+    /* Does what it says */  
+    print_date_and_time(fplog,cr->nodeid,"Finished mdrun",&runtime);
+
+    /* Close logfile already here if we were appending to it */
+    if (MASTER(cr) && (Flags & MD_APPENDFILES))
+    {
+        gmx_log_close(fplog);
+    }  
+
+    rc=(int)gmx_get_stop_condition();
+
+#ifdef GMX_THREAD_MPI
+    /* we need to join all threads. The sub-threads join when they
+       exit this function, but the master thread needs to be told to 
+       wait for that. */
+    if (PAR(cr) && MASTER(cr))
+    {
+        tMPI_Finalize();
+    }
+#endif
+
+    return rc;
+}
index bd7b9f93e5ea518b467c0181344bd1e5a92cfd85..e11367dc7e8451533c6a1aaf878ec9102c24b932 100644 (file)
@@ -2,7 +2,7 @@
  * This file is part of the GROMACS molecular simulation package.
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2010, The GROMACS development team,
+ * Copyright (c) 2001-2010,2012 The GROMACS development team,
  * check out http://www.gromacs.org for more information.
  * Copyright (c) 2012, by the GROMACS development team, led by
  * David van der Spoel, Berk Hess, Erik Lindahl, and including many
@@ -60,9 +60,6 @@
 #define TIMED_TESTS     MOD_20_32BIT | LOGIC_4_ITER_SHMEM | RANDOM_BLOCKS /*!< Bit flag with type of tests to
                                                                             run in time constrained memtest. */
 
-/*! Number of supported GPUs */
-#define NB_GPUS (sizeof(SupportedGPUs)/sizeof(SupportedGPUs[0]))
-
 static int cuda_max_device_count = 32; /*! Max number of devices supported by CUDA (for consistency checking).
                                            In reality it 16 with CUDA <=v5.0, but let's stay on the safe side. */
 
@@ -87,55 +84,6 @@ enum memtest_G80_test_types {
     LOGIC_4_ITER_SHMEM =        0x1000
 };
 
-// TODO put this list into an external file and include it so that the list is easily accessible
-/*! List of supported GPUs. */
-static const char * const SupportedGPUs[] = {
-    /* GT400 */
-    "Geforce GTX 480",
-    "Geforce GTX 470",
-    "Geforce GTX 465",
-    "Geforce GTX 460",
-
-    "Tesla C2070",
-    "Tesla C2050",
-    "Tesla S2070",
-    "Tesla S2050",
-    "Tesla M2070",
-    "Tesla M2050",
-
-    "Quadro 5000",
-    "Quadro 6000",
-
-    /* GT200 */
-    "Geforce GTX 295",
-    "Geforce GTX 285",
-    "Geforce GTX 280",
-    "Geforce GTX 275",
-    "Geforce GTX 260",
-    "GeForce GTS 250",
-    "GeForce GTS 150",
-
-    "GeForce GTX 285M",
-    "GeForce GTX 280M",
-
-    "Tesla S1070",
-    "Tesla C1060",
-    "Tesla M1060",
-
-    "Quadro FX 5800",
-    "Quadro FX 4800",
-    "Quadro CX",
-    "Quadro Plex 2200 D2",
-    "Quadro Plex 2200 S4",
-
-    /* G90 */
-    "GeForce 9800 G", /* GX2, GTX, GTX+, GT */
-    "GeForce 9800M GTX",
-
-    "Quadro FX 4700",
-    "Quadro Plex 2100 D4"
-};
-
 
 /*! 
   * \brief Runs GPU sanity checks.
@@ -243,48 +191,6 @@ static int do_sanity_checks(int dev_id, cudaDeviceProp *dev_prop)
 }
 
 
-/*! 
- * \brief Checks whether the GPU with the given name is supported in Gromacs-OpenMM.
- * 
- * \param[in] gpu_name  the name of the CUDA device
- * \returns             TRUE if the device is supported, otherwise FALSE
- */
-static bool is_gmx_openmm_supported_gpu_name(char *gpuName)
-{
-    size_t i;
-    for (i = 0; i < NB_GPUS; i++)
-    {
-        trim(gpuName);
-        if (gmx_strncasecmp(gpuName, SupportedGPUs[i], strlen(SupportedGPUs[i])) == 0)
-            return 1;
-    }
-    return 0;
-}
-
-/*! \brief Checks whether the GPU with the given device id is supported in Gromacs-OpenMM.
- *
- * \param[in] dev_id    the device id of the GPU or -1 if the device has already been selected
- * \param[out] gpu_name Set to contain the name of the CUDA device, if NULL passed, no device name is set. 
- * \returns             TRUE if the device is supported, otherwise FALSE
- * 
- */
-gmx_bool is_gmx_openmm_supported_gpu(int dev_id, char *gpu_name)
-{
-    cudaDeviceProp dev_prop;
-
-    if (debug) fprintf(debug, "Checking compatibility with device #%d, %s\n", dev_id, gpu_name);
-
-    if (do_sanity_checks(dev_id, &dev_prop) != 0)
-        return -1;
-
-    if (gpu_name != NULL)
-    { 
-        strcpy(gpu_name, dev_prop.name);
-    }
-    return is_gmx_openmm_supported_gpu_name(dev_prop.name);
-}
-
-
 /*!
  * \brief Runs a set of memory tests specified by the given bit-flags.
  * Tries to allocate and do the test on \p megs Mb memory or
index 459770df281e97c2fdd6f5012b66a616c79e7a28..5ab860dceeb40a62f34db6bf655331881861447a 100644 (file)
@@ -72,8 +72,7 @@ set(GMXPREPROCESS_SOURCES
 set(MDRUN_SOURCES 
     gctio.c    ionize.c runner.c
     do_gct.c     repl_ex.c  xutils.c pme_loadbal.c
-    md.c         mdrun.c    genalg.c membed.c
-    md_openmm.c)
+    md.c         mdrun.c    genalg.c membed.c)
 
 add_library(gmxpreprocess ${GMXPREPROCESS_SOURCES})
 target_link_libraries(gmxpreprocess md)
@@ -81,23 +80,18 @@ set_target_properties(gmxpreprocess PROPERTIES OUTPUT_NAME "gmxpreprocess${GMX_L
     COMPILE_FLAGS "${OpenMP_C_FLAGS}")
 
 
-if(GMX_GPU OR GMX_OPENMM) 
+if(GMX_GPU)
     include_directories(${CMAKE_SOURCE_DIR}/src/gmxlib/gpu_utils)
 endif()
 
-if(GMX_OPENMM) 
-    include_directories(${OpenMM_INCLUDE_DIR})
-    link_directories(${OpenMM_LIBRARY_DIR}) 
-    # with this define no evn.var. is needed with OPENMM_PLUGIN_DIR
-    # if the same OpenMM installation is used for running and building 
-    add_definitions( -DOPENMM_PLUGIN_DIR="${OpenMM_PLUGIN_DIR}" ) 
-    file(TO_CMAKE_PATH ${OpenMM_PLUGIN_DIR} _path)
-    add_library(openmm_api_wrapper STATIC openmm_wrapper.cpp)
-    target_link_libraries(openmm_api_wrapper ${OpenMM_LIBRARIES})
-    set(GMX_OPENMM_LIBRARIES openmm_api_wrapper ${OpenMM_LIBRARIES})   
+if(GMX_OPENMM)
+    # Even though the OpenMM build has "moved to contrib", many things
+    # have be be done from within the scope of the CMakeLists.txt that
+    # builds its mdrun, and that is here
+    include(../contrib/BuildMdrunOpenMM)
 endif(GMX_OPENMM)
 
-if(GMX_GPU OR GMX_OPENMM OR GMX_FORCE_CXX)
+if(GMX_GPU OR GMX_FORCE_CXX)
     set_source_files_properties(main.c PROPERTIES LANGUAGE CXX)
     if ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang")
         set_source_files_properties(main.c PROPERTIES COMPILE_FLAGS "-x c++")
@@ -128,16 +122,9 @@ endforeach()
 
 add_executable(mdrun ${MDRUN_SOURCES} main.c)
 gmx_add_man_page(mdrun)
-target_link_libraries(mdrun ${GMX_EXTRA_LIBRARIES} ${GMX_OPENMM_LIBRARIES})
+target_link_libraries(mdrun ${GMX_EXTRA_LIBRARIES})
 set_target_properties(mdrun PROPERTIES OUTPUT_NAME "mdrun${GMX_BINARY_SUFFIX}" COMPILE_FLAGS "${OpenMP_C_FLAGS}")
 
-# this is to circumvent the following MSVC error: 
-# warning LNK4098: defaultlib 'LIBCMT' conflicts with use of other libs
-# fatal error LNK1169: one or more multiply defined symbols found
-if(GMX_OPENMM AND MSVC)
-    set_target_properties(mdrun PROPERTIES LINK_FLAGS "/NODEFAULTLIB:LIBCMT")
-endif()
-
 # Construct component groups for installation; note that a component may
 # belong to only one group
 foreach(PROGRAM ${GMX_KERNEL_PROGRAMS})
@@ -180,7 +167,7 @@ endforeach()
 install(TARGETS gmxpreprocess DESTINATION ${LIB_INSTALL_DIR} COMPONENT libraries-gmxpreprocess)
 
 if (INSTALL_CUDART_LIB) #can be set manual by user
-    if (GMX_OPENMM OR GMX_GPU)
+    if (GMX_GPU)
         foreach(CUDA_LIB ${CUDA_LIBRARIES})
             string(REGEX MATCH "cudart" IS_CUDART ${CUDA_LIB})
             if(IS_CUDART) #libcuda should not be installed
@@ -191,7 +178,7 @@ if (INSTALL_CUDART_LIB) #can be set manual by user
             endif()
         endforeach()
     else()
-        message(WARNING "INSTALL_CUDART_LIB only makes sense with GMX_OPENMM or GMX_GPU")
+        message(WARNING "INSTALL_CUDART_LIB only makes sense with GMX_GPU")
     endif()
 endif ()
 endif(GMX_FAHCORE)
index 5d80e8422ba3c27f1c27e40fb201f9fe47b84f84..5c6a23b9c30e80b7e69335b2b1c633faea93f381 100644 (file)
 int cmain(int argc,char *argv[])
 {
   const char *desc[] = {
- #ifdef GMX_OPENMM
-    "This is an experimental release of GROMACS for accelerated",
-       "Molecular Dynamics simulations on GPU processors. Support is provided",
-       "by the OpenMM library (https://simtk.org/home/openmm).[PAR]",
-       "*Warning*[BR]",
-       "This release is targeted at developers and advanced users and",
-       "care should be taken before production use. The following should be",
-       "noted before using the program:[PAR]",
-       " * The current release runs only on modern nVidia GPU hardware with CUDA support.",
-       "Make sure that the necessary CUDA drivers and libraries for your operating system",
-       "are already installed. The CUDA SDK also should be installed in order to compile",
-       "the program from source (http://www.nvidia.com/object/cuda_home.html).[PAR]",
-       " * Multiple GPU cards are not supported.[PAR]",
-       " * Only a small subset of the GROMACS features and options are supported on the GPUs.",
-       "See below for a detailed list.[PAR]",
-       " * Consumer level GPU cards are known to often have problems with faulty memory.",
-       "It is recommended that a full memory check of the cards is done at least once",
-       "(for example, using the memtest=full option).",
-       "A partial memory check (for example, memtest=15) before and",
-       "after the simulation run would help spot",
-       "problems resulting from processor overheating.[PAR]",
-       " * The maximum size of the simulated systems depends on the available",
-       "GPU memory,for example, a GTX280 with 1GB memory has been tested with systems",
-       "of up to about 100,000 atoms.[PAR]",
-       " * In order to take a full advantage of the GPU platform features, many algorithms",
-       "have been implemented in a very different way than they are on the CPUs.",
-       "Therefore numercal correspondence between properties of the state of",
-       "simulated systems should not be expected. Moreover, the values will likely vary",
-       "when simulations are done on different GPU hardware.[PAR]",
-       " * Frequent retrieval of system state information such as",
-       "trajectory coordinates and energies can greatly influence the performance",
-       "of the program due to slow CPU<->GPU memory transfer speed.[PAR]",
-       " * MD algorithms are complex, and although the Gromacs code is highly tuned for them,",
-       "they often do not translate very well onto the streaming architetures.",
-       "Realistic expectations about the achievable speed-up from test with GTX280:",
-       "For small protein systems in implicit solvent using all-vs-all kernels the acceleration",
-       "can be as high as 20 times, but in most other setups involving cutoffs and PME the",
-       "acceleration is usually only ~4 times relative to a 3GHz CPU.[PAR]",
-       "Supported features:[PAR]",
-       " * Integrators: md/md-vv/md-vv-avek, sd/sd1 and bd.\n",
-       " * Long-range interactions (option coulombtype): Reaction-Field, Ewald, PME, and cut-off (for Implicit Solvent only)\n",
-       " * Temperature control: Supported only with the md/md-vv/md-vv-avek, sd/sd1 and bd integrators.\n",
-       " * Pressure control: Supported.\n",
-       " * Implicit solvent: Supported.\n",
-       "A detailed description can be found on the GROMACS website:\n",
-       "http://www.gromacs.org/gpu[PAR]",
-/* From the original mdrun documentaion */
-    "The [TT]mdrun[tt] program reads the run input file ([TT]-s[tt])",
-    "and distributes the topology over nodes if needed.",
-    "[TT]mdrun[tt] produces at least four output files.",
-    "A single log file ([TT]-g[tt]) is written, unless the option",
-    "[TT]-seppot[tt] is used, in which case each node writes a log file.",
-    "The trajectory file ([TT]-o[tt]), contains coordinates, velocities and",
-    "optionally forces.",
-    "The structure file ([TT]-c[tt]) contains the coordinates and",
-    "velocities of the last step.",
-    "The energy file ([TT]-e[tt]) contains energies, the temperature,",
-    "pressure, etc, a lot of these things are also printed in the log file.",
-    "Optionally coordinates can be written to a compressed trajectory file",
-    "([TT]-x[tt]).[PAR]",
-/* openmm specific information */
-       "Usage with OpenMM:[BR]",
-       "[TT]mdrun -device \"OpenMM:platform=Cuda,memtest=15,deviceid=0,force-device=no\"[tt][PAR]",
-       "Options:[PAR]",
-       "      [TT]platform[tt] = Cuda\t\t:\tThe only available value. OpenCL support will be available in future.\n",
-       "      [TT]memtest[tt] = 15\t\t:\tRun a partial, random GPU memory test for the given amount of seconds. A full test",
-       "(recommended!) can be run with \"memtest=full\". Memory testing can be disabled with \"memtest=off\".\n",
-       "      [TT]deviceid[tt] = 0\t\t:\tSpecify the target device when multiple cards are present.",
-       "Only one card can be used at any given time though.\n",
-       "      [TT]force-device[tt] = no\t\t:\tIf set to \"yes\" [TT]mdrun[tt]  will be forced to execute on",
-       "hardware that is not officially supported. GPU acceleration can also be achieved on older",
-       "but Cuda capable cards, although the simulation might be too slow, and the memory limits too strict.",
-#else
     "The [TT]mdrun[tt] program is the main computational chemistry engine",
     "within GROMACS. Obviously, it performs Molecular Dynamics simulations,",
     "but it can also perform Stochastic Dynamics, Energy Minimization,",
@@ -429,7 +356,6 @@ int cmain(int argc,char *argv[])
     "the [TT]mdrun[tt] process that is the parent of the others.",
     "[PAR]",
     "When [TT]mdrun[tt] is started with MPI, it does not run niced by default."
-#endif
   };
   t_commrec    *cr;
   t_filenm fnm[] = {
@@ -608,11 +534,6 @@ int cmain(int argc,char *argv[])
       "HIDDENReset cycle counters after these many time steps" },
     { "-resethway", FALSE, etBOOL, {&bResetCountersHalfWay},
       "HIDDENReset the cycle counters after half the number of steps or halfway [TT]-maxh[tt]" }
-#ifdef GMX_OPENMM
-    ,
-    { "-device",  FALSE, etSTR, {&deviceOptions},
-      "Device option string" }
-#endif
   };
   gmx_edsam_t  ed;
   unsigned long Flags, PCA_Flags;
index 5b381fd8e5c58e2fd932d6b501700a4d2bf0c94a..317b04433d92388d7781ad5bfb141802eee30b2a 100644 (file)
@@ -1041,10 +1041,6 @@ gmx_membed_t init_membed(FILE *fplog, int nfile, const t_filenm fnm[], gmx_mtop_
             gmx_input("Sorry, parallel g_membed is not yet fully functional.");
         }
 
-#ifdef GMX_OPENMM
-        gmx_input("Sorry, g_membed does not work with openmm.");
-#endif
-
         if(*cpt>=0)
         {
             fprintf(stderr,"\nSetting -cpt to -1, because embedding cannot be restarted from cpt-files.\n");
index d6cfe63e22f94524c067289d377ac39bd97ac8a9..8c640098d2c9f5ff2841f95d0944ed618822c225 100644 (file)
@@ -83,7 +83,6 @@
 #include "../mdlib/nbnxn_consts.h"
 #include "gmx_fatal_collective.h"
 #include "membed.h"
-#include "md_openmm.h"
 #include "gmx_omp.h"
 
 #include "thread_mpi/threads.h"
 #include "corewrap.h"
 #endif
 
-#ifdef GMX_OPENMM
-#include "md_openmm.h"
-#endif
-
 #include "gpu_utils.h"
 #include "nbnxn_cuda_data_mgmt.h"
 
@@ -111,11 +106,7 @@ typedef struct {
 } gmx_intp_t;
 
 /* The array should match the eI array in include/types/enums.h */
-#ifdef GMX_OPENMM  /* FIXME do_md_openmm needs fixing */
-const gmx_intp_t integrator[eiNR] = { {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm}, {do_md_openmm},{do_md_openmm}};
-#else
 const gmx_intp_t integrator[eiNR] = { {do_md}, {do_steep}, {do_cg}, {do_md}, {do_md}, {do_nm}, {do_lbfgs}, {do_tpi}, {do_tpi}, {do_md}, {do_md},{do_md}};
-#endif
 
 gmx_large_int_t     deform_init_init_step_tpx;
 matrix              deform_init_box_tpx;
@@ -1848,12 +1839,7 @@ int mdrunner(gmx_hw_opt_t *hw_opt,
     }
 
 
-    if (integrator[inputrec->eI].func == do_md
-#ifdef GMX_OPENMM
-        ||
-        integrator[inputrec->eI].func == do_md_openmm
-#endif
-        )
+    if (integrator[inputrec->eI].func == do_md)
     {
         /* Turn on signal handling on all nodes */
         /*
index fa95936b0564c45dcffa81f01b119e14c3a526d4..c2cd494e18423e886835be8f0779390a0f17052e 100644 (file)
@@ -181,6 +181,9 @@ t_mdebin *init_mdebin(ener_file_t fp_ene,
         md->bEInd[i]=FALSE;
     }
 
+    /* Even though the OpenMM build has moved to contrib, it's not
+     * practical to move/remove this code fragment, because of the
+     * fundamental mess that is the GROMACS library structure. */
 #ifndef GMX_OPENMM
     for(i=0; i<F_NRE; i++)
     {
index a158b5e0811a39dca01d6c77ddab89ffb716b37b..56800cb18db29e74cf32164ff1479e969dbb34c0 100644 (file)
@@ -100,7 +100,7 @@ set(GMX_TOOLS_PROGRAMS_NOT_FOR_INSTALLATION
 
 
 foreach(TOOL ${GMX_TOOLS_PROGRAMS} ${GMX_TOOLS_PROGRAMS_NOT_FOR_INSTALLATION})
-    if(GMX_GPU OR GMX_OPENMM OR GMX_FORCE_CXX)
+    if(GMX_GPU OR GMX_FORCE_CXX)
         set_source_files_properties(${TOOL}.c PROPERTIES LANGUAGE CXX)
         if ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang")
             set_source_files_properties(${TOOL}.c PROPERTIES COMPILE_FLAGS "-x c++")
index 25944666c5abea0f730c5109291ce9155daf3cf6..87813fe86219f032f5b5f94ccfc2d8ad371bf73a 100644 (file)
@@ -121,11 +121,6 @@ int gmx_membed(int argc,char *argv[])
     gmx_bool bVerbose=FALSE;
     char *mdrun_path=NULL;
 
-/* arguments relevant to OPENMM only*/
-#ifdef GMX_OPENMM
-    gmx_fatal(FARGS,"g_membed not implemented for openmm");
-#endif
-
     t_pargs pa[] = {
         { "-xyinit",   FALSE, etREAL,  {&xy_fac},
             "Resize factor for the protein in the xy dimension before starting embedding" },