Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

New warnings for METAD adn WALKERS_MPI #931

Merged
merged 15 commits into from
May 14, 2023
Merged
1 change: 1 addition & 0 deletions regtest/basic/rt-errormessages/Makefile
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
include ../../scripts/test.make
3 changes: 3 additions & 0 deletions regtest/basic/rt-errormessages/config
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
type=make
# this is needed because we are explicitly checking some exception message
export PLUMED_STACK_TRACE=no
1 change: 1 addition & 0 deletions regtest/basic/rt-errormessages/output.reference
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
METAD WALKER_MPI : Exception thrown, correct message
77 changes: 77 additions & 0 deletions regtest/basic/rt-errormessages/testWALKERS_MPI.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,77 @@
#include <fstream>

#include "plumed/wrapper/Plumed.h"
#include "plumed/tools/Communicator.h"

constexpr unsigned nat=10;
//it is a struct because we don't need getter/setters (for now)
struct plumedThrowChecker{
unsigned natoms=nat;
std::vector<double> positions=[](unsigned natoms){
std::vector<double> toret(natoms*3,0.0);
for(unsigned i=0; i<3*natoms; i++) toret[i]=i;
return toret;
}(nat);
std::vector<double> masses{nat,1.0};
std::vector<double> forces{3*nat,0.0};
std::vector<double> box{9,0.0};
std::vector<double> virial{9,0.0};
///This initializes a new plumed to test each throw in the cleanliest way possible
int checkThrow (std::ostream& out,std::string name, std::string cmd, std::string expectedMessage){
//GIVEN an initialized Plumed interface
PLMD::Plumed plumed;

plumed.cmd("setNatoms",&natoms);
int step=0;
plumed.cmd("setLogFile","test.log");
plumed.cmd("init");
plumed.cmd("setStep",&step);
plumed.cmd("setPositions",positions.data());
plumed.cmd("setBox",box.data());
plumed.cmd("setForces",forces.data());
plumed.cmd("setVirial",virial.data());
plumed.cmd("setMasses",masses.data());

plumed.cmd("readInputLine","d: DISTANCE ATOMS=1,2");
plumed.cmd("readInputLine","d1: DISTANCE ATOMS={1 2}");
///TODO: expand with a "readInputLines" to give the possibility to test in more situations
try {
//WHEN the user ask for the given input
plumed.cmd("readInputLine",cmd.c_str());
//THEN plumed should gracefully exit with a clear error message
} catch(PLMD::Plumed::ExceptionError &e) { //correct throw, we are happy
std::string exceptionText{e.what()};
out << name << " : ";
if (exceptionText.find(expectedMessage) != std::string::npos) {
out << "Exception thrown, correct message\n";
return 0;
}

out << "Exception thrown, wrong message: "
<< e.what() ;
out << "\tExpected message should contain: \""
<< expectedMessage << "\"\n";
return 1;
}
out << "Exception not thrown\n";
return 1;
}
};

int main(int, char**) {
std::ofstream out("output");
plumedThrowChecker ptc;
//When the user aks for a WALKERS_MPI in the METAD action,
//if MPI is not installed then the user must be informed
//if MPI is installed then the comunications must be already set up
//WHEN PLUMED is not compiled with MPI or the MPI routines are not initialized
std::string expectedMessage="WALKERS_MPI flag requires MPI compilation";
if (PLMD::Communicator::plumedHasMPI()) {
expectedMessage="WALKERS_MPI needs the communicator correctly initialized";
}
ptc.checkThrow(out,"METAD WALKER_MPI",
"METAD ARG=d,d1 SIGMA=0.1,0.2 HEIGHT=0.1 PACE=2 RESTART=YES WALKERS_MPI",
expectedMessage);

//add other throw messages checks
}
2 changes: 2 additions & 0 deletions regtest/scripts/run
Original file line number Diff line number Diff line change
Expand Up @@ -172,8 +172,10 @@ fi

if $plumed --is-installed ; then
export PLUMED_KERNEL="$root/../lib${PLUMED_PROGRAM_NAME:-plumed}Kernel.$($plumed info --soext)"
export PLUMED_LIB="$root/../lib${PLUMED_PROGRAM_NAME:-plumed}.$($plumed info --soext)"
else
export PLUMED_KERNEL="$root/src/lib/libplumedKernel.$($plumed info --soext)"
export PLUMED_LIB="$root/src/lib/libplumed.$($plumed info --soext)"
fi
if type -t plumed_regtest_before 1>/dev/null ; then
plumed_regtest_before
Expand Down
6 changes: 6 additions & 0 deletions src/bias/MetaD.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -878,6 +878,12 @@ MetaD::MetaD(const ActionOptions& ao):
// MPI version
parseFlag("WALKERS_MPI",walkers_mpi_);

//If this Action is not compiled with MPI the user is informed and we exit gracefully
if(walkers_mpi_) {
plumed_assert(Communicator::plumedHasMPI()) << "Invalid walkers configuration: WALKERS_MPI flag requires MPI compilation";
plumed_assert(Communicator::initialized()) << "Invalid walkers configuration: WALKERS_MPI needs the communicator correctly initialized.";
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Here I think you could just remove the first test.

Beside this, the way to raise the warnings is correct!!

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In my first draft i had a single error message with "MPI is not initialized or plumed is not compiled with MPI", but I think is more clear for the user to have two more specialized messages, and easier to test for the developers

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think you are right

}

// Flying Gaussian
parseFlag("FLYING_GAUSSIAN", flying_);

Expand Down
8 changes: 8 additions & 0 deletions src/tools/Communicator.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,14 @@

namespace PLMD {

bool Communicator::plumedHasMPI() {
#ifdef __PLUMED_HAS_MPI
return true;
#else
return false;
#endif
}

Communicator::Communicator()
#ifdef __PLUMED_HAS_MPI
: communicator(MPI_COMM_SELF)
Expand Down
3 changes: 3 additions & 0 deletions src/tools/Communicator.h
Original file line number Diff line number Diff line change
Expand Up @@ -121,6 +121,9 @@ class Communicator {
}
};
public:
///Runtime acces to the __PLUMED_HAS_MPI definition
static bool plumedHasMPI();

/// Wrapper class for MPI_Status
class Status {
int Get_count(MPI_Datatype)const;
Expand Down