Add Cray compiler support
authorDaniel Landau <dlandau@cray.com>
Mon, 24 Mar 2014 16:01:22 +0000 (11:01 -0500)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Wed, 26 Mar 2014 18:37:23 +0000 (19:37 +0100)
 * Add support for SIMD intstructions with the Cray compiler with
   -hgnu. This is a general GNU compatibility mode so there
   also has to be a guard in the atomics section in threadMPI.
 * Add support for atomics with the Cray compiler.
 * Add Cray compiler support for cycle counting.
 * OpenMP detection for the Cray compiler in CMake before
   version 3.0 is broken, so refuse to proceed with that
   combination.

Change-Id: I1209f3a61ed1d4d49e56dacbb315235088c5eea6

cmake/gmxManageOpenMP.cmake
cmake/gmxTestSimd.cmake
src/gromacs/legacyheaders/thread_mpi/atomic.h
src/gromacs/legacyheaders/thread_mpi/atomic/cce.h [new file with mode: 0644]
src/gromacs/legacyheaders/thread_mpi/atomic/cce_intrinsics.h [new file with mode: 0644]
src/gromacs/legacyheaders/thread_mpi/atomic/cce_spinlock.h [new file with mode: 0644]
src/gromacs/timing/cyclecounter.h

index 9993497ebfe37853effe73610541614005545146..0fa2faeb77a7970278b265c099f7854b1de44cfe 100644 (file)
@@ -45,6 +45,10 @@ if(GMX_OPENMP)
         message(STATUS "OpenMP multithreading not supported with gcc/llvm-gcc 4.2 on Mac OS X, disabled")
         set(GMX_OPENMP OFF CACHE BOOL
             "OpenMP multithreading not not supported with gcc/llvm-gcc 4.2 on Mac OS X, disabled!" FORCE)
+    elseif(CMAKE_C_COMPILER_ID MATCHES "Cray" AND CMAKE_VERSION VERSION_LESS 3)
+        message(STATUS "OpenMP multithreading is not detected correctly for the Cray compiler with CMake before version 3.0 (see http://public.kitware.com/Bug/view.php?id=14567)")
+        set(GMX_OPENMP OFF CACHE BOOL
+            "OpenMP multithreading is not detected correctly for the Cray compiler with CMake before version 3.0 (see http://public.kitware.com/Bug/view.php?id=14567)" FORCE)
     else()
         # We should do OpenMP detection if we get here
         # OpenMP check must come before other CFLAGS!
index bdf4e708e7509e7412b26b12f768fb32a6a91aa6..1137db92ecdd9769ab158929b7e49459e82dfdd0 100644 (file)
@@ -68,12 +68,12 @@ elseif(${GMX_SIMD} STREQUAL "SSE2")
                               "#include<xmmintrin.h>
                               int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_rsqrt_ps(x);return 0;}"
                               SIMD_C_FLAGS
-                              "-msse2" "/arch:SSE2")
+                              "-msse2" "/arch:SSE2" "-hgnu")
     gmx_find_cxxflag_for_source(CXXFLAGS_SSE2 "C++ compiler SSE2 flag"
                                 "#include<xmmintrin.h>
                                 int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_rsqrt_ps(x);return 0;}"
                                 SIMD_CXX_FLAGS
-                                "-msse2" "/arch:SSE2")
+                                "-msse2" "/arch:SSE2" "-hgnu")
 
     if(NOT CFLAGS_SSE2 OR NOT CXXFLAGS_SSE2)
         message(FATAL_ERROR "Cannot find SSE2 compiler flag. Use a newer compiler, or disable SIMD (slower).")
@@ -89,12 +89,12 @@ elseif(${GMX_SIMD} STREQUAL "SSE4.1")
                               "#include<smmintrin.h>
                               int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_dp_ps(x,x,0x77);return 0;}"
                               SIMD_C_FLAGS
-                              "-msse4.1" "/arch:SSE4.1" "/arch:SSE2")
+                              "-msse4.1" "/arch:SSE4.1" "/arch:SSE2" "-hgnu")
     gmx_find_cxxflag_for_source(CXXFLAGS_SSE4_1 "C++ compiler SSE4.1 flag"
                                 "#include<smmintrin.h>
                                 int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_dp_ps(x,x,0x77);return 0;}"
                                 SIMD_CXX_FLAGS
