diff --git a/CHANGES/v2.3.txt b/CHANGES/v2.3.txt
index d5aa9352d26b82e48f3bb88ad164a50f43d49a5e..b2586b9e3053fa92aa0801743ce246c9b56aa394 100644
--- a/CHANGES/v2.3.txt
+++ b/CHANGES/v2.3.txt
@@ -160,6 +160,7 @@ See branch \branch{v2.3} on git repository.
- Prevented users from causing segfaults by storing derivatives without LOWMEM flag. In these caess PLUMED crashes with meaningful errors.
- Fixed bug in \ref HISTOGRAM that causes nans when using KERNEL=DISCRETE option
- Fixed a bug in the parser related to braces, see \issue{229}
+- Fixed a bug that appeared when using \ref Q3, \ref Q4 and \ref Q6 with LOWEST or HIGHEST flag
For developers:
- plumedcheck validation has been made stricter. All the checks are now described in the developer manual.
diff --git a/regtest/crystallization/rt-q6/colv.reference b/regtest/crystallization/rt-q6/colv.reference
index 8d69c80ee2192c76108352c9c7255dee723a579b..320d6162236a8580640400b282a440620d62fb12 100644
--- a/regtest/crystallization/rt-q6/colv.reference
+++ b/regtest/crystallization/rt-q6/colv.reference
@@ -1,2 +1,2 @@
-#! FIELDS time q6.mean
- 0.000000 0.245096
+#! FIELDS time q6.mean q6.lowest
+ 0.000000 0.245096 0.150123
diff --git a/regtest/crystallization/rt-q6/plumed.dat b/regtest/crystallization/rt-q6/plumed.dat
index f09db50863178d3da093fdb70496de4401793b76..d107217876b4cb22783e894547bae7bb623d4e55 100644
--- a/regtest/crystallization/rt-q6/plumed.dat
+++ b/regtest/crystallization/rt-q6/plumed.dat
@@ -3,10 +3,10 @@ COORDINATIONNUMBER SPECIES=1-64 SWITCH={RATIONAL D_0=3.0 R_0=1.5} MEAN LABEL=c
PRINT ARG=c.* FILE=colv2
DUMPDERIVATIVES ARG=c.* FILE=deriv2 FMT=%8.4f
-Q6 SPECIES=1-64 D_0=3.0 R_0=1.5 MEAN LABEL=q6
+Q6 SPECIES=1-64 D_0=3.0 R_0=1.5 MEAN LOWEST LABEL=q6
# Q6 SPECIES=1-64 D_0=3.0 R_0=1.5 MEAN NUMERICAL_DERIVATIVES LABEL=q6n
PRINT ARG=q6.* FILE=colv
-DUMPDERIVATIVES ARG=q6.* FILE=deriv FMT=%8.4f
+DUMPDERIVATIVES ARG=q6.mean FILE=deriv FMT=%8.4f
Q4 SPECIES=1-64 D_0=3.0 R_0=1.5 MEAN LABEL=q4
# Q4 SPECIES=1-64 D_0=3.0 R_0=1.5 MEAN NUMERICAL_DERIVATIVES LABEL=q4n
diff --git a/src/multicolvar/CoordinationNumbers.cpp b/src/multicolvar/CoordinationNumbers.cpp
index e280931ed3a2a0d9e9dbd5f6b5d4b6ac7b5664b4..ea2723e1009d709a9df627da6b51f343121d8ff2 100644
--- a/src/multicolvar/CoordinationNumbers.cpp
+++ b/src/multicolvar/CoordinationNumbers.cpp
@@ -104,10 +104,10 @@ void CoordinationNumbers::registerKeywords( Keywords& keys ) {
keys.add("compulsory","D_0","0.0","The d_0 parameter of the switching function");
keys.add("compulsory","R_0","The r_0 parameter of the switching function");
keys.add("optional","R_POWER","Multiply the coordination number function by a power of r, "
- "as done in White and Voth (see note above, default: no)");
+ "as done in White and Voth (see note above, default: no)");
keys.add("optional","SWITCH","This keyword is used if you want to employ an alternative to the continuous swiching function defined above. "
- "The following provides information on the \\ref switchingfunction that are available. "
- "When this keyword is present you no longer need the NN, MM, D_0 and R_0 keywords.");
+ "The following provides information on the \\ref switchingfunction that are available. "
+ "When this keyword is present you no longer need the NN, MM, D_0 and R_0 keywords.");
// Use actionWithDistributionKeywords
keys.use("MEAN"); keys.use("MORE_THAN"); keys.use("LESS_THAN"); keys.use("MAX");
keys.use("MIN"); keys.use("BETWEEN"); keys.use("HISTOGRAM"); keys.use("MOMENTS");
@@ -144,7 +144,7 @@ CoordinationNumbers::CoordinationNumbers(const ActionOptions&ao):
log.printf(" Multiplying switching function by r^%d\n", r_power);
double offset = switchingFunction.calculate(rcut*0.9999, rcut2) * pow(rcut*0.9999, r_power);
log.printf(" You will have a discontinuous jump of %f to 0 near the cutoff of your switching function. "
- "Consider setting D_MAX or reducing R_POWER if this is large\n", offset);
+ "Consider setting D_MAX or reducing R_POWER if this is large\n", offset);
}
// Set the link cell cutoff
@@ -161,18 +161,18 @@ double CoordinationNumbers::compute( const unsigned& tindex, AtomValuePack& myat
for(unsigned i=1; i<myatoms.getNumberOfAtoms(); ++i) {
Vector& distance=myatoms.getPosition(i);
if ( (d2=distance[0]*distance[0])<rcut2 &&
- (d2+=distance[1]*distance[1])<rcut2 &&
- (d2+=distance[2]*distance[2])<rcut2) {
+ (d2+=distance[1]*distance[1])<rcut2 &&
+ (d2+=distance[2]*distance[2])<rcut2) {
sw = switchingFunction.calculateSqr( d2, dfunc );
if(r_power > 0) {
- d = sqrt(d2); raised = pow( d, r_power - 1 );
- accumulateSymmetryFunction( 1, i, sw * raised * d,
- (dfunc * d * raised + sw * r_power) * distance,
- (-dfunc * d * raised - sw * r_power) * Tensor(distance, distance),
- myatoms );
+ d = sqrt(d2); raised = pow( d, r_power - 1 );
+ accumulateSymmetryFunction( 1, i, sw * raised * d,
+ (dfunc * d * raised + sw * r_power) * distance,
+ (-dfunc * d * raised - sw * r_power) * Tensor(distance, distance),
+ myatoms );
} else {
- accumulateSymmetryFunction( 1, i, sw, (dfunc)*distance, (-dfunc)*Tensor(distance,distance), myatoms );
+ accumulateSymmetryFunction( 1, i, sw, (dfunc)*distance, (-dfunc)*Tensor(distance,distance), myatoms );
}
}
}
diff --git a/src/vesselbase/StoreDataVessel.cpp b/src/vesselbase/StoreDataVessel.cpp
index a89c50d1492a15febecde8267b38795300b40ce3..72c831290fd15d2d7241c10e7763702d5a4b5ebd 100644
--- a/src/vesselbase/StoreDataVessel.cpp
+++ b/src/vesselbase/StoreDataVessel.cpp
@@ -44,9 +44,6 @@ void StoreDataVessel::addActionThatUses( ActionWithVessel* actionThatUses ) {
}
void StoreDataVessel::resize() {
- vecsize=getAction()->getNumberOfQuantities();
- plumed_dbg_assert( vecsize>0 );
-
if( getAction()->lowmem || !getAction()->derivativesAreRequired() ) {
nspace = 1;
active_der.resize( max_lowmem_stash * ( 1 + getAction()->getNumberOfDerivatives() ) );
@@ -57,6 +54,8 @@ void StoreDataVessel::resize() {
nspace = 1 + getAction()->maxderivatives;
active_der.resize( getNumberOfStoredValues() * ( 1 + getAction()->maxderivatives ) );
}
+ vecsize=getAction()->getNumberOfQuantities();
+ plumed_dbg_assert( vecsize>0 );
resizeBuffer( getNumberOfStoredValues()*vecsize*nspace );
local_buffer.resize( getNumberOfStoredValues()*vecsize*nspace );
}
@@ -102,14 +101,14 @@ void StoreDataVessel::storeDerivatives( const unsigned& myelem, MultiValue& myva
}
void StoreDataVessel::retrieveSequentialValue( const unsigned& jelem, const bool& normed, std::vector<double>& values ) const {
- plumed_assert( values.size()==vecsize );
+ plumed_dbg_assert( values.size()==vecsize );
unsigned ibuf = jelem * vecsize * nspace;
for(unsigned i=0; i<vecsize; ++i) { values[i]=local_buffer[ibuf]; ibuf+=nspace; }
if( normed && values.size()>2 ) getAction()->normalizeVector( values );
}
void StoreDataVessel::retrieveValueWithIndex( const unsigned& myelem, const bool& normed, std::vector<double>& values ) const {
- plumed_assert( values.size()==vecsize );
+ plumed_dbg_assert( values.size()==vecsize );
unsigned jelem = getStoreIndex( myelem );
retrieveSequentialValue( jelem, normed, values );
}