diff --git a/src/adjmat/AdjacencyMatrixBase.cpp b/src/adjmat/AdjacencyMatrixBase.cpp
index 808725c6c796b4a2eef04d6fcc5074332e58f077..5af5c6751762e7960098fe6c8dbf46c03a562b4a 100644
--- a/src/adjmat/AdjacencyMatrixBase.cpp
+++ b/src/adjmat/AdjacencyMatrixBase.cpp
@@ -204,19 +204,6 @@ void AdjacencyMatrixBase::requestAtoms( const std::vector<AtomNumber>& atoms, co
setupMultiColvarBase( atoms );
}
-void AdjacencyMatrixBase::addAtomDerivatives( const unsigned& ival, const unsigned& iatom, const Vector& der, multicolvar::AtomValuePack& myatoms ) const {
- unsigned jatom=myatoms.getIndex(iatom);
-
- if( jatom>colvar_label.size() ){
- myatoms.addAtomsDerivatives( ival, jatom, der );
- } else {
- unsigned mmc=colvar_label[jatom];
- unsigned basen=0; for(unsigned i=0;i<mmc;++i) basen+=mybasemulticolvars[i]->getNumberOfAtoms();
- multicolvar::CatomPack atom0=mybasemulticolvars[mmc]->getCentralAtomPack( basen, convertToLocalIndex(jatom,mmc) );
- myatoms.addComDerivatives( ival, der, atom0 );
- }
-}
-
// Maybe put this back GAT to check that it is returning an atom number that is one of the nodes
// and not a hydrogen if we are doing HBPAMM
// AtomNumber AdjacencyMatrixBase::getAbsoluteIndexOfCentralAtom(const unsigned& i) const {
diff --git a/src/adjmat/AdjacencyMatrixBase.h b/src/adjmat/AdjacencyMatrixBase.h
index 8db15e752cdb7accd5501692462880b022f1eafd..4ff1a1d44a2936fa73541599258ef2319adcea53 100644
--- a/src/adjmat/AdjacencyMatrixBase.h
+++ b/src/adjmat/AdjacencyMatrixBase.h
@@ -66,8 +66,6 @@ protected:
void requestAtoms( const std::vector<AtomNumber>& atoms, const bool& symmetric, const unsigned& nrows );
/// Return the group this atom is a part of
unsigned getBaseColvarNumber( const unsigned& ) const ;
-/// Add some derivatives to a particular component of a particular atom
- void addAtomDerivatives( const unsigned& , const unsigned& , const Vector& , multicolvar::AtomValuePack& ) const ;
/// Add some derivatives to the relevant orientation
void addOrientationDerivatives( const unsigned&, const unsigned& , const std::vector<double>& , multicolvar::AtomValuePack& ) const ;
public:
diff --git a/src/crystallization/Fccubic.cpp b/src/crystallization/Fccubic.cpp
index 617ae7070bda514c1be1907ac83b991e9c036acd..47b95a6dadcc0d5ff1a6d0267dfd2b9cf2f58b1e 100644
--- a/src/crystallization/Fccubic.cpp
+++ b/src/crystallization/Fccubic.cpp
@@ -197,11 +197,11 @@ double Fccubic::compute( const unsigned& tindex, multicolvar::AtomValuePack& mya
fder = (+dfunc)*tmp*distance + sw*myder;
- myatoms.addAtomsDerivatives( 1, 0, -fder );
- myatoms.addAtomsDerivatives( 1, i, +fder );
+ addAtomDerivatives( 1, 0, -fder, myatoms );
+ addAtomDerivatives( 1, i, +fder, myatoms);
myatoms.addBoxDerivatives( 1, Tensor(distance,-fder) );
- myatoms.addAtomsDerivatives( 0, 0, (-dfunc)*distance );
- myatoms.addAtomsDerivatives( 0, i, (+dfunc)*distance );
+ addAtomDerivatives( 0, 0, (-dfunc)*distance, myatoms);
+ addAtomDerivatives( 0, i, (+dfunc)*distance, myatoms);
myatoms.addBoxDerivatives( 0, (-dfunc)*Tensor(distance,distance) );
}
}
diff --git a/src/crystallization/MoleculeOrientation.cpp b/src/crystallization/MoleculeOrientation.cpp
index ef892889762ebda83dd02fca3cb2fc9ccbb77886..4a7879ce4e2b61a88d72c7b9d8d3f17e904da641 100644
--- a/src/crystallization/MoleculeOrientation.cpp
+++ b/src/crystallization/MoleculeOrientation.cpp
@@ -89,18 +89,18 @@ VectorMultiColvar(ao)
void MoleculeOrientation::calculateVector( multicolvar::AtomValuePack& myatoms ) const {
Vector distance; distance=getSeparation( myatoms.getPosition(0), myatoms.getPosition(1) );
- myatoms.addAtomsDerivatives( 2, 0, Vector(-1.0,0,0) );
- myatoms.addAtomsDerivatives( 2, 1, Vector(+1.0,0,0) );
+ addAtomDerivatives( 2, 0, Vector(-1.0,0,0), myatoms );
+ addAtomDerivatives( 2, 1, Vector(+1.0,0,0), myatoms );
myatoms.addBoxDerivatives( 2, Tensor(distance,Vector(-1.0,0,0)) );
myatoms.addValue( 2, distance[0] );
- myatoms.addAtomsDerivatives( 3, 0, Vector(0,-1.0,0) );
- myatoms.addAtomsDerivatives( 3, 1, Vector(0,+1.0,0) );
+ addAtomDerivatives( 3, 0, Vector(0,-1.0,0), myatoms );
+ addAtomDerivatives( 3, 1, Vector(0,+1.0,0), myatoms );
myatoms.addBoxDerivatives( 3, Tensor(distance,Vector(0,-1.0,0)) );
myatoms.addValue( 3, distance[1] );
- myatoms.addAtomsDerivatives( 4, 0, Vector(0,0,-1.0) );
- myatoms.addAtomsDerivatives( 4, 1, Vector(0,0,+1.0) );
+ addAtomDerivatives( 4, 0, Vector(0,0,-1.0), myatoms );
+ addAtomDerivatives( 4, 1, Vector(0,0,+1.0), myatoms );
myatoms.addBoxDerivatives( 4, Tensor(distance,Vector(0,0,-1.0)) );
myatoms.addValue( 4, distance[2] );
}
diff --git a/src/crystallization/MoleculePlane.cpp b/src/crystallization/MoleculePlane.cpp
index ae259d325cb2763bce9a0ebb210b7ffc6db4d54b..cb4f953dcc25c2dfe467bbc1257c02d050e1f0c4 100644
--- a/src/crystallization/MoleculePlane.cpp
+++ b/src/crystallization/MoleculePlane.cpp
@@ -81,38 +81,38 @@ void MoleculePlane::calculateVector( multicolvar::AtomValuePack& myatoms ) const
}
cp = crossProduct( d1, d2 );
- myatoms.addAtomsDerivatives( 2, 0, crossProduct( Vector(-1.0,0,0), d2 ) );
+ addAtomDerivatives( 2, 0, crossProduct( Vector(-1.0,0,0), d2 ), myatoms );
if( myatoms.getNumberOfAtoms()==3 ){
- myatoms.addAtomsDerivatives( 2, 1, crossProduct( Vector(+1.0,0,0), d2 ) + crossProduct( Vector(-1.0,0,0), d1 ) );
- myatoms.addAtomsDerivatives( 2, 2, crossProduct( Vector(+1.0,0,0), d1 ) );
+ addAtomDerivatives( 2, 1, crossProduct( Vector(+1.0,0,0), d2 ) + crossProduct( Vector(-1.0,0,0), d1 ), myatoms );
+ addAtomDerivatives( 2, 2, crossProduct( Vector(+1.0,0,0), d1 ), myatoms );
} else {
- myatoms.addAtomsDerivatives( 2, 1, crossProduct( Vector(+1.0,0,0), d2 ) );
- myatoms.addAtomsDerivatives( 2, 2, crossProduct( Vector(-1.0,0,0), d1 ) );
- myatoms.addAtomsDerivatives( 2, 3, crossProduct( Vector(+1.0,0,0), d1 ) );
+ addAtomDerivatives( 2, 1, crossProduct( Vector(+1.0,0,0), d2 ), myatoms );
+ addAtomDerivatives( 2, 2, crossProduct( Vector(-1.0,0,0), d1 ), myatoms );
+ addAtomDerivatives( 2, 3, crossProduct( Vector(+1.0,0,0), d1 ), myatoms );
}
myatoms.addBoxDerivatives( 2, Tensor(d1,crossProduct(Vector(+1.0,0,0), d2)) + Tensor( d2, crossProduct(Vector(-1.0,0,0), d1)) );
myatoms.addValue( 2, cp[0] );
- myatoms.addAtomsDerivatives( 3, 0, crossProduct( Vector(0,-1.0,0), d2 ) );
+ addAtomDerivatives( 3, 0, crossProduct( Vector(0,-1.0,0), d2 ), myatoms );
if( myatoms.getNumberOfAtoms()==3 ){
- myatoms.addAtomsDerivatives( 3, 1, crossProduct( Vector(0,+1.0,0), d2 ) + crossProduct( Vector(0,-1.0,0), d1 ) );
- myatoms.addAtomsDerivatives( 3, 2, crossProduct( Vector(0,+1.0,0), d1 ) );
+ addAtomDerivatives( 3, 1, crossProduct( Vector(0,+1.0,0), d2 ) + crossProduct( Vector(0,-1.0,0), d1 ), myatoms );
+ addAtomDerivatives( 3, 2, crossProduct( Vector(0,+1.0,0), d1 ), myatoms );
} else {
- myatoms.addAtomsDerivatives( 3, 1, crossProduct( Vector(0,+1.0,0), d2 ) );
- myatoms.addAtomsDerivatives( 3, 2, crossProduct( Vector(0,-1.0,0), d1 ) );
- myatoms.addAtomsDerivatives( 3, 3, crossProduct( Vector(0,+1.0,0), d1 ) );
+ addAtomDerivatives( 3, 1, crossProduct( Vector(0,+1.0,0), d2 ), myatoms );
+ addAtomDerivatives( 3, 2, crossProduct( Vector(0,-1.0,0), d1 ), myatoms );
+ addAtomDerivatives( 3, 3, crossProduct( Vector(0,+1.0,0), d1 ), myatoms );
}
myatoms.addBoxDerivatives( 3, Tensor(d1,crossProduct(Vector(0,+1.0,0), d2)) + Tensor( d2, crossProduct(Vector(0,-1.0,0), d1)) );
myatoms.addValue( 3, cp[1] );
- myatoms.addAtomsDerivatives( 4, 0, crossProduct( Vector(0,0,-1.0), d2 ) );
+ addAtomDerivatives( 4, 0, crossProduct( Vector(0,0,-1.0), d2 ), myatoms );
if( myatoms.getNumberOfAtoms()==3 ){
- myatoms.addAtomsDerivatives( 4, 1, crossProduct( Vector(0,0,+1.0), d2 ) + crossProduct( Vector(0,0,-1.0), d1 ) );
- myatoms.addAtomsDerivatives( 4, 2, crossProduct( Vector(0,0,+1.0), d1 ) );
+ addAtomDerivatives( 4, 1, crossProduct( Vector(0,0,+1.0), d2 ) + crossProduct( Vector(0,0,-1.0), d1 ), myatoms );
+ addAtomDerivatives( 4, 2, crossProduct( Vector(0,0,+1.0), d1 ), myatoms);
} else {
- myatoms.addAtomsDerivatives( 4, 1, crossProduct( Vector(0,0,-1.0), d2 ) );
- myatoms.addAtomsDerivatives( 4, 2, crossProduct( Vector(0,0,-1.0), d1 ) );
- myatoms.addAtomsDerivatives( 4, 3, crossProduct( Vector(0,0,+1.0), d1 ) );
+ addAtomDerivatives( 4, 1, crossProduct( Vector(0,0,-1.0), d2 ), myatoms);
+ addAtomDerivatives( 4, 2, crossProduct( Vector(0,0,-1.0), d1 ), myatoms);
+ addAtomDerivatives( 4, 3, crossProduct( Vector(0,0,+1.0), d1 ), myatoms);
}
myatoms.addBoxDerivatives( 4, Tensor(d1,crossProduct(Vector(0,0,+1.0), d2)) + Tensor( d2, crossProduct(Vector(0,0,-1.0), d1)) );
myatoms.addValue( 4, cp[2] );
diff --git a/src/crystallization/OrientationSphere.cpp b/src/crystallization/OrientationSphere.cpp
index 319de6dcb5871939d0224fa8e0398f4513b1b510..2ea92618c5b755049fbe5ce8352338433a0b8eb2 100644
--- a/src/crystallization/OrientationSphere.cpp
+++ b/src/crystallization/OrientationSphere.cpp
@@ -108,13 +108,12 @@ double OrientationSphere::compute( const unsigned& tindex, multicolvar::AtomValu
mergeVectorDerivatives( 1, 2, this_orient.size(), myatoms.getIndex(0), this_orient, myder0, myatoms );
mergeVectorDerivatives( 1, 2, catom_der.size(), myatoms.getIndex(i), catom_der, myder1, myatoms );
myatoms.addComDerivatives( 1, f_dot*(-dfunc)*distance, atom0 );
- multicolvar::CatomPack atom1=getCentralAtomPackFromInput( myatoms.getIndex(i) );
- myatoms.addComDerivatives( 1, f_dot*(dfunc)*distance, atom1 );
+ addAtomDerivatives( 1, i, f_dot*(dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 1, (-dfunc)*f_dot*Tensor(distance,distance) );
myder1.clearAll();
myatoms.addComDerivatives( 0, (-dfunc)*distance, atom0 );
- myatoms.addComDerivatives( 0, (dfunc)*distance, atom1 );
+ addAtomDerivatives( 0, i, (dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 0, (-dfunc)*Tensor(distance,distance) );
}
diff --git a/src/crystallization/SimpleCubic.cpp b/src/crystallization/SimpleCubic.cpp
index 3bb23d476da60cc10636d0b8f2d0bdadeb0ffe06..4165bf0b53aba129005d624fc35746820a1a378d 100644
--- a/src/crystallization/SimpleCubic.cpp
+++ b/src/crystallization/SimpleCubic.cpp
@@ -148,14 +148,14 @@ double SimpleCubic::compute( const unsigned& tindex, multicolvar::AtomValuePack&
myder[2] = 4*z3/t2-4*distance[2]*t3;
value += sw*tmp; fder = (+dfunc)*tmp*distance + sw*myder;
- myatoms.addAtomsDerivatives( 1, 0, -fder );
- myatoms.addAtomsDerivatives( 1, i, +fder );
+ addAtomDerivatives( 1, 0, -fder, myatoms );
+ addAtomDerivatives( 1, i, +fder, myatoms );
// Tens is a constructor that you build by taking the vector product of two vectors (remember the scalars!)
myatoms.addBoxDerivatives( 1, Tensor(distance,-fder) );
norm += sw;
- myatoms.addAtomsDerivatives( 0, 0, (-dfunc)*distance );
- myatoms.addAtomsDerivatives( 0, i, (+dfunc)*distance );
+ addAtomDerivatives( 0, 0, (-dfunc)*distance, myatoms );
+ addAtomDerivatives( 0, i, (+dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 0, (-dfunc)*Tensor(distance,distance) );
}
}
diff --git a/src/crystallization/Steinhardt.cpp b/src/crystallization/Steinhardt.cpp
index 3ce6c8be6d60991e018c3b1ec878f1e575ccf666..6041e0dd6597722552b4a8bd00462d181ac20e9d 100644
--- a/src/crystallization/Steinhardt.cpp
+++ b/src/crystallization/Steinhardt.cpp
@@ -86,8 +86,8 @@ void Steinhardt::calculateVector( multicolvar::AtomValuePack& myatoms ) const {
double dlen3 = d2*dlen;
// Store derivatives of weight
- myatoms.addAtomsDerivatives( 1, 0, (-dfunc)*distance );
- myatoms.addAtomsDerivatives( 1, i, (+dfunc)*distance );
+ addAtomDerivatives( 1, 0, (-dfunc)*distance, myatoms );
+ addAtomDerivatives( 1, i, (+dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 1, (-dfunc)*Tensor( distance,distance ) );
// Do stuff for m=0
@@ -97,8 +97,8 @@ void Steinhardt::calculateVector( multicolvar::AtomValuePack& myatoms ) const {
// Derivative wrt to the vector connecting the two atoms
myrealvec = (+sw)*dpoly_ass*dz + poly_ass*(+dfunc)*distance;
// Accumulate the derivatives
- myatoms.addAtomsDerivatives( 2 + tmom, 0, -myrealvec );
- myatoms.addAtomsDerivatives( 2 + tmom, i, myrealvec );
+ addAtomDerivatives( 2 + tmom, 0, -myrealvec, myatoms );
+ addAtomDerivatives( 2 + tmom, i, myrealvec, myatoms);
myatoms.addBoxDerivatives( 2 + tmom, Tensor( -myrealvec,distance ) );
// And store the vector function
myatoms.addValue( 2 + tmom, sw*poly_ass );
@@ -134,13 +134,13 @@ void Steinhardt::calculateVector( multicolvar::AtomValuePack& myatoms ) const {
// Real part
myatoms.addValue( 2+tmom+m, sw*tq6 );
- myatoms.addAtomsDerivatives( 2+tmom+m, 0, -myrealvec );
- myatoms.addAtomsDerivatives( 2+tmom+m, i, myrealvec );
+ addAtomDerivatives( 2+tmom+m, 0, -myrealvec, myatoms );
+ addAtomDerivatives( 2+tmom+m, i, myrealvec, myatoms );
myatoms.addBoxDerivatives( 2+tmom+m, Tensor( -myrealvec,distance ) );
// Imaginary part
myatoms.addValue( 2+ncomp+tmom+m, sw*itq6 );
- myatoms.addAtomsDerivatives( 2+ncomp+tmom+m, 0, -myimagvec );
- myatoms.addAtomsDerivatives( 2+ncomp+tmom+m, i, myimagvec );
+ addAtomDerivatives( 2+ncomp+tmom+m, 0, -myimagvec, myatoms );
+ addAtomDerivatives( 2+ncomp+tmom+m, i, myimagvec, myatoms );
myatoms.addBoxDerivatives( 2+ncomp+tmom+m, Tensor( -myimagvec,distance ) );
// Store -m part of vector
double pref=pow(-1.0,m);
@@ -148,13 +148,13 @@ void Steinhardt::calculateVector( multicolvar::AtomValuePack& myatoms ) const {
// conjugate of Legendre polynomial
// Real part
myatoms.addValue( 2+tmom-m, pref*sw*tq6 );
- myatoms.addAtomsDerivatives( 2+tmom-m, 0, -pref*myrealvec );
- myatoms.addAtomsDerivatives( 2+tmom-m, i, pref*myrealvec );
+ addAtomDerivatives( 2+tmom-m, 0, -pref*myrealvec, myatoms );
+ addAtomDerivatives( 2+tmom-m, i, pref*myrealvec, myatoms );
myatoms.addBoxDerivatives( 2+tmom-m, pref*Tensor( -myrealvec,distance ) );
// Imaginary part
myatoms.addValue( 2+ncomp+tmom-m, -pref*sw*itq6 );
- myatoms.addAtomsDerivatives( 2+ncomp+tmom-m, 0, pref*myimagvec );
- myatoms.addAtomsDerivatives( 2+ncomp+tmom-m, i, -pref*myimagvec );
+ addAtomDerivatives( 2+ncomp+tmom-m, 0, pref*myimagvec, myatoms );
+ addAtomDerivatives( 2+ncomp+tmom-m, i, -pref*myimagvec, myatoms );
myatoms.addBoxDerivatives( 2+ncomp+tmom-m, pref*Tensor( myimagvec,distance ) );
}
}
diff --git a/src/crystallization/Tetrahedral.cpp b/src/crystallization/Tetrahedral.cpp
index cb1f09d4084cb98dd87bd7cc22eacadbaed81043..a70da37acd24b240202172626c6de0616ae74e77 100644
--- a/src/crystallization/Tetrahedral.cpp
+++ b/src/crystallization/Tetrahedral.cpp
@@ -174,14 +174,14 @@ double Tetrahedral::compute( const unsigned& tindex, multicolvar::AtomValuePack&
myder[2] = (tt1-(distance[2]*t1)) + (-tt2-(distance[2]*t2)) + (-tt3-(distance[2]*t3)) + (tt4-(distance[2]*t4));
value += sw*tmp; fder = (+dfunc)*tmp*distance + sw*myder;
- myatoms.addAtomsDerivatives( 1, 0, -fder );
- myatoms.addAtomsDerivatives( 1, i, +fder );
+ addAtomDerivatives( 1, 0, -fder, myatoms );
+ addAtomDerivatives( 1, i, +fder, myatoms );
// Tens is a constructor that you build by taking the vector product of two vectors (remember the scalars!)
myatoms.addBoxDerivatives( 1, Tensor(distance,-fder) );
norm += sw;
- myatoms.addAtomsDerivatives( 0, 0, (-dfunc)*distance );
- myatoms.addAtomsDerivatives( 0, i, (+dfunc)*distance );
+ addAtomDerivatives( 0, 0, (-dfunc)*distance, myatoms );
+ addAtomDerivatives( 0, i, (+dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 0, (-dfunc)*Tensor(distance,distance) );
}
}
diff --git a/src/multicolvar/AlphaBeta.cpp b/src/multicolvar/AlphaBeta.cpp
index 307a3b26e2b4dd6e4ca4cd6eeafc8f52c0e4ec3c..0005d6951fb96ced13491ad7200531b945fd01a8 100644
--- a/src/multicolvar/AlphaBeta.cpp
+++ b/src/multicolvar/AlphaBeta.cpp
@@ -165,10 +165,10 @@ double AlphaBeta::compute( const unsigned& tindex, AtomValuePack& myatoms ) cons
dd2 *= svalue;
value = 0.5*cvalue;
- myatoms.addAtomsDerivatives(1, 0,dd0);
- myatoms.addAtomsDerivatives(1, 1,dd1-dd0);
- myatoms.addAtomsDerivatives(1, 2,dd2-dd1);
- myatoms.addAtomsDerivatives(1, 3,-dd2);
+ addAtomDerivatives(1, 0,dd0,myatoms);
+ addAtomDerivatives(1, 1,dd1-dd0,myatoms);
+ addAtomDerivatives(1, 2,dd2-dd1,myatoms);
+ addAtomDerivatives(1, 3,-dd2,myatoms);
myatoms.addBoxDerivatives(1, -(extProduct(d0,dd0)+extProduct(d1,dd1)+extProduct(d2,dd2)));
diff --git a/src/multicolvar/Angles.cpp b/src/multicolvar/Angles.cpp
index cf108fd5e92e85de7f31a999f64f436d483c79d0..714dfade11f2077e63bbe71b9f7e005f97d698cc 100644
--- a/src/multicolvar/Angles.cpp
+++ b/src/multicolvar/Angles.cpp
@@ -185,9 +185,9 @@ void Angles::calculateWeight( AtomValuePack& myatoms ) const {
wtot=w1*w2; dw1*=w2; dw2*=w1;
myatoms.setValue( 0, wtot );
- myatoms.addAtomsDerivatives( 0, 1, dw2*dik );
- myatoms.addAtomsDerivatives( 0, 0, -dw1*dij - dw2*dik );
- myatoms.addAtomsDerivatives( 0, 2, dw1*dij );
+ addAtomDerivatives( 0, 1, dw2*dik, myatoms );
+ addAtomDerivatives( 0, 0, -dw1*dij - dw2*dik, myatoms );
+ addAtomDerivatives( 0, 2, dw1*dij, myatoms );
myatoms.addBoxDerivatives( 0, (-dw1)*Tensor(dij,dij) + (-dw2)*Tensor(dik,dik) );
}
@@ -199,9 +199,9 @@ double Angles::compute( const unsigned& tindex, AtomValuePack& myatoms ) const {
double angle=a.compute(dij,dik,ddij,ddik);
// And finish the calculation
- myatoms.addAtomsDerivatives( 1, 1, ddik );
- myatoms.addAtomsDerivatives( 1, 0, - ddik - ddij );
- myatoms.addAtomsDerivatives( 1, 2, ddij );
+ addAtomDerivatives( 1, 1, ddik, myatoms );
+ addAtomDerivatives( 1, 0, - ddik - ddij, myatoms );
+ addAtomDerivatives( 1, 2, ddij, myatoms );
myatoms.addBoxDerivatives( 1, -(Tensor(dij,ddij)+Tensor(dik,ddik)) );
return angle;
diff --git a/src/multicolvar/AtomValuePack.h b/src/multicolvar/AtomValuePack.h
index ce2a81daac46d48c3b880e27e8fe2d6b4188a28e..be587009acdf2ef8b389eb19f10ac0c9a0534aae 100644
--- a/src/multicolvar/AtomValuePack.h
+++ b/src/multicolvar/AtomValuePack.h
@@ -34,6 +34,8 @@ namespace multicolvar {
class CatomPack;
class AtomValuePack {
+friend class MultiColvarBase;
+friend class LocalAverage;
private:
/// Copy of the values that we are adding to
MultiValue& myvals;
@@ -47,6 +49,10 @@ private:
std::vector<unsigned>& sort_vector;
/// This holds atom positions
std::vector<Vector>& myatoms;
+///
+ void addDerivative( const unsigned& , const unsigned& , const double& );
+///
+ void addAtomsDerivatives( const unsigned& , const unsigned& , const Vector& );
public:
AtomValuePack( MultiValue& vals, MultiColvarBase const * mcolv );
/// Set the number of atoms
@@ -71,10 +77,6 @@ public:
void addValue( const unsigned& ival, const double& vv );
///
double getValue( const unsigned& ) const ;
-///
- void addDerivative( const unsigned& , const unsigned& , const double& );
-///
- void addAtomsDerivatives( const unsigned& , const unsigned& , const Vector& );
///
void addBoxDerivatives( const unsigned& , const Tensor& );
///
diff --git a/src/multicolvar/Bridge.cpp b/src/multicolvar/Bridge.cpp
index e8352e31a08c37fcd69f4d83b84e962558ce24e0..11911b6fecb68c7722c93aacd64faafde5a3abf0 100644
--- a/src/multicolvar/Bridge.cpp
+++ b/src/multicolvar/Bridge.cpp
@@ -144,8 +144,8 @@ void Bridge::calculateWeight( AtomValuePack& myatoms ) const {
double dw, w=sf2.calculateSqr( ldij, dw );
myatoms.setValue( 0, w );
- myatoms.addAtomsDerivatives( 0, 0, -dw*dij );
- myatoms.addAtomsDerivatives( 0, 2, dw*dij );
+ addAtomDerivatives( 0, 0, -dw*dij, myatoms );
+ addAtomDerivatives( 0, 2, dw*dij, myatoms );
myatoms.addBoxDerivatives( 0, (-dw)*Tensor(dij,dij) );
}
@@ -154,8 +154,8 @@ double Bridge::compute( const unsigned& tindex, AtomValuePack& myatoms ) const {
double dw, w=sf1.calculateSqr( dik.modulo2(), dw );
// And finish the calculation
- myatoms.addAtomsDerivatives( 1, 0, -dw*dik );
- myatoms.addAtomsDerivatives( 1, 1, dw*dik );
+ addAtomDerivatives( 1, 0, -dw*dik, myatoms );
+ addAtomDerivatives( 1, 1, dw*dik, myatoms );
myatoms.addBoxDerivatives( 1, (-dw)*Tensor(dik,dik) );
return w;
}
diff --git a/src/multicolvar/CoordinationNumbers.cpp b/src/multicolvar/CoordinationNumbers.cpp
index ae618b6d5e357649b3a09240921d07a4f7e2005f..fae71b1c13a45c2b372f4575afb7cb0794f483c0 100644
--- a/src/multicolvar/CoordinationNumbers.cpp
+++ b/src/multicolvar/CoordinationNumbers.cpp
@@ -134,8 +134,8 @@ double CoordinationNumbers::compute( const unsigned& tindex, AtomValuePack& myat
sw = switchingFunction.calculateSqr( d2, dfunc );
value += sw;
- myatoms.addAtomsDerivatives( 1, 0, (-dfunc)*distance );
- myatoms.addAtomsDerivatives( 1, i, (dfunc)*distance );
+ addAtomDerivatives( 1, 0, (-dfunc)*distance, myatoms );
+ addAtomDerivatives( 1, i, (dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 1, (-dfunc)*Tensor(distance,distance) );
}
}
diff --git a/src/multicolvar/DihedralCorrelation.cpp b/src/multicolvar/DihedralCorrelation.cpp
index f5d1b7aa42a336ec7f64e7a0f96ef5697c9fdb16..5134aadfc70a8637fc0bdb554945d778b84f4dd7 100644
--- a/src/multicolvar/DihedralCorrelation.cpp
+++ b/src/multicolvar/DihedralCorrelation.cpp
@@ -141,20 +141,20 @@ double DihedralCorrelation::compute( const unsigned& tindex, AtomValuePack& myat
dd11 *= 0.5*sin( phi2 - phi1 );
dd12 *= 0.5*sin( phi2 - phi1 );
// And add
- myatoms.addAtomsDerivatives(1, 0, dd10);
- myatoms.addAtomsDerivatives(1, 1, dd11-dd10);
- myatoms.addAtomsDerivatives(1, 2, dd12-dd11);
- myatoms.addAtomsDerivatives(1, 3, -dd12);
+ addAtomDerivatives(1, 0, dd10, myatoms );
+ addAtomDerivatives(1, 1, dd11-dd10, myatoms );
+ addAtomDerivatives(1, 2, dd12-dd11, myatoms );
+ addAtomDerivatives(1, 3, -dd12, myatoms );
myatoms.addBoxDerivatives (1, -(extProduct(d10,dd10)+extProduct(d11,dd11)+extProduct(d12,dd12)));
// Derivative wrt phi2
dd20 *= -0.5*sin( phi2 - phi1 );
dd21 *= -0.5*sin( phi2 - phi1 );
dd22 *= -0.5*sin( phi2 - phi1 );
// And add
- myatoms.addAtomsDerivatives(1, 4, dd20);
- myatoms.addAtomsDerivatives(1, 5, dd21-dd20);
- myatoms.addAtomsDerivatives(1, 6, dd22-dd21);
- myatoms.addAtomsDerivatives(1, 7, -dd22);
+ addAtomDerivatives(1, 4, dd20, myatoms );
+ addAtomDerivatives(1, 5, dd21-dd20, myatoms );
+ addAtomDerivatives(1, 6, dd22-dd21, myatoms );
+ addAtomDerivatives(1, 7, -dd22, myatoms );
myatoms.addBoxDerivatives(1, -(extProduct(d20,dd20)+extProduct(d21,dd21)+extProduct(d22,dd22)));
return value;
diff --git a/src/multicolvar/Distances.cpp b/src/multicolvar/Distances.cpp
index 2bf9839d9b2cfaedfba90139f13ea6c6b636859d..83bc7e8f364759b761f6360afdf7c8c50042da3a 100644
--- a/src/multicolvar/Distances.cpp
+++ b/src/multicolvar/Distances.cpp
@@ -191,8 +191,8 @@ double Distances::compute( const unsigned& tindex, AtomValuePack& myatoms ) cons
const double invvalue=1.0/value;
// And finish the calculation
- myatoms.addAtomsDerivatives( 1, 0,-invvalue*distance );
- myatoms.addAtomsDerivatives( 1, 1, invvalue*distance );
+ addAtomDerivatives( 1, 0,-invvalue*distance, myatoms );
+ addAtomDerivatives( 1, 1, invvalue*distance, myatoms );
myatoms.addBoxDerivatives( 1, -invvalue*Tensor(distance,distance) );
return value;
}
diff --git a/src/multicolvar/InPlaneDistances.cpp b/src/multicolvar/InPlaneDistances.cpp
index 0fb720ea0eb1250fec427f3a6e1bb33daf1c1f65..1b59ac7d875abb594a9f0aa252c61282d085de8d 100644
--- a/src/multicolvar/InPlaneDistances.cpp
+++ b/src/multicolvar/InPlaneDistances.cpp
@@ -115,9 +115,9 @@ double InPlaneDistances::compute( const unsigned& tindex, AtomValuePack& myatoms
double sangle=sin(angle), cangle=cos(angle);
double dd=dir.modulo(), invdd=1.0/dd, val=dd*sangle;
- myatoms.addAtomsDerivatives( 1, 0, dd*cangle*ddik + sangle*invdd*dir );
- myatoms.addAtomsDerivatives( 1, 1, -dd*cangle*(ddik+ddij) - sangle*invdd*dir );
- myatoms.addAtomsDerivatives( 1, 2, dd*cangle*ddij );
+ addAtomDerivatives( 1, 0, dd*cangle*ddik + sangle*invdd*dir, myatoms );
+ addAtomDerivatives( 1, 1, -dd*cangle*(ddik+ddij) - sangle*invdd*dir, myatoms );
+ addAtomDerivatives( 1, 2, dd*cangle*ddij, myatoms );
myatoms.addBoxDerivatives( 1, -dd*cangle*(Tensor(normal,ddij)+Tensor(dir,ddik)) - sangle*invdd*Tensor(dir,dir) );
return val;
diff --git a/src/multicolvar/MultiColvarBase.cpp b/src/multicolvar/MultiColvarBase.cpp
index 63bdf985bff7b1b24dba5786b82ae5f7051141a7..39eb6039c3fb656353bfb2b4e534f9de1d9b56b1 100644
--- a/src/multicolvar/MultiColvarBase.cpp
+++ b/src/multicolvar/MultiColvarBase.cpp
@@ -388,6 +388,19 @@ void MultiColvarBase::calculate(){
runAllTasks();
}
+void MultiColvarBase::addAtomDerivatives( const unsigned& ival, const unsigned& iatom, const Vector& der, multicolvar::AtomValuePack& myatoms ) const {
+ unsigned jatom=myatoms.getIndex(iatom);
+
+ if( jatom<colvar_label.size() ){
+ unsigned mmc=colvar_label[jatom];
+ unsigned basen=0; for(unsigned i=0;i<mmc;++i) basen+=mybasemulticolvars[i]->getNumberOfAtoms();
+ multicolvar::CatomPack atom0=mybasemulticolvars[mmc]->getCentralAtomPack( basen, convertToLocalIndex(jatom,mmc) );
+ myatoms.addComDerivatives( ival, der, atom0 );
+ } else {
+ myatoms.addAtomsDerivatives( ival, iatom, der );
+ }
+}
+
void MultiColvarBase::performTask( const unsigned& task_index, const unsigned& current, MultiValue& myvals ) const {
AtomValuePack myatoms( myvals, this );
diff --git a/src/multicolvar/MultiColvarBase.h b/src/multicolvar/MultiColvarBase.h
index fd3c5268ce4ee49b651ca788ae088c3d37d22da1..0958b74d4b35501b8c5aabe1343ed741fef8d3d8 100644
--- a/src/multicolvar/MultiColvarBase.h
+++ b/src/multicolvar/MultiColvarBase.h
@@ -87,6 +87,8 @@ protected:
void addTaskToList( const unsigned& taskCode );
/// Finish setting up the multicolvar base
void setupMultiColvarBase( const std::vector<AtomNumber>& atoms );
+/// Add some derivatives to a particular component of a particular atom
+ void addAtomDerivatives( const unsigned& , const unsigned& , const Vector& , multicolvar::AtomValuePack& ) const ;
/// Set which atoms are to be used to calculate the central atom position
void setAtomsForCentralAtom( const std::vector<bool>& catom_ind );
/// Set the value of the cutoff for the link cells
@@ -179,7 +181,7 @@ Vector MultiColvarBase::getPositionOfAtomForLinkCells( const unsigned& iatom ) c
return mybasemulticolvars[mmc]->getCentralAtomPos( convertToLocalIndex(iatom,mmc) );
}
return ActionAtomistic::getPosition( iatom );
-}
+}
inline
unsigned MultiColvarBase::getNumberOfDerivatives(){
diff --git a/src/multicolvar/NumberOfLinks.cpp b/src/multicolvar/NumberOfLinks.cpp
index a810206ce72af816720cb98a7293f75ed6394951..011cf185909dc67dfe8589cd0ad3f3d367db1319 100644
--- a/src/multicolvar/NumberOfLinks.cpp
+++ b/src/multicolvar/NumberOfLinks.cpp
@@ -136,10 +136,8 @@ void NumberOfLinks::calculateWeight( AtomValuePack& myatoms ) const {
myatoms.setValue(0,sw);
if( !doNotCalculateDerivatives() ){
- CatomPack atom0=getCentralAtomPackFromInput( myatoms.getIndex(0) );
- myatoms.addComDerivatives( 0, (-dfunc)*distance, atom0 );
- CatomPack atom1=getCentralAtomPackFromInput( myatoms.getIndex(1) );
- myatoms.addComDerivatives( 0, (dfunc)*distance, atom1 );
+ addAtomDerivatives( 0, 0, (-dfunc)*distance, myatoms );
+ addAtomDerivatives( 0, 1, (dfunc)*distance, myatoms );
myatoms.addBoxDerivatives( 0, (-dfunc)*Tensor(distance,distance) );
}
}
diff --git a/src/multicolvar/Torsions.cpp b/src/multicolvar/Torsions.cpp
index 235a34283ba15eb65bd4d24d2adafb4040f69b26..93864c264a0f5b8206a96a89744bed9c0d99b41f 100644
--- a/src/multicolvar/Torsions.cpp
+++ b/src/multicolvar/Torsions.cpp
@@ -110,10 +110,10 @@ double Torsions::compute( const unsigned& tindex, AtomValuePack& myatoms ) const
Vector dd0,dd1,dd2; PLMD::Torsion t;
double value = t.compute(d0,d1,d2,dd0,dd1,dd2);
- myatoms.addAtomsDerivatives(1,0,dd0);
- myatoms.addAtomsDerivatives(1,1,dd1-dd0);
- myatoms.addAtomsDerivatives(1,2,dd2-dd1);
- myatoms.addAtomsDerivatives(1,3,-dd2);
+ addAtomDerivatives(1,0,dd0,myatoms);
+ addAtomDerivatives(1,1,dd1-dd0,myatoms);
+ addAtomDerivatives(1,2,dd2-dd1,myatoms);
+ addAtomDerivatives(1,3,-dd2,myatoms);
myatoms.addBoxDerivatives (1, -(extProduct(d0,dd0)+extProduct(d1,dd1)+extProduct(d2,dd2)));
diff --git a/src/multicolvar/XDistances.cpp b/src/multicolvar/XDistances.cpp
index 0333cf952a7b3694f29d3daec68ad7a7e4b06d1c..e34859880c2a1e3eea0badfec3cedea4dfc02e3e 100644
--- a/src/multicolvar/XDistances.cpp
+++ b/src/multicolvar/XDistances.cpp
@@ -158,8 +158,8 @@ double XDistances::compute( const unsigned& tindex, AtomValuePack& myatoms ) con
Vector myvec; myvec.zero();
// And finish the calculation
- myvec[myc]=+1; myatoms.addAtomsDerivatives( 1, 1, myvec );
- myvec[myc]=-1; myatoms.addAtomsDerivatives( 1, 0, myvec );
+ myvec[myc]=+1; addAtomDerivatives( 1, 1, myvec, myatoms );
+ myvec[myc]=-1; addAtomDerivatives( 1, 0, myvec, myatoms );
myatoms.addBoxDerivatives( 1, Tensor(distance,myvec) );
return value;
}
diff --git a/src/multicolvar/XYDistances.cpp b/src/multicolvar/XYDistances.cpp
index 9836c0d105aeb9908fb66861f1383b7938ecd3d2..2919b39c2a2b3fd6cd635467b9cd70be73d4b995 100644
--- a/src/multicolvar/XYDistances.cpp
+++ b/src/multicolvar/XYDistances.cpp
@@ -136,8 +136,8 @@ double XYDistances::compute( const unsigned& tindex, AtomValuePack& myatoms ) co
Vector myvec; myvec.zero();
// And finish the calculation
- myvec[myc1]=+invvalue*distance[myc1]; myvec[myc2]=+invvalue*distance[myc2]; myatoms.addAtomsDerivatives( 1, 1, myvec );
- myvec[myc1]=-invvalue*distance[myc1]; myvec[myc2]=-invvalue*distance[myc2]; myatoms.addAtomsDerivatives( 1, 0, myvec );
+ myvec[myc1]=+invvalue*distance[myc1]; myvec[myc2]=+invvalue*distance[myc2]; addAtomDerivatives( 1, 1, myvec, myatoms );
+ myvec[myc1]=-invvalue*distance[myc1]; myvec[myc2]=-invvalue*distance[myc2]; addAtomDerivatives( 1, 0, myvec, myatoms );
myatoms.addBoxDerivatives( 1, Tensor(distance,myvec) );
return value;
}