-                                "-msse4.1" "/arch:SSE4.1" "/arch:SSE2")
+                                "-msse4.1" "/arch:SSE4.1" "/arch:SSE2" "-hgnu")
 
     if(NOT CFLAGS_SSE4_1 OR NOT CXXFLAGS_SSE4_1)
         message(FATAL_ERROR "Cannot find SSE4.1 compiler flag. "
@@ -122,12 +122,12 @@ elseif(${GMX_SIMD} STREQUAL "AVX_128_FMA")
                               "#include<immintrin.h>
                               int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_permute_ps(x,1);return 0;}"
                               SIMD_C_FLAGS
-                              "-mavx" "/arch:AVX")
+                              "-mavx" "/arch:AVX" "-hgnu")
     gmx_find_cxxflag_for_source(CXXFLAGS_AVX_128 "C++ compiler AVX (128 bit) flag"
                                 "#include<immintrin.h>
                                 int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_permute_ps(x,1);return 0;}"
                                 SIMD_CXX_FLAGS
-                                "-mavx" "/arch:AVX")
+                                "-mavx" "/arch:AVX" "-hgnu")
 
     ### STAGE 2: Find the fused-multiply add flag.
     # GCC requires x86intrin.h for FMA support. MSVC 2010 requires intrin.h for FMA support.
@@ -146,14 +146,14 @@ ${INCLUDE_X86INTRIN_H}
 ${INCLUDE_INTRIN_H}
 int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_macc_ps(x,x,x);return 0;}"
                               SIMD_C_FLAGS
-                              "-mfma4")
+                              "-mfma4" "-hgnu")
     gmx_find_cxxflag_for_source(CXXFLAGS_AVX_128_FMA "C++ compiler AVX (128 bit) FMA4 flag"
 "#include<immintrin.h>
 ${INCLUDE_X86INTRIN_H}
 ${INCLUDE_INTRIN_H}
 int main(){__m128 x=_mm_set1_ps(0.5);x=_mm_macc_ps(x,x,x);return 0;}"
                                 SIMD_CXX_FLAGS
-                                "-mfma4")
+                                "-mfma4" "-hgnu")
 
     # We only need to check the last (FMA) test; that will always fail if the basic AVX128 test failed
     if(NOT CFLAGS_AVX_128_FMA OR NOT CXXFLAGS_AVX_128_FMA)
@@ -209,12 +209,12 @@ elseif(${GMX_SIMD} STREQUAL "AVX_256")
                               "#include<immintrin.h>
                               int main(){__m256 x=_mm256_set1_ps(0.5);x=_mm256_add_ps(x,x);return 0;}"
                               SIMD_C_FLAGS
-                              "-mavx" "/arch:AVX")
+                              "-mavx" "/arch:AVX" "-hgnu")
     gmx_find_cxxflag_for_source(CXXFLAGS_AVX "C++ compiler AVX (256 bit) flag"
                                 "#include<immintrin.h>
                                 int main(){__m256 x=_mm256_set1_ps(0.5);x=_mm256_add_ps(x,x);return 0;}"
                                 SIMD_CXX_FLAGS
-                                "-mavx" "/arch:AVX")
+                                "-mavx" "/arch:AVX" "-hgnu")
 
     if(NOT CFLAGS_AVX OR NOT CXXFLAGS_AVX)
         message(FATAL_ERROR "Cannot find AVX compiler flag. Use a newer compiler, or choose SSE4.1 SIMD (slower).")
@@ -233,12 +233,12 @@ elseif(${GMX_SIMD} STREQUAL "AVX2_256")
                               "#include<immintrin.h>
                               int main(){__m256 x=_mm256_set1_ps(0.5);x=_mm256_fmadd_ps(x,x,x);return 0;}"
                               SIMD_C_FLAGS
