From 996c72871c65e4c2d002fa026eb5947a97357cd9 Mon Sep 17 00:00:00 2001
From: Gareth Tribello <gt@eider.phy.qub.ac.uk>
Date: Wed, 4 Nov 2015 09:41:25 +0000
Subject: [PATCH] A few optimisations in multicolvar and adjacency matrix
Also if we run DFS on an adjacency matrix without an input colvar
then the values used in DFS clustering are the sums of the rows
of the adjacency matrix
---
src/adjmat/ActionWithInputMatrix.cpp | 59 ++++++++++++++++++++++------
src/adjmat/ActionWithInputMatrix.h | 1 +
src/adjmat/AdjacencyMatrixBase.cpp | 7 +++-
src/adjmat/AdjacencyMatrixBase.h | 12 ++++++
src/adjmat/DFSBase.cpp | 4 --
src/adjmat/DFSBase.h | 2 -
src/adjmat/MatrixColumnSums.cpp | 16 ++------
src/adjmat/MatrixRowSums.cpp | 16 ++------
src/adjmat/MatrixSummationBase.cpp | 23 ++++++++++-
src/adjmat/MatrixSummationBase.h | 4 ++
src/multicolvar/AtomValuePack.cpp | 6 ++-
src/multicolvar/MultiColvarBase.h | 1 +
12 files changed, 107 insertions(+), 44 deletions(-)
diff --git a/src/adjmat/ActionWithInputMatrix.cpp b/src/adjmat/ActionWithInputMatrix.cpp
index 9c80388e3..b06919827 100644
--- a/src/adjmat/ActionWithInputMatrix.cpp
+++ b/src/adjmat/ActionWithInputMatrix.cpp
@@ -105,21 +105,52 @@ bool ActionWithInputMatrix::isCurrentlyActive( const unsigned& ind ) const {
}
void ActionWithInputMatrix::getVectorForTask( const unsigned& ind, const bool& normed, std::vector<double>& orient0 ) const {
- plumed_dbg_assert( isCurrentlyActive( ind ) );
- plumed_dbg_assert( ind<(mymatrix->function)->colvar_label.size() ); unsigned mmc=(mymatrix->function)->colvar_label[ind];
- plumed_dbg_assert( ((mymatrix->function)->mybasedata[mmc])->storedValueIsActive( (mymatrix->function)->convertToLocalIndex(ind,mmc) ) );
- ((mymatrix->function)->mybasedata[mmc])->retrieveValue( (mymatrix->function)->convertToLocalIndex(ind,mmc), normed, orient0 );
+ if( (mymatrix->function)->colvar_label.size()==0 ){
+ double df, sum=0.0; std::vector<double> tvals( mymatrix->getNumberOfComponents() );
+ unsigned vin; unsigned ncols = mymatrix->getNumberOfColumns(); orient0.assign(orient0.size(),0);
+ for(unsigned i=0;i<ncols;++i){
+ if( mymatrix->isSymmetric() && ind==i ) continue;
+ unsigned myelem = mymatrix->getStoreIndexFromMatrixIndices( ind, i );
+ if( !mymatrix->storedValueIsActive( myelem ) ) continue ;
+ mymatrix->retrieveValue( myelem, false, tvals );
+ orient0[1]+=(mymatrix->function)->transformStoredValues( tvals, vin, df);
+ }
+ orient0[0]=1.0;
+ } else {
+ plumed_dbg_assert( isCurrentlyActive( ind ) );
+ plumed_dbg_assert( ind<(mymatrix->function)->colvar_label.size() ); unsigned mmc=(mymatrix->function)->colvar_label[ind];
+ plumed_dbg_assert( ((mymatrix->function)->mybasedata[mmc])->storedValueIsActive( (mymatrix->function)->convertToLocalIndex(ind,mmc) ) );
+ ((mymatrix->function)->mybasedata[mmc])->retrieveValue( (mymatrix->function)->convertToLocalIndex(ind,mmc), normed, orient0 );
+ }
}
void ActionWithInputMatrix::getVectorDerivatives( const unsigned& ind, const bool& normed, MultiValue& myder ) const {
- plumed_dbg_assert( isCurrentlyActive( ind ) );
- plumed_dbg_assert( ind<(mymatrix->function)->colvar_label.size() ); unsigned mmc=(mymatrix->function)->colvar_label[ind];
- plumed_dbg_assert( ((mymatrix->function)->mybasedata[mmc])->storedValueIsActive( (mymatrix->function)->convertToLocalIndex(ind,mmc) ) );
- if( myder.getNumberOfValues()!=(mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfQuantities() ||
- myder.getNumberOfDerivatives()!=(mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfDerivatives() ){
- myder.resize( (mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfQuantities(), (mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfDerivatives() );
+ if( (mymatrix->function)->colvar_label.size()==0 ){
+ MultiValue myvals( 2, myder.getNumberOfDerivatives() );
+ double df, sum=0.0; std::vector<double> tvals( mymatrix->getNumberOfComponents() );
+ unsigned vin; unsigned ncols = mymatrix->getNumberOfColumns();
+ for(unsigned i=0;i<ncols;++i){
+ if( mymatrix->isSymmetric() && ind==i ) continue;
+ unsigned myelem = mymatrix->getStoreIndexFromMatrixIndices( ind, i );
+ if( !mymatrix->storedValueIsActive( myelem ) ) continue ;
+ mymatrix->retrieveValue( myelem, false, tvals );
+ double dum=(mymatrix->function)->transformStoredValues( tvals, vin, df);
+ mymatrix->retrieveDerivatives( myelem, false, myvals );
+ for(unsigned jd=0;jd<myvals.getNumberActive();++jd){
+ unsigned ider=myvals.getActiveIndex(jd);
+ myder.addDerivative( 1, ider, df*myvals.getDerivative( vin, ider ) );
+ }
+ }
+ } else {
+ plumed_dbg_assert( isCurrentlyActive( ind ) );
+ plumed_dbg_assert( ind<(mymatrix->function)->colvar_label.size() ); unsigned mmc=(mymatrix->function)->colvar_label[ind];
+ plumed_dbg_assert( ((mymatrix->function)->mybasedata[mmc])->storedValueIsActive( (mymatrix->function)->convertToLocalIndex(ind,mmc) ) );
+ if( myder.getNumberOfValues()!=(mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfQuantities() ||
+ myder.getNumberOfDerivatives()!=(mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfDerivatives() ){
+ myder.resize( (mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfQuantities(), (mymatrix->function)->mybasemulticolvars[mmc]->getNumberOfDerivatives() );
+ }
+ (mymatrix->function)->mybasedata[mmc]->retrieveDerivatives( (mymatrix->function)->convertToLocalIndex(ind,mmc), normed, myder );
}
- (mymatrix->function)->mybasedata[mmc]->retrieveDerivatives( (mymatrix->function)->convertToLocalIndex(ind,mmc), normed, myder );
}
Vector ActionWithInputMatrix::getSeparation( const Vector& vec1, const Vector& vec2 ) const {
@@ -133,12 +164,18 @@ unsigned ActionWithInputMatrix::getNumberOfNodeTypes() const {
return size;
}
+unsigned ActionWithInputMatrix::getNumberOfQuantities(){
+ if( (mymatrix->function)->colvar_label.size()==0 ) return 2;
+ else return (mymatrix->function)->mybasemulticolvars[0]->getNumberOfQuantities();
+}
+
unsigned ActionWithInputMatrix::getNumberOfAtomsInGroup( const unsigned& igrp ) const {
plumed_dbg_assert( igrp<(mymatrix->function)->mybasemulticolvars.size() );
return (mymatrix->function)->mybasemulticolvars[igrp]->getFullNumberOfTasks();
}
multicolvar::MultiColvarBase* ActionWithInputMatrix::getBaseMultiColvar( const unsigned& igrp ) const {
+ plumed_dbg_assert( igrp<(mymatrix->function)->mybasemulticolvars.size() );
return (mymatrix->function)->mybasemulticolvars[igrp];
}
diff --git a/src/adjmat/ActionWithInputMatrix.h b/src/adjmat/ActionWithInputMatrix.h
index b5770d8b0..f1c874568 100644
--- a/src/adjmat/ActionWithInputMatrix.h
+++ b/src/adjmat/ActionWithInputMatrix.h
@@ -72,6 +72,7 @@ public:
// Turn on the derivatives
virtual void turnOnDerivatives();
bool isPeriodic(){ return false; }
+ unsigned getNumberOfQuantities();
void apply();
///
AtomNumber getAbsoluteIndexOfCentralAtom(const unsigned& i) const ;
diff --git a/src/adjmat/AdjacencyMatrixBase.cpp b/src/adjmat/AdjacencyMatrixBase.cpp
index 1fc3357cd..f07d5bb0c 100644
--- a/src/adjmat/AdjacencyMatrixBase.cpp
+++ b/src/adjmat/AdjacencyMatrixBase.cpp
@@ -53,7 +53,12 @@ bool AdjacencyMatrixBase::parseAtomList(const std::string& key, const int& num,
if( mlabs.size()==0 ) return false;
bool found_acts=interpretInputMultiColvars(mlabs,0.0);
- if( !found_acts ) ActionAtomistic::interpretAtomList( mlabs, t );
+ if( !found_acts ){
+ ActionAtomistic::interpretAtomList( mlabs, t );
+ log.printf(" involving atoms ");
+ for(unsigned i=0;i<t.size();++i) log.printf("%d ",t[i].serial() );
+ log.printf("\n");
+ }
return true;
}
diff --git a/src/adjmat/AdjacencyMatrixBase.h b/src/adjmat/AdjacencyMatrixBase.h
index 3b269b60d..c11ee64df 100644
--- a/src/adjmat/AdjacencyMatrixBase.h
+++ b/src/adjmat/AdjacencyMatrixBase.h
@@ -82,6 +82,8 @@ public:
virtual double transformStoredValues( const std::vector<double>& myvals, unsigned& vout, double& df ) const ;
/// Used to check for connections between atoms
virtual bool checkForConnection( const std::vector<double>& myvals ) const=0;
+/// Get the atom number
+ AtomNumber getAbsoluteIndexOfCentralAtom( const unsigned& i ) const ;
};
inline
@@ -95,6 +97,16 @@ unsigned AdjacencyMatrixBase::getBaseColvarNumber( const unsigned& inum ) const
return 0;
}
+inline
+AtomNumber AdjacencyMatrixBase::getAbsoluteIndexOfCentralAtom( const unsigned& iatom ) const {
+ if( iatom<colvar_label.size() ){
+ unsigned mmc=colvar_label[ iatom ];
+ return mybasemulticolvars[mmc]->getAbsoluteIndexOfCentralAtom( convertToLocalIndex(iatom,mmc) );
+ }
+ return ActionAtomistic::getAbsoluteIndex( iatom );
+}
+
+
inline
void AdjacencyMatrixBase::getOrientationVector( const unsigned& ind, const bool& normed, std::vector<double>& orient ) const {
plumed_dbg_assert( ind<colvar_label.size() ); unsigned mmc=colvar_label[ind];
diff --git a/src/adjmat/DFSBase.cpp b/src/adjmat/DFSBase.cpp
index df91c3a3d..2578bcee6 100644
--- a/src/adjmat/DFSBase.cpp
+++ b/src/adjmat/DFSBase.cpp
@@ -68,10 +68,6 @@ void DFSBase::turnOnDerivatives(){
ActionWithInputMatrix::turnOnDerivatives();
}
-unsigned DFSBase::getNumberOfQuantities(){
- return getBaseMultiColvar(0)->getNumberOfQuantities();
-}
-
void DFSBase::performClustering(){
// All the clusters have zero size initially
for(unsigned i=0;i<cluster_sizes.size();++i){ cluster_sizes[i].first=0; cluster_sizes[i].second=i; }
diff --git a/src/adjmat/DFSBase.h b/src/adjmat/DFSBase.h
index b3b776061..1a968ab86 100644
--- a/src/adjmat/DFSBase.h
+++ b/src/adjmat/DFSBase.h
@@ -62,8 +62,6 @@ public:
static void registerKeywords( Keywords& keys );
/// Constructor
explicit DFSBase(const ActionOptions&);
-/// Required as we have to be able to deal with vectors
- unsigned getNumberOfQuantities();
/// This checks whether derivatives can be computed given the base multicolvar
void turnOnDerivatives();
};
diff --git a/src/adjmat/MatrixColumnSums.cpp b/src/adjmat/MatrixColumnSums.cpp
index 7f119baff..a84bb300b 100644
--- a/src/adjmat/MatrixColumnSums.cpp
+++ b/src/adjmat/MatrixColumnSums.cpp
@@ -73,27 +73,19 @@ MatrixSummationBase(ao)
}
double MatrixColumnSums::compute( const unsigned& tinded, multicolvar::AtomValuePack& myatoms ) const {
- double sum=0.0; std::vector<double> tvals(2);
+ double sum=0.0; std::vector<double> tvals( mymatrix->getNumberOfComponents() );
unsigned nrows = mymatrix->getNumberOfRows();
for(unsigned i=0;i<nrows;++i){
if( mymatrix->isSymmetric() && tinded==i ) continue;
- unsigned myelem = mymatrix->getStoreIndexFromMatrixIndices( i, tinded );
- mymatrix->retrieveValue( myelem, false, tvals );
- sum+=tvals[1];
+ sum+=retrieveConnectionValue( i, tinded, tvals );
}
if( !doNotCalculateDerivatives() ){
- MultiValue myvals( 2, myatoms.getNumberOfDerivatives() );
+ MultiValue myvals( mymatrix->getNumberOfComponents(), myatoms.getNumberOfDerivatives() );
MultiValue& myvout=myatoms.getUnderlyingMultiValue();
for(unsigned i=0;i<nrows;++i){
if( mymatrix->isSymmetric() && tinded==i ) continue ;
- unsigned myelem = mymatrix->getStoreIndexFromMatrixIndices( i, tinded );
- if( !mymatrix->storedValueIsActive( myelem ) ) continue ;
- mymatrix->retrieveDerivatives( myelem, false, myvals );
- for(unsigned jd=0;jd<myvals.getNumberActive();++jd){
- unsigned ider=myvals.getActiveIndex(jd);
- myvout.addDerivative( 1, ider, myvals.getDerivative( 1, ider ) );
- }
+ addConnectionDerivatives( i, tinded, tvals, myvals, myvout );
}
}
return sum;
diff --git a/src/adjmat/MatrixRowSums.cpp b/src/adjmat/MatrixRowSums.cpp
index 7807c2002..356649bcd 100644
--- a/src/adjmat/MatrixRowSums.cpp
+++ b/src/adjmat/MatrixRowSums.cpp
@@ -67,27 +67,19 @@ MatrixSummationBase(ao)
}
double MatrixRowSums::compute( const unsigned& tinded, multicolvar::AtomValuePack& myatoms ) const {
- double sum=0.0; std::vector<double> tvals(2);
+ double sum=0.0; std::vector<double> tvals( mymatrix->getNumberOfComponents() );
unsigned ncols = mymatrix->getNumberOfColumns();
for(unsigned i=0;i<ncols;++i){
if( mymatrix->isSymmetric() && tinded==i ) continue;
- unsigned myelem = mymatrix->getStoreIndexFromMatrixIndices( tinded, i );
- mymatrix->retrieveValue( myelem, false, tvals );
- sum+=tvals[1];
+ sum+=retrieveConnectionValue( tinded, i, tvals );
}
if( !doNotCalculateDerivatives() ){
- MultiValue myvals( 2, myatoms.getNumberOfDerivatives() );
+ MultiValue myvals( mymatrix->getNumberOfComponents(), myatoms.getNumberOfDerivatives() );
MultiValue& myvout=myatoms.getUnderlyingMultiValue();
for(unsigned i=0;i<ncols;++i){
if( mymatrix->isSymmetric() && tinded==i ) continue ;
- unsigned myelem = mymatrix->getStoreIndexFromMatrixIndices( tinded, i );
- if( !mymatrix->storedValueIsActive( myelem ) ) continue ;
- mymatrix->retrieveDerivatives( myelem, false, myvals );
- for(unsigned jd=0;jd<myvals.getNumberActive();++jd){
- unsigned ider=myvals.getActiveIndex(jd);
- myvout.addDerivative( 1, ider, myvals.getDerivative( 1, ider ) );
- }
+ addConnectionDerivatives( tinded, i, tvals, myvals, myvout );
}
}
return sum;
diff --git a/src/adjmat/MatrixSummationBase.cpp b/src/adjmat/MatrixSummationBase.cpp
index 49d704cc5..40d47e27b 100644
--- a/src/adjmat/MatrixSummationBase.cpp
+++ b/src/adjmat/MatrixSummationBase.cpp
@@ -72,7 +72,28 @@ Vector MatrixSummationBase::getPositionOfAtomForLinkCells( const unsigned& iatom
AtomNumber MatrixSummationBase::getAbsoluteIndexOfCentralAtom(const unsigned& i) const {
return (mymatrix->function)->getAbsoluteIndexOfCentralAtom(i);
-}
+}
+
+double MatrixSummationBase::retrieveConnectionValue( const unsigned& i, const unsigned& j, std::vector<double>& vals ) const {
+ unsigned vi, myelem = mymatrix->getStoreIndexFromMatrixIndices( i, j );
+ if( !mymatrix->storedValueIsActive( myelem ) ) return 0;
+
+ mymatrix->retrieveValue( myelem, false, vals ); double df;
+ return (mymatrix->function)->transformStoredValues( vals, vi, df );
+}
+
+void MatrixSummationBase::addConnectionDerivatives( const unsigned& i, const unsigned& j, std::vector<double>& vals, MultiValue& myvals, MultiValue& myvout ) const {
+ unsigned vi, myelem = mymatrix->getStoreIndexFromMatrixIndices( i, j );
+ if( !mymatrix->storedValueIsActive( myelem ) ) return;
+
+ mymatrix->retrieveValue( myelem, false, vals ); double df;
+ double vv = (mymatrix->function)->transformStoredValues( vals, vi, df );
+ mymatrix->retrieveDerivatives( myelem, false, myvals );
+ for(unsigned jd=0;jd<myvals.getNumberActive();++jd){
+ unsigned ider=myvals.getActiveIndex(jd);
+ myvout.addDerivative( 1, ider, df*myvals.getDerivative( vi, ider ) );
+ }
+}
}
}
diff --git a/src/adjmat/MatrixSummationBase.h b/src/adjmat/MatrixSummationBase.h
index bba8d5698..a22c58b68 100644
--- a/src/adjmat/MatrixSummationBase.h
+++ b/src/adjmat/MatrixSummationBase.h
@@ -33,6 +33,10 @@ friend class ActionWithInputMatrix;
protected:
/// The vessel that holds the adjacency matrix
AdjacencyMatrixVessel* mymatrix;
+/// Get the value of a connection
+ double retrieveConnectionValue( const unsigned& i, const unsigned& j, std::vector<double>& vals ) const ;
+/// Add the derivatives
+ void addConnectionDerivatives( const unsigned& i, const unsigned& j, std::vector<double>& vals, MultiValue& myvals, MultiValue& myvout ) const ;
public:
static void registerKeywords( Keywords& keys );
explicit MatrixSummationBase(const ActionOptions&);
diff --git a/src/multicolvar/AtomValuePack.cpp b/src/multicolvar/AtomValuePack.cpp
index 66f7e4cbf..da752f5aa 100644
--- a/src/multicolvar/AtomValuePack.cpp
+++ b/src/multicolvar/AtomValuePack.cpp
@@ -44,7 +44,11 @@ myatoms( vals.getAtomVector() )
unsigned AtomValuePack::setupAtomsFromLinkCells( const std::vector<unsigned>& cind, const Vector& cpos, const LinkCells& linkcells ){
natoms=cind.size(); for(unsigned i=0;i<natoms;++i) indices[i]=cind[i];
linkcells.retrieveNeighboringAtoms( cpos, natoms, indices );
- for(unsigned i=0;i<natoms;++i) myatoms[i]=mycolv->getPositionOfAtomForLinkCells( indices[i] ) - cpos;
+ if( mycolv->colvar_label.size()==0 ){
+ for(unsigned i=0;i<natoms;++i) myatoms[i]=mycolv->getPosition( indices[i] ) - cpos;
+ } else {
+ for(unsigned i=0;i<natoms;++i) myatoms[i]=mycolv->getPositionOfAtomForLinkCells( indices[i] ) - cpos;
+ }
if( mycolv->usesPbc() ) mycolv->applyPbc( myatoms, natoms );
return natoms;
}
diff --git a/src/multicolvar/MultiColvarBase.h b/src/multicolvar/MultiColvarBase.h
index b3eaed4bb..ff0852732 100644
--- a/src/multicolvar/MultiColvarBase.h
+++ b/src/multicolvar/MultiColvarBase.h
@@ -45,6 +45,7 @@ class MultiColvarBase :
friend class BridgedMultiColvarFunction;
friend class VolumeGradientBase;
friend class MultiColvarFilter;
+ friend class AtomValuePack;
private:
/// Use periodic boundary conditions
bool usepbc;
--
GitLab