From 29cf9aea11893398a5a9d691ca299cbc353cd357 Mon Sep 17 00:00:00 2001 From: Berk Hess Date: Mon, 24 May 2021 07:44:50 +0000 Subject: [PATCH] Fix incorrectly sized ga2la hash table --- docs/release-notes/2021/2021.3.rst | 8 ++++++++ src/gromacs/domdec/domdec_constraints.cpp | 4 ++-- src/gromacs/domdec/domdec_vsite.cpp | 4 ++-- src/gromacs/domdec/ga2la.h | 13 ++++++++++--- src/gromacs/domdec/hashedmap.h | 22 ++++++++++++++-------- src/gromacs/domdec/partition.cpp | 4 ++-- src/gromacs/domdec/tests/hashedmap.cpp | 8 ++++---- src/gromacs/mdlib/md_support.cpp | 2 +- src/gromacs/mdlib/updategroupscog.cpp | 4 ++-- 9 files changed, 45 insertions(+), 24 deletions(-) diff --git a/docs/release-notes/2021/2021.3.rst b/docs/release-notes/2021/2021.3.rst index b324ac639d..54cce59905 100644 --- a/docs/release-notes/2021/2021.3.rst +++ b/docs/release-notes/2021/2021.3.rst @@ -32,3 +32,11 @@ The source code validation could otherwise fail a build with cryptic errors. Miscellaneous ^^^^^^^^^^^^^ +Removed performance loss in the mdrun domain decomposition +"""""""""""""""""""""""""""""""""""""""""""""""""""""""""" + +With 16 or more so-called PP MPI ranks, the domain decomposition +repartitioning could incur large performance overheads due to a sub-optimally +sized hash table. This has now been fixed. + +:issue:`4054` diff --git a/src/gromacs/domdec/domdec_constraints.cpp b/src/gromacs/domdec/domdec_constraints.cpp index b11ef6bdb3..501c0b2c80 100644 --- a/src/gromacs/domdec/domdec_constraints.cpp +++ b/src/gromacs/domdec/domdec_constraints.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 2006,2007,2008,2009,2010 by the GROMACS development team. * Copyright (c) 2012,2013,2014,2015,2016 by the GROMACS development team. - * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2017,2018,2019,2020,2021, 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. @@ -138,7 +138,7 @@ void dd_clear_local_constraint_indices(gmx_domdec_t* dd) if (dd->constraint_comm) { - dc->ga2la->clear(); + dc->ga2la->clearAndResizeHashTable(); } } diff --git a/src/gromacs/domdec/domdec_vsite.cpp b/src/gromacs/domdec/domdec_vsite.cpp index 0fe374ee45..957249bc72 100644 --- a/src/gromacs/domdec/domdec_vsite.cpp +++ b/src/gromacs/domdec/domdec_vsite.cpp @@ -3,7 +3,7 @@ * * Copyright (c) 2006,2007,2008,2009,2010 by the GROMACS development team. * Copyright (c) 2012,2013,2014,2015,2016 by the GROMACS development team. - * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2017,2018,2019,2020,2021, 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. @@ -97,7 +97,7 @@ void dd_clear_local_vsite_indices(gmx_domdec_t* dd) { if (dd->vsite_comm) { - dd->ga2la_vsite->clear(); + dd->ga2la_vsite->clearAndResizeHashTable(); } } diff --git a/src/gromacs/domdec/ga2la.h b/src/gromacs/domdec/ga2la.h index e4b972afbf..6c054ebd19 100644 --- a/src/gromacs/domdec/ga2la.h +++ b/src/gromacs/domdec/ga2la.h @@ -4,7 +4,7 @@ * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. * Copyright (c) 2010,2014,2015,2017,2018 by the GROMACS development team. - * Copyright (c) 2019,2020, by the GROMACS development team, led by + * Copyright (c) 2019,2020,2021, 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. @@ -149,8 +149,11 @@ public: } } - //! Clear all the entries in the list. - void clear() + /*! \brief Clear all the entries in the list. + * + * \param[in] resizeHashTable When true the hash table is optimized based on the current number of entries stored + */ + void clear(const bool resizeHashTable) { if (usingDirect_) { @@ -159,6 +162,10 @@ public: entry.cell = -1; } } + else if (resizeHashTable) + { + data_.hashed.clearAndResizeHashTable(); + } else { data_.hashed.clear(); diff --git a/src/gromacs/domdec/hashedmap.h b/src/gromacs/domdec/hashedmap.h index 00f37fc81a..6877b7d82f 100644 --- a/src/gromacs/domdec/hashedmap.h +++ b/src/gromacs/domdec/hashedmap.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2020,2021, 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. @@ -285,15 +285,9 @@ public: return nullptr; } - /*! \brief Clear all the entries in the list - * - * Also optimizes the size of the table based on the current - * number of elements stored. - */ + //! Clear all the entries in the list void clear() { - const int oldNumElements = numElements_; - for (hashEntry& entry : table_) { entry.key = -1; @@ -301,6 +295,18 @@ public: } startIndexForSpaceForListEntry_ = bucket_count(); numElements_ = 0; + } + + /*! \brief Clear all the entries in the list and resizes the hash table + * + * Optimizes the size of the hash table based on the current + * number of elements stored. + */ + void clearAndResizeHashTable() + { + const int oldNumElements = numElements_; + + clear(); /* Resize the hash table when the occupation is far from optimal. * Do not resize with 0 elements to avoid minimal size when clear() diff --git a/src/gromacs/domdec/partition.cpp b/src/gromacs/domdec/partition.cpp index 67a79d8e85..44453ac334 100644 --- a/src/gromacs/domdec/partition.cpp +++ b/src/gromacs/domdec/partition.cpp @@ -606,7 +606,7 @@ static void clearDDStateIndices(gmx_domdec_t* dd, const bool keepLocalAtomIndice if (!keepLocalAtomIndices) { /* Clear the whole list without the overhead of searching */ - ga2la.clear(); + ga2la.clear(true); } else { @@ -3019,7 +3019,7 @@ void dd_partition_system(FILE* fplog, state_change_natoms(state_local, comm->atomRanges.numHomeAtoms()); /* Rebuild all the indices */ - dd->ga2la->clear(); + dd->ga2la->clear(false); ncgindex_set = 0; wallcycle_sub_stop(wcycle, ewcsDD_GRID); diff --git a/src/gromacs/domdec/tests/hashedmap.cpp b/src/gromacs/domdec/tests/hashedmap.cpp index 7157633e1d..f13a17618e 100644 --- a/src/gromacs/domdec/tests/hashedmap.cpp +++ b/src/gromacs/domdec/tests/hashedmap.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2018,2019, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2021, 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. @@ -192,18 +192,18 @@ TEST(HashedMap, ResizesTable) EXPECT_LT(map.bucket_count(), 128); // Check that the table size is double #elements after clear() - map.clear(); + map.clearAndResizeHashTable(); EXPECT_EQ(map.bucket_count(), 128); // Check that calling clear() a second time does not resize - map.clear(); + map.clearAndResizeHashTable(); EXPECT_EQ(map.bucket_count(), 128); map.insert(2, 'b'); EXPECT_EQ(map.bucket_count(), 128); // Check that calling clear with 1 elements sizes down - map.clear(); + map.clearAndResizeHashTable(); EXPECT_LT(map.bucket_count(), 128); } diff --git a/src/gromacs/mdlib/md_support.cpp b/src/gromacs/mdlib/md_support.cpp index cecbc3ca0c..5b72136114 100644 --- a/src/gromacs/mdlib/md_support.cpp +++ b/src/gromacs/mdlib/md_support.cpp @@ -4,7 +4,7 @@ * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. * Copyright (c) 2013,2014,2015,2016,2017 The GROMACS development team. - * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2020,2021, 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. diff --git a/src/gromacs/mdlib/updategroupscog.cpp b/src/gromacs/mdlib/updategroupscog.cpp index 2739a398bc..a929f5e92f 100644 --- a/src/gromacs/mdlib/updategroupscog.cpp +++ b/src/gromacs/mdlib/updategroupscog.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2018,2019, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2021, 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. @@ -138,7 +138,7 @@ void UpdateGroupsCog::clear() cogIndices_.clear(); cogs_.clear(); numAtomsPerCog_.clear(); - globalToLocalMap_.clear(); + globalToLocalMap_.clearAndResizeHashTable(); } } // namespace gmx -- 2.22.0