From d15c12f00d2b98f103e6bb43bc0035de43795095 Mon Sep 17 00:00:00 2001
From: Gareth Aneurin Tribello <garethtribello@Gareths-MacBook-Pro.local>
Date: Wed, 10 Jul 2024 10:43:59 +0100
Subject: [PATCH] Optimisations in adjacency matrix calculator

---
 src/adjmat/AdjacencyMatrixBase.cpp | 15 +++++----------
 src/adjmat/AdjacencyMatrixBase.h   |  6 ++++++
 2 files changed, 11 insertions(+), 10 deletions(-)

diff --git a/src/adjmat/AdjacencyMatrixBase.cpp b/src/adjmat/AdjacencyMatrixBase.cpp
index 7f6fad4f75..8363d04965 100644
--- a/src/adjmat/AdjacencyMatrixBase.cpp
+++ b/src/adjmat/AdjacencyMatrixBase.cpp
@@ -228,6 +228,11 @@ void AdjacencyMatrixBase::updateNeighbourList() {
     }
     threecells.buildCellLists( ltmp_pos2, threeblocks, getPbc() );
   }
+  // And finally work out maximum number of columns to use
+  maxcol = nlist[0];
+  for(unsigned i=1; i<getConstPntrToComponent(0)->getShape()[0]; ++i) {
+    if( nlist[i]>maxcol ) maxcol = nlist[i];
+  }
 }
 
 void AdjacencyMatrixBase::getAdditionalTasksRequired( ActionWithVector* action, std::vector<unsigned>& atasks ) {
@@ -249,16 +254,6 @@ void AdjacencyMatrixBase::getAdditionalTasksRequired( ActionWithVector* action,
   }
 }
 
-unsigned AdjacencyMatrixBase::getNumberOfColumns() const {
-  unsigned maxcol=nlist[0];
-  for(unsigned i=1; i<getConstPntrToComponent(0)->getShape()[0]; ++i) {
-    if( nlist[i]>maxcol ) maxcol = nlist[i];
-  }
-  // This triggers appropriate storage of symmetric matrices which can only have
-  // a maximumum of shape[1]-1 non-zero elements
-  return maxcol;
-}
-
 unsigned AdjacencyMatrixBase::retrieveNeighbours( const unsigned& current, std::vector<unsigned> & indices ) const {
   unsigned natoms=nlist[current]; indices[0]=current;
   unsigned lstart = getConstPntrToComponent(0)->getShape()[0] + current*(1+natoms_per_list); plumed_dbg_assert( nlist[lstart]==current );
diff --git a/src/adjmat/AdjacencyMatrixBase.h b/src/adjmat/AdjacencyMatrixBase.h
index 3d500cf9c1..e529688d8d 100644
--- a/src/adjmat/AdjacencyMatrixBase.h
+++ b/src/adjmat/AdjacencyMatrixBase.h
@@ -36,6 +36,7 @@ class AdjacencyMatrixBase : public ActionWithMatrix {
   LinkCells linkcells, threecells;
   std::vector<unsigned> ablocks, threeblocks;
   double nl_cut, nl_cut2;
+  unsigned maxcol;
   unsigned nl_stride;
   unsigned natoms_per_list;
   std::vector<unsigned> nlist;
@@ -106,6 +107,11 @@ void AdjacencyMatrixBase::addBoxDerivatives( const Tensor& vir, MultiValue& myva
   myvals.addDerivative( w_index, nbase+8, vir(2,2) );
 }
 
+inline
+unsigned AdjacencyMatrixBase::getNumberOfColumns() const {
+  return maxcol;
+}
+
 
 }
 }