diff --git a/regtest/rt19-mpi/plumed.dat b/regtest/rt19-mpi/plumed.dat
index 255e12a62a0ae0d8585f95bdc00fd1c1a5cdc5eb..fa64e0c134e5dd2714674bca20c009a6748f0fef 100644
--- a/regtest/rt19-mpi/plumed.dat
+++ b/regtest/rt19-mpi/plumed.dat
@@ -1,5 +1,5 @@
-DISTANCES ATOMS1=1,2 ATOMS2=5,4 HISTOGRAM={NBINS=3 LOWER=0.5 UPPER=2.0} LABEL=d1
-DISTANCES ATOMS1=1,2 ATOMS2=5,4 HISTOGRAM={NBINS=3 LOWER=0.5 UPPER=2.0} NUMERICAL_DERIVATIVES LABEL=d1n
+DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={NBINS=3 LOWER=0.5 UPPER=2.0} LABEL=d1
+DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={NBINS=3 LOWER=0.5 UPPER=2.0} NUMERICAL_DERIVATIVES LABEL=d1n
DUMPDERIVATIVES ARG=d1.*,d1n.* FILE=derivatives1 FMT=%8.5f
DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={LOWER=0.5 UPPER=1.5} LABEL=d2
DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={LOWER=0.5 UPPER=1.5} NUMERICAL_DERIVATIVES LABEL=d2n
diff --git a/regtest/rt19/plumed.dat b/regtest/rt19/plumed.dat
index 255e12a62a0ae0d8585f95bdc00fd1c1a5cdc5eb..fa64e0c134e5dd2714674bca20c009a6748f0fef 100644
--- a/regtest/rt19/plumed.dat
+++ b/regtest/rt19/plumed.dat
@@ -1,5 +1,5 @@
-DISTANCES ATOMS1=1,2 ATOMS2=5,4 HISTOGRAM={NBINS=3 LOWER=0.5 UPPER=2.0} LABEL=d1
-DISTANCES ATOMS1=1,2 ATOMS2=5,4 HISTOGRAM={NBINS=3 LOWER=0.5 UPPER=2.0} NUMERICAL_DERIVATIVES LABEL=d1n
+DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={NBINS=3 LOWER=0.5 UPPER=2.0} LABEL=d1
+DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={NBINS=3 LOWER=0.5 UPPER=2.0} NUMERICAL_DERIVATIVES LABEL=d1n
DUMPDERIVATIVES ARG=d1.*,d1n.* FILE=derivatives1 FMT=%8.5f
DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={LOWER=0.5 UPPER=1.5} LABEL=d2
DISTANCES ATOMS1=1,2 ATOMS2=5,4 WITHIN={LOWER=0.5 UPPER=1.5} NUMERICAL_DERIVATIVES LABEL=d2n
diff --git a/src/DistributionHistogram.cpp b/src/DistributionHistogram.cpp
deleted file mode 100644
index 70ee83ca650ce84e360c8e1ae78598952abd9daa..0000000000000000000000000000000000000000
--- a/src/DistributionHistogram.cpp
+++ /dev/null
@@ -1,101 +0,0 @@
-/* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
- Copyright (c) 2012 The plumed team
- (see the PEOPLE file at the root of the distribution for a list of names)
-
- See http://www.plumed-code.org for more information.
-
- This file is part of plumed, version 2.0.
-
- plumed is free software: you can redistribute it and/or modify
- it under the terms of the GNU Lesser General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- plumed is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- GNU Lesser General Public License for more details.
-
- You should have received a copy of the GNU Lesser General Public License
- along with plumed. If not, see <http://www.gnu.org/licenses/>.
-+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
-#include "FunctionVessel.h"
-#include "MultiColvar.h"
-#include "HistogramBead.h"
-
-namespace PLMD {
-
-class histogram : public NormedSumVessel {
-private:
- MultiColvar* mycolv;
- std::vector<HistogramBead> hist;
-public:
- static void reserveKeyword( Keywords& keys );
- histogram( const VesselOptions& da );
- void getWeight( const unsigned& i, Value& weight );
- void compute( const unsigned& i, const unsigned& j, Value& theval );
- void printKeywords();
-};
-
-PLUMED_REGISTER_VESSEL(histogram,"HISTOGRAM")
-
-void histogram::reserveKeyword( Keywords& keys ){
- keys.reserve("optional","HISTOGRAM", "create a discretized histogram of the distribution of collective variables. " + HistogramBead::histodocs() );
-}
-
-histogram::histogram( const VesselOptions& da ) :
-NormedSumVessel(da)
-{
-
- mycolv=dynamic_cast<MultiColvar*>( getAction() );
- plumed_massert( mycolv, "histogram is used to calculate functions of multi colvars");
-
- std::string errormsg;
- std::vector<std::string> data=Tools::getWords(da.parameters);
- bool usenorm=false; Tools::parseFlag(data,"NORM",usenorm);
- if(usenorm) useNorm();
-
- bool isPeriodic=getAction()->isPeriodic();
- double min, max;
- if( isPeriodic ) getAction()->retrieveDomain( min, max );
-
- std::vector<std::string> bins; HistogramBead::generateBins( da.parameters, "", bins );
- hist.resize( bins.size() );
- for(unsigned i=0;i<hist.size();++i){
- hist[i].set( bins[i], "", errormsg );
- if( !isPeriodic ) hist[i].isNotPeriodic();
- else hist[i].isPeriodic( min, max );
- if( errormsg.size()!=0 ) error( errormsg );
-
- std::string lb, ub;
- Tools::convert( hist[i].getlowb(), lb );
- Tools::convert( hist[i].getbigb(), ub );
- addOutput( "between" + lb + "&" + ub );
- log.printf(" value %s.between%s&%s contains the ",(getAction()->getLabel()).c_str(),lb.c_str(),ub.c_str());
- if(usenorm) log.printf("fraction of values %s\n",(hist[i].description()).c_str());
- else log.printf("number of values %s\n",(hist[i].description()).c_str());
- }
-}
-
-void histogram::printKeywords(){
- Keywords keys;
- keys.add("compulsory","NBINS","the number of bins in the histogram");
- keys.add("compulsory","LOWER","the lower bound for the histogram");
- keys.add("compulsory","UPPER","the upper boudn for the histogram");
- keys.add("compulsory","SMEAR","0.5","the ammount to smear the values by to smooth the histogram");
- keys.addFlag("NORM",false,"normalize the histogram");
- keys.print(log);
-}
-
-void histogram::compute( const unsigned& icv, const unsigned& jfunc, Value& theval ){
- plumed_assert( jfunc<hist.size() );
- mycolv->retreiveLastCalculatedValue( theval );
- double df, f; f=hist[jfunc].calculate( theval.get() , df );
- theval.chainRule(df); theval.set(f);
-}
-
-void histogram::getWeight( const unsigned& i, Value& weight ){
- mycolv->retrieveColvarWeight( i, weight );
-}
-
-}
diff --git a/src/DistributionWithin.cpp b/src/DistributionWithin.cpp
index 6fddea3359aa985b5f076bfe094aba7c8c13f5ab..8b8f016e744c20c4dd8df20636f6d642663e30c6 100644
--- a/src/DistributionWithin.cpp
+++ b/src/DistributionWithin.cpp
@@ -1,25 +1,3 @@
-/* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
- Copyright (c) 2012 The plumed team
- (see the PEOPLE file at the root of the distribution for a list of names)
-
- See http://www.plumed-code.org for more information.
-
- This file is part of plumed, version 2.0.
-
- plumed is free software: you can redistribute it and/or modify
- it under the terms of the GNU Lesser General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- plumed is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- GNU Lesser General Public License for more details.
-
- You should have received a copy of the GNU Lesser General Public License
- along with plumed. If not, see <http://www.gnu.org/licenses/>.
-+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
-
#include "FunctionVessel.h"
#include "MultiColvar.h"
#include "HistogramBead.h"
@@ -29,7 +7,7 @@ namespace PLMD {
class within : public NormedSumVessel {
private:
MultiColvar* mycolv;
- HistogramBead hist;
+ std::vector<HistogramBead> hist;
public:
static void reserveKeyword( Keywords& keys );
within( const VesselOptions& da );
@@ -41,10 +19,11 @@ public:
PLUMED_REGISTER_VESSEL(within,"WITHIN")
void within::reserveKeyword( Keywords& keys ){
- keys.reserve("numbered", "WITHIN", "calculate the number variabels that are within a certain range and store it in a value called between<lowerbound>&<upperbound>. "
- "To make this quantity continuous it is calculated using " + HistogramBead::documentation(false) + ". "
- "Adding the NORM flag allows you to calculate the fraction of colvars in the particular range rather than the total number "
- "(N.B. this option should probably not used if you are using neighbor lists and the derivatives of the WITHIN value).");
+ keys.reserve("numbered","WITHIN", "calculate the number of variables that are within a certain range and store it in a value called between<lowerbound>&<upperbound> "
+ "or create a discretized histogram of the distribution for a particular range. To make these quantities continuous they are "
+ "calculated using " + HistogramBead::documentation(false) + " If you add the NBINS keyword the range between your upper and "
+ "lower bounds is divided into a discrete number of bins. Adding the NORM flag will calculate the fraction of colvars in "
+ "range of interest rather than the total number");
}
within::within( const VesselOptions& da ) :
@@ -54,38 +33,54 @@ NormedSumVessel(da)
mycolv=dynamic_cast<MultiColvar*>( getAction() );
plumed_massert( mycolv, "within is used to calculate functions of multi colvars");
- std::string errormsg;
- std::vector<std::string> data=Tools::getWords(da.parameters);
- bool usenorm=false; Tools::parseFlag(data,"NORM",usenorm);
- if(usenorm) useNorm();
-
bool isPeriodic=getAction()->isPeriodic();
double min, max;
if( isPeriodic ) getAction()->retrieveDomain( min, max );
- hist.set( da.parameters, "", errormsg );
- if( !isPeriodic ) hist.isNotPeriodic();
- else hist.isPeriodic( min, max );
-
- if( errormsg.size()!=0 ) error( errormsg );
-
- std::string lb, ub;
- Tools::convert( hist.getlowb(), lb );
- Tools::convert( hist.getbigb(), ub );
- addOutput( "between" + lb + "&" + ub );
- log.printf(" value %s.between%s&%s contains the ",(getAction()->getLabel()).c_str(),lb.c_str(),ub.c_str());
- if(usenorm) log.printf("fraction of values %s\n",(hist.description()).c_str());
- else log.printf("number of values %s\n",(hist.description()).c_str());
+ std::string errormsg;
+ std::vector<std::string> data=Tools::getWords(da.parameters);
+ bool usenorm=false; Tools::parseFlag(data,"NORM",usenorm);
+ if(usenorm) useNorm();
+ bool hasbins=false; unsigned nbins=1;
+ hasbins=Tools::parse(data,"NBINS",nbins);
+ if(!hasbins){
+ hist.resize(1);
+ hist[0].set( da.parameters,"",errormsg );
+
+ } else {
+ std::vector<std::string> bins; HistogramBead::generateBins( da.parameters, "", bins );
+ hist.resize( bins.size() );
+ for(unsigned i=0;i<hist.size();++i) hist[i].set( bins[i], "", errormsg );
+ }
+ for(unsigned i=0;i<hist.size();++i){
+ if( !isPeriodic ) hist[i].isNotPeriodic();
+ else hist[i].isPeriodic( min, max );
+ if( errormsg.size()!=0 ) error( errormsg );
+
+ std::string lb, ub;
+ Tools::convert( hist[i].getlowb(), lb );
+ Tools::convert( hist[i].getbigb(), ub );
+ addOutput( "between" + lb + "&" + ub );
+ log.printf(" value %s.between%s&%s contains the ",(getAction()->getLabel()).c_str(),lb.c_str(),ub.c_str());
+ if(usenorm) log.printf("fraction of values %s\n",(hist[i].description()).c_str());
+ else log.printf("number of values %s\n",(hist[i].description()).c_str());
+ }
}
void within::printKeywords(){
- hist.printKeywords( log );
+ Keywords keys;
+ keys.add("compulsory","NBINS","1","the number of bins you wish to divide the range into");
+ keys.add("compulsory","LOWER","the lower bound");
+ keys.add("compulsory","UPPER","the upper bound");
+ keys.add("compulsory","SMEAR","0.5","the ammount to smear the values by to smooth the histogram");
+ keys.addFlag("NORM",false,"normalize the histogram");
+ keys.print(log);
}
-void within::compute( const unsigned& i, const unsigned& j, Value& theval ){
- plumed_assert( j==0 );
+void within::compute( const unsigned& icv, const unsigned& jfunc, Value& theval ){
+ plumed_assert( jfunc<hist.size() );
mycolv->retreiveLastCalculatedValue( theval );
- double df, f; f=hist.calculate( theval.get() , df );
+ double df, f; f=hist[jfunc].calculate( theval.get() , df );
theval.chainRule(df); theval.set(f);
}
diff --git a/src/HistogramBead.cpp b/src/HistogramBead.cpp
index 02359bb159fd75e4787a5b825f9c271f39b9f19c..55de76715c0dcb34fb178a0ce013eb8e2edd19c8 100644
--- a/src/HistogramBead.cpp
+++ b/src/HistogramBead.cpp
@@ -32,11 +32,12 @@ std::string HistogramBead::documentation( bool dir ) {
if(dir){
ostr<<"\\f$ w(x)=\\int_a^b \\frac{1}{\\sqrt{2\\pi}\\sigma_d} \\exp\\left( -\\frac{(x'-x)^2}{2\\sigma_x^2} \\right) \\textrm{d}x' \\f$";
ostr<<"where \\f$ \\sigma_x=(b_x-a_x)k_x \\f$. The parameters of the functions are specifed in fractional coordinates using ";
- ostr<<"(XLOWER=\\f$a_x\\f$ XUPPER=\\f$b_x\\f$ XSMEAR=\\f$k_x\\f$ YLOWER=\\f$a_y\\f$ YUPPER=\\f$b_y\\f$ YSMEAR=\\f$k_y\\f$ ZLOWER=\\f$a_z\\f$ ZUPPER=\\f$b_z\\f$ ZSMEAR=\\f$k_z\\f$).";
+ ostr<<"{XLOWER=\\f$a_x\\f$ XUPPER=\\f$b_x\\f$ XSMEAR=\\f$k_x\\f$ YLOWER=\\f$a_y\\f$ YUPPER=\\f$b_y\\f$ YSMEAR=\\f$k_y\\f$ ZLOWER=\\f$a_z\\f$ ZUPPER=\\f$b_z\\f$ ZSMEAR=\\f$k_z\\f$}.";
ostr<<"If any of the SMEAR keywords are not present then the default \\f$k_x=0.5\\f$ is used in that direction. ";
} else {
ostr<<"\\f$ w(r)=\\int_a^b \\frac{1}{\\sqrt{2\\pi}\\sigma} \\exp\\left( -\\frac{(r - r')^2}{2\\sigma^2} \\right) \\textrm{d}r' \\f$";
- ostr<<"where \\f$ \\sigma=(b-a)k \\f$. The parameters of the function are specified using (LOWER=\\f$a\\f$ UPPER=\\f$b\\f$ SMEAR=\\f$k\\f$). ";
+ ostr<<"where \\f$ \\sigma=k\\frac{b-a}{n} \\f$ and \\f$n\\f$ is the number of bins you are dividing the range into. ";
+ ostr<<"The parameters of the function are specified using {LOWER=\\f$a\\f$ UPPER=\\f$b\\f$ SMEAR=\\f$k\\f$}. ";
ostr<<"If the SMEAR keyword is not present then by default \\f$k=0.5\\f$.";
}
return ostr.str();
@@ -48,16 +49,6 @@ std::string HistogramBead::description() const {
return ostr.str();
}
-std::string HistogramBead::histodocs() {
- std::ostringstream ostr;
- ostr<<"The range is divided into a discete number of bins and the number of values that fall within each bin is calculated using ";
- ostr<<"\\f$ w(r)=\\int_a^b \\frac{1}{\\sqrt{2\\pi}\\sigma} \\exp\\left( -\\frac{(r'-r)^2}{2\\sigma^2} \\right) \\textrm{d}r' \\f$";
- ostr<<"where \\f$ \\sigma=(b-a)k \\f$. The particular range of interest and number of bins are specified using ";
- ostr<<"(NBINS=\\f$n\\f$ LOWER=\\f$a\\f$ UPPER=\\f$b\\f$ SMEAR=\\f$x\\f$). If the SMEAR keyword is not present then by default \\f$k=0.5\\f$. ";
- ostr<<"You can calculate a normalized histogram using the NORM flag (N.B. Don't use this if you are using derivatives of the histogram and neighbor lists)";
- return ostr.str();
-}
-
void HistogramBead::generateBins( const std::string& params, const std::string& dd, std::vector<std::string>& bins ){
if( dd.size()!=0 && params.find(dd)==std::string::npos) return;
std::vector<std::string> data=Tools::getWords(params);
diff --git a/src/HistogramBead.h b/src/HistogramBead.h
index dc1c5ff9dc4ee79626e917d22b8286e26b299bc0..35daa40150921dfa6bbbd18a6c792547f2facecc 100644
--- a/src/HistogramBead.h
+++ b/src/HistogramBead.h
@@ -49,7 +49,6 @@ private:
double difference( const double& d1, const double& d2 ) const ;
public:
static std::string documentation( bool dir );
- static std::string histodocs();
static void generateBins( const std::string& params, const std::string& dd, std::vector<std::string>& bins );
HistogramBead();
std::string description() const ;
diff --git a/src/MultiColvarCoordination.cpp b/src/MultiColvarCoordination.cpp
index dfa53f5405034a8a80fa738cdac4f8a3a2d886b4..9ff837806bc93a2ecc183f3c19ac593ff4e72ce8 100644
--- a/src/MultiColvarCoordination.cpp
+++ b/src/MultiColvarCoordination.cpp
@@ -84,7 +84,7 @@ void MultiColvarCoordination::registerKeywords( Keywords& keys ){
keys.add("optional","R_0","The r_0 parameter of the switching function");
// Use actionWithDistributionKeywords
keys.use("AVERAGE"); keys.use("MORE_THAN"); keys.use("LESS_THAN");
- keys.use("MIN"); keys.use("WITHIN"); keys.use("HISTOGRAM"); keys.use("MOMENTS");
+ keys.use("MIN"); keys.use("WITHIN"); keys.use("MOMENTS");
keys.use("SUBCELL");
}
diff --git a/src/MultiColvarDistance.cpp b/src/MultiColvarDistance.cpp
index 81cbd71872768598fa4b6cc3287059874d010ce4..4cbdd9032785b3fcb7d50cbc591c51e7a81d0ba1 100644
--- a/src/MultiColvarDistance.cpp
+++ b/src/MultiColvarDistance.cpp
@@ -84,7 +84,7 @@ void MultiColvarDistance::registerKeywords( Keywords& keys ){
ActionWithDistribution::autoParallelize( keys );
keys.use("ATOMS"); keys.use("GROUP"); keys.use("GROUPA"); keys.use("GROUPB");
keys.use("AVERAGE"); keys.use("MIN"); keys.use("LESS_THAN");
- keys.use("MORE_THAN"); keys.use("WITHIN"); keys.use("HISTOGRAM"); keys.use("MOMENTS");
+ keys.use("MORE_THAN"); keys.use("WITHIN"); keys.use("MOMENTS");
}
MultiColvarDistance::MultiColvarDistance(const ActionOptions&ao):