From f31c4619c4fab9f13740e133e8bdabebf508a40a Mon Sep 17 00:00:00 2001 From: Roland Schulz Date: Mon, 25 Nov 2013 16:06:26 -0500 Subject: [PATCH] Add bitmask type Adds support for arbirtrary number of OpenMP threads Fixes #1386 Change-Id: I3747b4c7cd4b45d5901f710fcab47984c8913bc6 --- cmake/gmxManageOpenMP.cmake | 8 + src/config.h.cmakein | 3 + .../gmock-1.7.0/gtest/include/gtest/gtest.h | 2 + src/gromacs/gmxlib/copyrite.cpp | 2 +- src/gromacs/mdlib/bitmask.h | 196 ++++++++++++++++++ src/gromacs/mdlib/clincs.cpp | 64 +++--- src/gromacs/mdlib/nbnxn_atomdata.c | 12 +- .../mdlib/nbnxn_kernels/nbnxn_kernel_common.c | 6 +- src/gromacs/mdlib/nbnxn_pairlist.h | 18 +- src/gromacs/mdlib/nbnxn_search.c | 24 ++- src/gromacs/mdlib/tests/CMakeLists.txt | 3 +- src/gromacs/mdlib/tests/bitmask.h | 132 ++++++++++++ src/gromacs/mdlib/tests/bitmask128.cpp | 40 ++++ src/gromacs/mdlib/tests/bitmask32.cpp | 40 ++++ src/gromacs/mdlib/tests/bitmask64.cpp | 40 ++++ 15 files changed, 521 insertions(+), 69 deletions(-) create mode 100644 src/gromacs/mdlib/bitmask.h create mode 100644 src/gromacs/mdlib/tests/bitmask.h create mode 100644 src/gromacs/mdlib/tests/bitmask128.cpp create mode 100644 src/gromacs/mdlib/tests/bitmask32.cpp create mode 100644 src/gromacs/mdlib/tests/bitmask64.cpp diff --git a/cmake/gmxManageOpenMP.cmake b/cmake/gmxManageOpenMP.cmake index d1a8e4dca5..ddcf7c4c9f 100644 --- a/cmake/gmxManageOpenMP.cmake +++ b/cmake/gmxManageOpenMP.cmake @@ -81,3 +81,11 @@ if(GMX_OPENMP) endif() endif() endif() +gmx_dependent_cache_variable(GMX_OPENMP_MAX_THREADS + "Maximum number of OpenMP Threads supported. Has to be 32 or a multiple of 64." + STRING 32 GMX_OPENMP) +mark_as_advanced(GMX_OPENMP_MAX_THREADS) +math(EXPR MAX_THREAD_MOD "${GMX_OPENMP_MAX_THREADS} % 64") +if (NOT GMX_OPENMP_MAX_THREADS EQUAL 32 AND NOT ${MAX_THREAD_MOD} EQUAL 0) + message(FATAL_ERROR "Only 32 or multiples of 64 supported for GMX_OPENMP_MAX_THREADS.") +endif() diff --git a/src/config.h.cmakein b/src/config.h.cmakein index 39b643332f..2000bbb8ef 100644 --- a/src/config.h.cmakein +++ b/src/config.h.cmakein @@ -172,6 +172,9 @@ /* Can and should use nice(3) to set priority */ #cmakedefine GMX_USE_NICE +/* Maximum number of OpenMP threads supported */ +#define GMX_OPENMP_MAX_THREADS @GMX_OPENMP_MAX_THREADS@ + /* Use if can't rename checkpoints */ #cmakedefine GMX_NO_RENAME diff --git a/src/external/gmock-1.7.0/gtest/include/gtest/gtest.h b/src/external/gmock-1.7.0/gtest/include/gtest/gtest.h index 6fa0a3925e..f76d8bea1e 100644 --- a/src/external/gmock-1.7.0/gtest/include/gtest/gtest.h +++ b/src/external/gmock-1.7.0/gtest/include/gtest/gtest.h @@ -1443,6 +1443,8 @@ AssertionResult CmpHelperEQ(const char* expected_expression, # pragma warning(push) // Saves the current warning state. # pragma warning(disable:4389) // Temporarily disables warning on // signed/unsigned mismatch. +# pragma warning(disable:4805) // Temporarily disables warning on + // unsafe mix of types #endif if (expected == actual) { diff --git a/src/gromacs/gmxlib/copyrite.cpp b/src/gromacs/gmxlib/copyrite.cpp index fe7045a22e..be67e2ec48 100644 --- a/src/gromacs/gmxlib/copyrite.cpp +++ b/src/gromacs/gmxlib/copyrite.cpp @@ -667,7 +667,7 @@ static void gmx_print_version_info(FILE *fp) fprintf(fp, "MPI library: none\n"); #endif #ifdef GMX_OPENMP - fprintf(fp, "OpenMP support: enabled\n"); + fprintf(fp, "OpenMP support: enabled (GMX_OPENMP_MAX_THREADS = %d)\n", GMX_OPENMP_MAX_THREADS); #else fprintf(fp, "OpenMP support: disabled\n"); #endif diff --git a/src/gromacs/mdlib/bitmask.h b/src/gromacs/mdlib/bitmask.h new file mode 100644 index 0000000000..4b9c812513 --- /dev/null +++ b/src/gromacs/mdlib/bitmask.h @@ -0,0 +1,196 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2014, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * 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. + */ +/*! \libinternal \file + * \brief + * Declares gmx_bitmask_t and associated functions + * + * \author Roland Schulz + * \inlibraryapi + */ + +#ifndef GMX_MDLIB_BITMASK_H +#define GMX_MDLIB_BITMASK_H + +#include "config.h" /* for GMX_MAX_OPENMP_THREADS */ + +#include + +#include "gromacs/utility/basedefinitions.h" + +/*! \brief Size of bitmask. Has to be 32 or multiple of 64. */ +#ifndef BITMASK_SIZE +#define BITMASK_SIZE GMX_OPENMP_MAX_THREADS +#endif + +#if BITMASK_SIZE != 32 && BITMASK_SIZE%64 != 0 +#error BITMASK_SIZE has to be 32 or a multiple of 64. +#endif + +#if BITMASK_SIZE <= 64 || defined DOXYGEN +#if BITMASK_SIZE == 32 +typedef gmx_uint32_t gmx_bitmask_t; +#else +typedef gmx_uint64_t gmx_bitmask_t; /**< bitmask type */ +#endif + +/*! \brief Initialize all bits to 0 */ +gmx_inline static void bitmask_clear(gmx_bitmask_t* m) +{ + *m = 0; +} + +/*! \brief Set bit at position b to 1. */ +gmx_inline static void bitmask_set_bit(gmx_bitmask_t* m, int b) +{ + *m |= ((gmx_bitmask_t)1 << b); +} + +/*! \brief Initialize all bits: bit b to 1, others to 0 */ +gmx_inline static void bitmask_init_bit(gmx_bitmask_t* m, int b) +{ + *m = ((gmx_bitmask_t)1 << b); +} + +/*! \brief Initialize all bits: all bits below b to 1, others to 0 */ +gmx_inline static void bitmask_init_low_bits(gmx_bitmask_t* m, int b) +{ + *m = ((gmx_bitmask_t)1 << b) - 1; +} + +/*! \brief Test if bit b is set */ +gmx_inline static gmx_bool bitmask_is_set(gmx_bitmask_t m, int b) +{ + return (m & ((gmx_bitmask_t)1 << b)) != 0; +} + +/*! \brief Test if both bitmasks have no common bits enabled */ +gmx_inline static gmx_bool bitmask_is_disjoint(gmx_bitmask_t a, gmx_bitmask_t b) +{ + return !(a & b); +} + +/*! \brief Test if both bitmasks are equal */ +gmx_inline static gmx_bool bitmask_is_equal(gmx_bitmask_t a, gmx_bitmask_t b) +{ + return a == b; +} + +/*! \brief Test if bitmask has no enabled bits */ +gmx_inline static gmx_bool bitmask_is_zero(gmx_bitmask_t m) +{ + return !m; +} + +/*! \brief Set all bits enabled in either mask and write into a */ +gmx_inline static void bitmask_union(gmx_bitmask_t* a, gmx_bitmask_t b) +{ + *a |= b; +} +#else +#define BITMASK_ALEN (BITMASK_SIZE/64) +typedef gmx_uint64_t gmx_bitmask_t[BITMASK_ALEN]; + +gmx_inline static void bitmask_clear(gmx_bitmask_t* m) +{ + memset(*m, 0, BITMASK_SIZE/8); +} + +gmx_inline static void bitmask_set_bit(gmx_bitmask_t* m, int b) +{ + (*m)[b/64] |= ((gmx_uint64_t)1 << (b%64)); +} + +gmx_inline static void bitmask_init_bit(gmx_bitmask_t* m, int b) +{ + bitmask_clear(m); + (*m)[b/64] = ((gmx_uint64_t)1 << (b%64)); +} + +gmx_inline static void bitmask_init_low_bits(gmx_bitmask_t* m, int b) +{ + memset(*m, 255, b/64*8); + (*m)[b/64] = ((gmx_uint64_t)1 << (b%64)) - 1; + memset(&(*m)[b/64+1], 0, (BITMASK_ALEN-b/64-1)*8); +} + +gmx_inline static gmx_bool bitmask_is_set(gmx_bitmask_t m, int b) +{ + return (m[b/64] & ((gmx_uint64_t)1 << (b%64))) != 0; +} + +gmx_inline static gmx_bool bitmask_is_disjoint(gmx_bitmask_t a, gmx_bitmask_t b) +{ + int i; + gmx_bool r = 1; + for (i = 0; i < BITMASK_ALEN; i++) + { + r = r && !(a[i] & b[i]); + } + return r; +} + +gmx_inline static gmx_bool bitmask_is_equal(gmx_bitmask_t a, gmx_bitmask_t b) +{ + int i; + gmx_bool r = 1; + for (i = 0; i < BITMASK_ALEN; i++) + { + r = r && (a[i] == b[i]); + } + return r; +} + +gmx_inline static gmx_bool bitmask_is_zero(gmx_bitmask_t m) +{ + int i; + gmx_bool r = 1; + for (i = 0; i < BITMASK_ALEN; i++) + { + r = r && !m[i]; + } + return r; +} + +gmx_inline static void bitmask_union(gmx_bitmask_t* a, gmx_bitmask_t b) +{ + int i; + for (i = 0; i < BITMASK_ALEN; i++) + { + (*a)[i] |= b[i]; + } +} +#endif + +#endif diff --git a/src/gromacs/mdlib/clincs.cpp b/src/gromacs/mdlib/clincs.cpp index 32bf8fc72f..b0c700047c 100644 --- a/src/gromacs/mdlib/clincs.cpp +++ b/src/gromacs/mdlib/clincs.cpp @@ -50,6 +50,7 @@ #include "gromacs/legacyheaders/types/commrec.h" #include "gromacs/math/units.h" #include "gromacs/math/vec.h" +#include "gromacs/mdlib/bitmask.h" #include "gromacs/pbcutil/pbc.h" #include "gromacs/topology/block.h" #include "gromacs/topology/mtop_util.h" @@ -96,7 +97,7 @@ typedef struct gmx_lincsdata { real *bllen; /* the reference bond length */ int nth; /* The number of threads doing LINCS */ lincs_thread_t *th; /* LINCS thread division */ - unsigned *atf; /* atom flags for thread parallelization */ + gmx_bitmask_t *atf; /* atom flags for thread parallelization */ int atf_nalloc; /* allocation size of atf */ /* arrays for temporary storage in the LINCS algorithm */ rvec *tmpv; @@ -1016,7 +1017,7 @@ static void lincs_thread_setup(struct gmx_lincsdata *li, int natoms) { lincs_thread_t *li_m; int th; - unsigned *atf; + gmx_bitmask_t *atf; int a; if (natoms > li->atf_nalloc) @@ -1029,7 +1030,12 @@ static void lincs_thread_setup(struct gmx_lincsdata *li, int natoms) /* Clear the atom flags */ for (a = 0; a < natoms; a++) { - atf[a] = 0; + bitmask_clear(&atf[a]); + } + + if (li->nth > BITMASK_SIZE) + { + gmx_fatal(FARGS, "More than %d threads is not supported for LINCS.", BITMASK_SIZE); } for (th = 0; th < li->nth; th++) @@ -1043,14 +1049,11 @@ static void lincs_thread_setup(struct gmx_lincsdata *li, int natoms) li_th->b0 = (li->nc* th )/li->nth; li_th->b1 = (li->nc*(th+1))/li->nth; - if (th < static_cast(sizeof(*atf)*8)) + /* For each atom set a flag for constraints from each */ + for (b = li_th->b0; b < li_th->b1; b++) { - /* For each atom set a flag for constraints from each */ - for (b = li_th->b0; b < li_th->b1; b++) - { - atf[li->bla[b*2] ] |= (1U<bla[b*2+1]] |= (1U<bla[b*2]], th); + bitmask_set_bit(&atf[li->bla[b*2+1]], th); } } @@ -1058,7 +1061,7 @@ static void lincs_thread_setup(struct gmx_lincsdata *li, int natoms) for (th = 0; th < li->nth; th++) { lincs_thread_t *li_th; - unsigned mask; + gmx_bitmask_t mask; int b; li_th = &li->th[th]; @@ -1070,35 +1073,24 @@ static void lincs_thread_setup(struct gmx_lincsdata *li, int natoms) srenew(li_th->ind_r, li_th->ind_nalloc); } - if (th < static_cast(sizeof(*atf)*8)) - { - mask = (1U<nind = 0; - li_th->nind_r = 0; - for (b = li_th->b0; b < li_th->b1; b++) + li_th->nind = 0; + li_th->nind_r = 0; + for (b = li_th->b0; b < li_th->b1; b++) + { + /* We let the constraint with the lowest thread index + * operate on atoms with constraints from multiple threads. + */ + if (bitmask_is_disjoint(atf[li->bla[b*2]], mask) && + bitmask_is_disjoint(atf[li->bla[b*2+1]], mask)) { - /* We let the constraint with the lowest thread index - * operate on atoms with constraints from multiple threads. - */ - if (((atf[li->bla[b*2]] & mask) == 0) && - ((atf[li->bla[b*2+1]] & mask) == 0)) - { - /* Add the constraint to the local atom update index */ - li_th->ind[li_th->nind++] = b; - } - else - { - /* Add the constraint to the rest block */ - li_th->ind_r[li_th->nind_r++] = b; - } + /* Add the constraint to the local atom update index */ + li_th->ind[li_th->nind++] = b; } - } - else - { - /* We are out of bits, assign all constraints to rest */ - for (b = li_th->b0; b < li_th->b1; b++) + else { + /* Add the constraint to the rest block */ li_th->ind_r[li_th->nind_r++] = b; } } diff --git a/src/gromacs/mdlib/nbnxn_atomdata.c b/src/gromacs/mdlib/nbnxn_atomdata.c index ca0f4999bc..c20dfb2e54 100644 --- a/src/gromacs/mdlib/nbnxn_atomdata.c +++ b/src/gromacs/mdlib/nbnxn_atomdata.c @@ -1524,7 +1524,7 @@ static void nbnxn_atomdata_add_nbat_f_to_f_treereduce(const nbnxn_atomdata_t *nb i0 = b *NBNXN_BUFFERFLAG_SIZE*nbat->fstride; i1 = (b+1)*NBNXN_BUFFERFLAG_SIZE*nbat->fstride; - if ((flags->flag[b] & (1ULL< 2) + if (bitmask_is_set(flags->flag[b], index[1]) || group_size > 2) { #ifdef GMX_NBNXN_SIMD nbnxn_atomdata_reduce_reals_simd @@ -1532,11 +1532,11 @@ static void nbnxn_atomdata_add_nbat_f_to_f_treereduce(const nbnxn_atomdata_t *nb nbnxn_atomdata_reduce_reals #endif (nbat->out[index[0]].f, - (flags->flag[b] & (1ULL< 2, + bitmask_is_set(flags->flag[b], index[0]) || group_size > 2, &(nbat->out[index[1]].f), 1, i0, i1); } - else if (!(flags->flag[b] & (1ULL<flag[b], index[0])) { nbnxn_atomdata_clear_reals(nbat->out[index[0]].f, i0, i1); @@ -1576,7 +1576,7 @@ static void nbnxn_atomdata_add_nbat_f_to_f_stdreduce(const nbnxn_atomdata_t *nba nfptr = 0; for (out = 1; out < nbat->nout; out++) { - if (flags->flag[b] & (1U<flag[b], out)) { fptr[nfptr++] = nbat->out[out].f; } @@ -1589,11 +1589,11 @@ static void nbnxn_atomdata_add_nbat_f_to_f_stdreduce(const nbnxn_atomdata_t *nba nbnxn_atomdata_reduce_reals #endif (nbat->out[0].f, - flags->flag[b] & (1U<<0), + bitmask_is_set(flags->flag[b], 0), fptr, nfptr, i0, i1); } - else if (!(flags->flag[b] & (1U<<0))) + else if (!bitmask_is_set(flags->flag[b], 0)) { nbnxn_atomdata_clear_reals(nbat->out[0].f, i0, i1); diff --git a/src/gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_common.c b/src/gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_common.c index 2fded68209..9663861948 100644 --- a/src/gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_common.c +++ b/src/gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_common.c @@ -53,16 +53,16 @@ static void clear_f_flagged(const nbnxn_atomdata_t *nbat, int output_index, real *f) { const nbnxn_buffer_flags_t *flags; - unsigned our_flag; + gmx_bitmask_t our_flag; int g, b, a0, a1, i; flags = &nbat->buffer_flags; - our_flag = (1U << output_index); + bitmask_init_bit(&our_flag, output_index); for (b = 0; b < flags->nflag; b++) { - if (flags->flag[b] & our_flag) + if (!bitmask_is_disjoint(flags->flag[b], our_flag)) { a0 = b*NBNXN_BUFFERFLAG_SIZE; a1 = a0 + NBNXN_BUFFERFLAG_SIZE; diff --git a/src/gromacs/mdlib/nbnxn_pairlist.h b/src/gromacs/mdlib/nbnxn_pairlist.h index 413af35801..4558266bc3 100644 --- a/src/gromacs/mdlib/nbnxn_pairlist.h +++ b/src/gromacs/mdlib/nbnxn_pairlist.h @@ -42,6 +42,7 @@ #include "gromacs/legacyheaders/types/nblist.h" #include "gromacs/math/vectypes.h" +#include "gromacs/mdlib/bitmask.h" #include "gromacs/utility/basedefinitions.h" #include "gromacs/utility/real.h" @@ -207,21 +208,16 @@ typedef struct { #define NBNXN_BUFFERFLAG_SIZE 16 #endif -/* We currently store the reduction flags as bits in an unsigned int. - * In most cases this limits the number of flags to 32. - * The reduction will automatically disable the flagging and do a full - * reduction when the flags won't fit, but this will lead to very slow - * reduction. As we anyhow don't expect reasonable performance with - * more than 32 threads, we put in this hard limit. - * You can increase this number, but the reduction will be very slow. +/* We store the reduction flags as gmx_bitmask_t. + * This limits the number of flags to BITMASK_SIZE. */ -#define NBNXN_BUFFERFLAG_MAX_THREADS 32 +#define NBNXN_BUFFERFLAG_MAX_THREADS (BITMASK_SIZE) /* Flags for telling if threads write to force output buffers */ typedef struct { - int nflag; /* The number of flag blocks */ - unsigned int *flag; /* Bit i is set when thread i writes to a cell-block */ - int flag_nalloc; /* Allocation size of cxy_flag */ + int nflag; /* The number of flag blocks */ + gmx_bitmask_t *flag; /* Bit i is set when thread i writes to a cell-block */ + int flag_nalloc; /* Allocation size of cxy_flag */ } nbnxn_buffer_flags_t; /* LJ combination rules: geometric, Lorentz-Berthelot, none */ diff --git a/src/gromacs/mdlib/nbnxn_search.c b/src/gromacs/mdlib/nbnxn_search.c index 56e616967c..5c368b2ddf 100644 --- a/src/gromacs/mdlib/nbnxn_search.c +++ b/src/gromacs/mdlib/nbnxn_search.c @@ -1767,7 +1767,7 @@ static void init_buffer_flags(nbnxn_buffer_flags_t *flags, } for (b = 0; b < flags->nflag; b++) { - flags->flag[b] = 0; + bitmask_clear(&(flags->flag[b])); } } @@ -4846,7 +4846,7 @@ static void nbnxn_make_pairlist_part(const nbnxn_search_t nbs, int ndistc; int ncpcheck; int gridi_flag_shift = 0, gridj_flag_shift = 0; - unsigned int *gridj_flag = NULL; + gmx_bitmask_t *gridj_flag = NULL; int ncj_old_i, ncj_old_j; nbs_cycle_start(&work->cc[enbsCCsearch]); @@ -5352,7 +5352,7 @@ static void nbnxn_make_pairlist_part(const nbnxn_search_t nbs, cbl = nbl->cj[nbl->ncj-1].cj >> gridj_flag_shift; for (cb = cbf; cb <= cbl; cb++) { - gridj_flag[cb] = 1U<ncj > ncj_old_i) { - work->buffer_flags.flag[(gridi->cell0+ci)>>gridi_flag_shift] = 1U<buffer_flags.flag[(gridi->cell0+ci)>>gridi_flag_shift]), th); } } @@ -5451,8 +5451,8 @@ static void reduce_buffer_flags(const nbnxn_search_t nbs, int nsrc, const nbnxn_buffer_flags_t *dest) { - int s, b; - const unsigned int *flag; + int s, b; + gmx_bitmask_t *flag; for (s = 0; s < nsrc; s++) { @@ -5460,33 +5460,35 @@ static void reduce_buffer_flags(const nbnxn_search_t nbs, for (b = 0; b < dest->nflag; b++) { - dest->flag[b] |= flag[b]; + bitmask_union(&(dest->flag[b]), flag[b]); } } } static void print_reduction_cost(const nbnxn_buffer_flags_t *flags, int nout) { - int nelem, nkeep, ncopy, nred, b, c, out; + int nelem, nkeep, ncopy, nred, b, c, out; + gmx_bitmask_t mask_0; nelem = 0; nkeep = 0; ncopy = 0; nred = 0; + bitmask_init_bit(&mask_0, 0); for (b = 0; b < flags->nflag; b++) { - if (flags->flag[b] == 1) + if (bitmask_is_equal(flags->flag[b], mask_0)) { /* Only flag 0 is set, no copy of reduction required */ nelem++; nkeep++; } - else if (flags->flag[b] > 0) + else if (!bitmask_is_zero(flags->flag[b])) { c = 0; for (out = 0; out < nout; out++) { - if (flags->flag[b] & (1U<flag[b], out)) { c++; } diff --git a/src/gromacs/mdlib/tests/CMakeLists.txt b/src/gromacs/mdlib/tests/CMakeLists.txt index 73e80a0bd3..4e00c4717a 100644 --- a/src/gromacs/mdlib/tests/CMakeLists.txt +++ b/src/gromacs/mdlib/tests/CMakeLists.txt @@ -32,5 +32,6 @@ # To help us fund GROMACS development, we humbly ask that you cite # the research papers on the package. Check out http://www.gromacs.org. -gmx_add_unit_test(ShakeUnitTests shake-test +gmx_add_unit_test(MdlibUnitTest mdlib-test + bitmask32.cpp bitmask64.cpp bitmask128.cpp shake.cpp) diff --git a/src/gromacs/mdlib/tests/bitmask.h b/src/gromacs/mdlib/tests/bitmask.h new file mode 100644 index 0000000000..7041c770c2 --- /dev/null +++ b/src/gromacs/mdlib/tests/bitmask.h @@ -0,0 +1,132 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2014, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * 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. + */ +/*! \internal \file + * \brief + * Tests for bitmask functionality. + * + * These tests check the functionality of bitmask.h + + * \author Roland Schulz + */ +#include + +#include "gromacs/mdlib/bitmask.h" + +//! Implemenation of BITMASK_CLASSNAME +#define BITMASK_CLASSNAME_(S) BitmaskTest ## S +//! Returns name of Bitmask test fixture class +#define BITMASK_CLASSNAME(S) BITMASK_CLASSNAME_(S) +//! Implementation of BITMASK_TEST_P +#define BITMASK_TEST_P_(C, T) TEST_P(C, T) +//! Defines a parameterized bitmask test +#define BITMASK_TEST_P(T) BITMASK_TEST_P_(BITMASK_CLASSNAME(BITMASK_SIZE), T) + +class BITMASK_CLASSNAME(BITMASK_SIZE) : public ::testing::TestWithParam +{ +}; + +BITMASK_TEST_P(SetAndClear) +{ + gmx_bitmask_t m; + int i = GetParam(); + bitmask_clear(&m); + EXPECT_TRUE(bitmask_is_zero(m)); + EXPECT_FALSE(bitmask_is_set(m, i)); + bitmask_set_bit(&m, i); + for (int j = 0; j < BITMASK_SIZE; j++) + { + EXPECT_EQ(bitmask_is_set(m, j), j == i); + } + bitmask_clear(&m); + EXPECT_TRUE(bitmask_is_zero(m)); +} + +BITMASK_TEST_P(InitBit) +{ + gmx_bitmask_t m1, m2; + int i = GetParam(); + bitmask_init_bit(&m1, i); + bitmask_clear(&m2); + EXPECT_FALSE(bitmask_is_equal(m1, m2)); + bitmask_set_bit(&m2, i); + EXPECT_TRUE(bitmask_is_equal(m1, m2)); +} + +BITMASK_TEST_P(InitLowBits) +{ + gmx_bitmask_t m; + int i = GetParam(); + bitmask_init_low_bits(&m, i); + for (int j = 0; j < BITMASK_SIZE; j++) + { + EXPECT_EQ(bitmask_is_set(m, j), j < i); + } +} + +BITMASK_TEST_P(Disjoint) +{ + gmx_bitmask_t m1, m2; + int i = GetParam(); + bitmask_init_bit(&m1, i); + bitmask_init_bit(&m2, i); + EXPECT_FALSE(bitmask_is_disjoint(m1, m2)); + bitmask_init_low_bits(&m2, i); + EXPECT_TRUE(bitmask_is_disjoint(m1, m2)); +} + +BITMASK_TEST_P(Union) +{ + gmx_bitmask_t m1, m2; + int i = GetParam(); + int j = (i + BITMASK_SIZE/2)%BITMASK_SIZE; + bitmask_init_bit(&m1, i); + bitmask_init_bit(&m2, j); + bitmask_union(&m1, m2); + for (int k = 0; k < BITMASK_SIZE; k++) + { + EXPECT_EQ(bitmask_is_set(m1, k), k == i || k == j); + } + + bitmask_init_bit(&m1, i); + bitmask_clear(&m2); + bitmask_union(&m1, m2); + bitmask_init_bit(&m2, i); + EXPECT_TRUE(bitmask_is_equal(m1, m2)); + + bitmask_clear(&m1); + bitmask_init_bit(&m2, i); + bitmask_union(&m1, m2); + EXPECT_TRUE(bitmask_is_equal(m1, m2)); +} diff --git a/src/gromacs/mdlib/tests/bitmask128.cpp b/src/gromacs/mdlib/tests/bitmask128.cpp new file mode 100644 index 0000000000..b6024dbf8d --- /dev/null +++ b/src/gromacs/mdlib/tests/bitmask128.cpp @@ -0,0 +1,40 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2014, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * 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 "gmxpre.h" + +#define BITMASK_SIZE 128 +#include "bitmask.h" + +INSTANTIATE_TEST_CASE_P(BitmaskTest128_9_78, BitmaskTest128, ::testing::Values(9, 78)); diff --git a/src/gromacs/mdlib/tests/bitmask32.cpp b/src/gromacs/mdlib/tests/bitmask32.cpp new file mode 100644 index 0000000000..a2bddef98a --- /dev/null +++ b/src/gromacs/mdlib/tests/bitmask32.cpp @@ -0,0 +1,40 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2014, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * 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 "gmxpre.h" + +#define BITMASK_SIZE 32 +#include "bitmask.h" + +INSTANTIATE_TEST_CASE_P(BitmaskTest32_11, BitmaskTest32, ::testing::Values(11)); diff --git a/src/gromacs/mdlib/tests/bitmask64.cpp b/src/gromacs/mdlib/tests/bitmask64.cpp new file mode 100644 index 0000000000..1807efd42f --- /dev/null +++ b/src/gromacs/mdlib/tests/bitmask64.cpp @@ -0,0 +1,40 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2014, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * 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 "gmxpre.h" + +#define BITMASK_SIZE 64 +#include "bitmask.h" + +INSTANTIATE_TEST_CASE_P(BitmaskTest64_10_42, BitmaskTest64, ::testing::Values(10, 42)); -- 2.22.0