-                              "-march=core-avx2" "-mavx2" "/arch:AVX") # no AVX2-specific flag for MSVC yet
+                              "-march=core-avx2" "-mavx2" "/arch:AVX" "-hgnu") # no AVX2-specific flag for MSVC yet
     gmx_find_cxxflag_for_source(CXXFLAGS_AVX2 "C++ compiler AVX2 flag"
                                 "#include<immintrin.h>
                                 int main(){__m256 x=_mm256_set1_ps(0.5);x=_mm256_fmadd_ps(x,x,x);return 0;}"
                                 SIMD_CXX_FLAGS
-                                "-march=core-avx2" "-mavx2" "/arch:AVX") # no AVX2-specific flag for MSVC yet
+                                "-march=core-avx2" "-mavx2" "/arch:AVX" "-hgnu") # no AVX2-specific flag for MSVC yet
 
     if(NOT CFLAGS_AVX2 OR NOT CXXFLAGS_AVX2)
         message(FATAL_ERROR "Cannot find AVX2 compiler flag. Use a newer compiler, or choose AVX SIMD (slower).")
index 0185bcbb20790189598ee09def0049f9320fc4db..5aa81630f68651bba660b80dcb74f5b94c31644e 100644 (file)
@@ -107,7 +107,7 @@ extern "C"
    Some compatible compilers, like icc on linux+mac will take this path,
    too */
 #if ( (defined(__GNUC__) || defined(__PATHSCALE__) || defined(__PGI)) && \
-    (!defined(__xlc__)) && (!defined(TMPI_TEST_NO_ATOMICS)) )
+    (!defined(__xlc__)) && (!defined(_CRAYC)) && (!defined(TMPI_TEST_NO_ATOMICS)) )
 
 #ifdef __GNUC__
 #define TMPI_GCC_VERSION (__GNUC__ * 10000 \
@@ -171,7 +171,10 @@ extern "C"
 /* Fujitsu FX10 SPARC compiler requires gcc compatibility with -Xg */
 #error Atomics support for Fujitsu FX10 compiler requires -Xg (gcc compatibility)
 
+#elif defined(_CRAYC)
 
+/* Cray compiler */
+#include "atomic/cce.h"
 #else
 
 #ifndef DOXYGEN
diff --git a/src/gromacs/legacyheaders/thread_mpi/atomic/cce.h b/src/gromacs/legacyheaders/thread_mpi/atomic/cce.h
new file mode 100644 (file)
index 0000000..9030730
--- /dev/null
@@ -0,0 +1,69 @@
+/*
+   This source code file is part of thread_mpi.
+   Original for gcc written by Sander Pronk, Erik Lindahl, and possibly
+   others. Modified for the Cray compiler by Daniel Landau.
+
+   Copyright (c) 2009, Sander Pronk, Erik Lindahl.
+   Copyright 2014, Cray Inc.
+   All rights reserved.
+
+   Redistribution and use in source and binary forms, with or without
+   modification, are permitted provided that the following conditions are met:
+   1) Redistributions of source code must retain the above copyright
+   notice, this list of conditions and the following disclaimer.
+   2) Redistributions in binary form must reproduce the above copyright
+   notice, this list of conditions and the following disclaimer in the
+   documentation and/or other materials provided with the distribution.
+   3) Neither the name of the copyright holders nor the
+   names of its contributors may be used to endorse or promote products
+   derived from this software without specific prior written permission.
+
+   THIS SOFTWARE IS PROVIDED BY US ''AS IS'' AND ANY
+   EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
+   WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+   DISCLAIMED. IN NO EVENT SHALL WE BE LIABLE FOR ANY
+   DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
+   (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
+   LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
+   ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+   SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+   If you want to redistribute modifications, 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 should not
+   be called official thread_mpi. Details are found in the README & COPYING
+   files.
+ */
+
+#include <intrinsics.h>
+
+#define tMPI_Atomic_memory_barrier() __builtin_ia32_mfence()
+
+
+typedef struct tMPI_Atomic
+{
+    volatile long value;
+}
+tMPI_Atomic_t;
+
+typedef struct tMPI_Atomic_ptr
+{
+    volatile void* value;
+}
+tMPI_Atomic_ptr_t;
+
+
+/* these are guaranteed to be  atomic on x86 and x86_64 */
+#define tMPI_Atomic_get(a)  ((int)( (a)->value) )
+#define tMPI_Atomic_set(a, i)  (((a)->value) = (i))
+
+
+#define tMPI_Atomic_ptr_get(a)  ((void*)((a)->value) )
+#define tMPI_Atomic_ptr_set(a, i)  (((a)->value) = (void*)(i))
+
+
+#include "cce_intrinsics.h"
+
+#include "cce_spinlock.h"
diff --git a/src/gromacs/legacyheaders/thread_mpi/atomic/cce_intrinsics.h b/src/gromacs/legacyheaders/thread_mpi/atomic/cce_intrinsics.h
new file mode 100644 (file)
index 0000000..40934e5
--- /dev/null
@@ -0,0 +1,71 @@
+/*
+   This source code file is part of thread_mpi.
+   Original for gcc written by Sander Pronk, Erik Lindahl, and possibly
+   others. Modified for the Cray compiler by Daniel Landau.
+
+
+   Copyright (c) 2009, Sander Pronk, Erik Lindahl.
+   Copyright 2014, Cray Inc.
+   All rights reserved.
+
+   Redistribution and use in source and binary forms, with or without
+   modification, are permitted provided that the following conditions are met:
+   1) Redistributions of source code must retain the above copyright
+   notice, this list of conditions and the following disclaimer.
+   2) Redistributions in binary form must reproduce the above copyright
+   notice, this list of conditions and the following disclaimer in the
+   documentation and/or other materials provided with the distribution.
+   3) Neither the name of the copyright holders nor the
+   names of its contributors may be used to endorse or promote products
+   derived from this software without specific prior written permission.
+
+   THIS SOFTWARE IS PROVIDED BY US ''AS IS'' AND ANY
+   EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
+   WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+   DISCLAIMED. IN NO EVENT SHALL WE BE LIABLE FOR ANY
+   DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
+   (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
+   LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
+   ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+   SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+   If you want to redistribute modifications, 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 should not
+   be called official thread_mpi. Details are found in the README & COPYING
+   files.
+ */
+
+#include <intrinsics.h>
+
+#define tMPI_Atomic_memory_barrier() __builtin_ia32_mfence()
+
+TMPI_EXPORT
+static inline int tMPI_Atomic_cas(tMPI_Atomic_t *a, int oldval, int newval)
+{
+    return __sync_val_compare_and_swap(&(a->value), oldval, newval) == oldval;
+}
+
+TMPI_EXPORT
+static inline int tMPI_Atomic_ptr_cas(tMPI_Atomic_ptr_t* a, void *oldval,
+                                      void *newval)
+{
+    return __sync_val_compare_and_swap((size_t*)&(a->value), (size_t)oldval, (size_t)newval) == (size_t)oldval;
+}
+
+TMPI_EXPORT
+static inline int tMPI_Atomic_add_return(tMPI_Atomic_t *a, volatile int i)
+{
+    return __sync_add_and_fetch( &(a->value), i);
+}
+#define TMPI_ATOMIC_HAVE_NATIVE_ADD_RETURN
+
+
+TMPI_EXPORT
+static inline int tMPI_Atomic_fetch_add(tMPI_Atomic_t *a, volatile int i)
+{
+    return __sync_fetch_and_add( &(a->value), i);
+}
+#define TMPI_ATOMIC_HAVE_NATIVE_FETCH_ADD
diff --git a/src/gromacs/legacyheaders/thread_mpi/atomic/cce_spinlock.h b/src/gromacs/legacyheaders/thread_mpi/atomic/cce_spinlock.h
new file mode 100644 (file)
index 0000000..31de26f
--- /dev/null
@@ -0,0 +1,93 @@
+/*
+   This source code file is part of thread_mpi.
+   Original for gcc written by Sander Pronk, Erik Lindahl, and possibly
+   others. Modified for the Cray compiler by Daniel Landau.
+
+   Copyright (c) 2009, Sander Pronk, Erik Lindahl.
+   Copyright 2014, Cray Inc.
+   All rights reserved.
+
+   Redistribution and use in source and binary forms, with or without
+   modification, are permitted provided that the following conditions are met:
+   1) Redistributions of source code must retain the above copyright
+   notice, this list of conditions and the following disclaimer.
+   2) Redistributions in binary form must reproduce the above copyright
+   notice, this list of conditions and the following disclaimer in the
+   documentation and/or other materials provided with the distribution.
+   3) Neither the name of the copyright holders nor the
+   names of its contributors may be used to endorse or promote products
+   derived from this software without specific prior written permission.
+
+   THIS SOFTWARE IS PROVIDED BY US ''AS IS'' AND ANY
+   EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
+   WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+   DISCLAIMED. IN NO EVENT SHALL WE BE LIABLE FOR ANY
+   DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
+   (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
+   LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
+   ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+   SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+   If you want to redistribute modifications, 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 should not
+   be called official thread_mpi. Details are found in the README & COPYING
+   files.
+ */
+
+#include <intrinsics.h>
+
+typedef struct tMPI_Spinlock
+{
+    volatile long lock /*__attribute__ ((aligned(64)))*/;
+} tMPI_Spinlock_t;
+
+#define TMPI_SPINLOCK_INITIALIZER   { 0 }
+
+#define TMPI_ATOMIC_HAVE_NATIVE_SPINLOCK
+
+
+
+static inline void tMPI_Spinlock_init(tMPI_Spinlock_t *x)
+{
+    x->lock = 0;
+}
+
+
+static inline void tMPI_Spinlock_lock(tMPI_Spinlock_t *x)
+{
+    while (__sync_lock_test_and_set(&(x->lock), 1) == 1)
+    {
+        /* this is nicer on the system bus: */
+        while (x->lock == 1)
+        {
+        }
+    }
+}
+
+
+static inline int tMPI_Spinlock_trylock(tMPI_Spinlock_t *x)
+{
+    return __sync_lock_test_and_set(&(x->lock), 1);
+}
+
+
+static inline void tMPI_Spinlock_unlock(tMPI_Spinlock_t *x)
+{
+    x->lock = 0;
+}
+
+static inline int tMPI_Spinlock_islocked(const tMPI_Spinlock_t *x)
+{
+    return ( x->lock == 1 );
+}
+
+static inline void tMPI_Spinlock_wait(tMPI_Spinlock_t *x)
+{
+    do
+    {
+    }
+    while (x->lock == 1);
+}
index b551009476bd69e42cd2bcf785ffe37cb196611e..87d19b19a1c3865653b314e789343fd9d810af41 100644 (file)
@@ -203,7 +203,7 @@ typedef long
  *       one when later linking to the library it might happen that the
  *       library supports cyclecounters but not the headers, or vice versa.
  */
-#if ((defined(__GNUC__) || defined(__INTEL_COMPILER) || defined(__PATHSCALE__) || defined(__PGIC__)) && \
+#if ((defined(__GNUC__) || defined(__INTEL_COMPILER) || defined(__PATHSCALE__) || defined(__PGIC__) || defined(_CRAYC)) && \
     (defined(__i386__) || defined(__x86_64__)))
 static __inline__ int gmx_cycles_have_counter(void)
 {
@@ -333,7 +333,7 @@ static int gmx_cycles_have_counter(void)
  *  routine.
  */
 #if ((defined(__GNUC__) || defined(__INTEL_COMPILER) || defined(__PATHSCALE__) || defined(__PGIC__)) && \
-    (defined(__i386__) || defined(__x86_64__)))
+    (defined(__i386__) || defined(__x86_64__)) && !defined(_CRAYC))
 static __inline__ gmx_cycles_t gmx_cycles_read(void)
 {
     /* x86 with GCC inline assembly - pentium TSC register */
@@ -505,6 +505,13 @@ static __inline__ gmx_cycles_t gmx_cycles_read(void)
     return ret;
 }
 
+#elif defined(_CRAYC)
+#include <intrinsics.h>
+
+static __inline gmx_cycles_t gmx_cycles_read(void)
+{
+    return _rtc();
+}
 #else
 static gmx_cycles_t gmx_cycles_read(void)
 {