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

Macro for fasterMulticolvar+named modes #1158

Draft
wants to merge 10 commits into
base: derivatives-from-backpropegation
Choose a base branch
from
2 changes: 1 addition & 1 deletion regtest/crystdistrib/rt-dops-shortcut/config
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
type=driver
plumed_modules=crystdistrib
plumed_modules="crystdistrib adjmat"
arg="--plumed plumed.dat --mf_xtc em.xtc --mc mcfile --timestep=0.001 --trajectory-stride=25"
extra_files="../../trajectories/em.xtc"
1 change: 1 addition & 0 deletions regtest/multicolvar/rt-input/Makefile
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
include ../../scripts/test.make
1 change: 1 addition & 0 deletions regtest/multicolvar/rt-input/config
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
type=make
201 changes: 201 additions & 0 deletions regtest/multicolvar/rt-input/main.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,201 @@
#include "plumed/colvar/MultiColvarTemplate.h"
#include "plumed/tools/Vector.h"

#include <fstream>
#include <iostream>
#include <variant>
#include <vector>

using namespace PLMD;
using colvar::multiColvars::Input;

template <typename T> void pv(std::ofstream &stream, const std::vector<T> &v) {
// todo: put this is a common header
std::string div = "";
stream << "(" << v.size() << ") ";
for (const auto &t : v) {
stream << div << t;
div = ", ";
}
}

void print(std::ofstream &stream, const std::vector<double> &masses,
const std::vector<double> &charges,
const std::vector<Vector> positions) {
stream << "Masses: ";
pv(stream, masses);
stream << "\nCharges: ";
pv(stream, charges);
stream << "\nPositions: ";
pv(stream, positions);
stream << "\n";
}

void inputTest(std::ofstream &stream, Input input) {
try {
const auto &masses = input.masses();
stream << " Masses: ";
pv(stream, masses);
stream << "\n";
} catch (std::bad_variant_access &) {
stream << " There are no masses\n";
}

try {
const auto &charges = input.charges();
stream << " Charges: ";
pv(stream, charges);
stream << "\n";
} catch (std::bad_variant_access &) {
stream << " There are no charges\n";
}
try {
const auto &positions = input.positions();
stream << " Positions: ";
pv(stream, positions);
stream << "\n";
} catch (std::bad_variant_access &) {
stream << " There are no positions\n";
}
stream << std::endl;
}

void workOnMass(Input input) {
auto &masses = input.var_masses();
for (auto &m : masses) {
m *= 2;
}
}

void workOnCharges(Input input) {
auto &charges = input.var_charges();
for (auto &c : charges) {
c *= 2;
}
}

void workOnPositions(Input input) {
auto &pos = input.var_positions();
for (auto &p : pos) {
p = p * 2.0;
}
}

class testData {
using vd = std::vector<double>;
using vv = std::vector<Vector>;
vv p = {{0, 0, 0}, {1, 1, 1}};
vd c = {-1, 1};
vd m = {3, 4};

public:
testData() = default;
vd &masses() { return m; }
vd &charges() { return c; }
vv &positions() { return p; }
const vd &c_masses() const { return m; }
const vd &c_charges() const { return c; }
const vv &c_positions() const { return p; }
};

class test {
using vd = std::vector<double>;
using vv = std::vector<Vector>;

public:
test() = default;
void dotestWithReference(std::ofstream &output) {
testData td;
output << "*dotestWithReference:" << std::endl;
output << "*Passing full input:\n";
inputTest(output, Input(td.masses(), td.charges(), td.positions()));
output << "*Passing only positions:\n";
inputTest(output, Input().positions(td.positions()));
output << "*Passing only masses:\n";
inputTest(output, Input().masses(td.masses()));
output << "*Passing only charges:\n";
inputTest(output, Input().charges(td.charges()));
}
void dotestWithConstReference(std::ofstream &output) {
testData td;
output << "*dotestWithConstReference:" << std::endl;
output << "*Passing full input:\n";
inputTest(output, Input(td.c_masses(), td.c_charges(), td.c_positions()));
output << "*Passing only positions:\n";
inputTest(output, Input().positions(td.c_positions()));
output << "*Passing only masses:\n";
inputTest(output, Input().masses(td.c_masses()));
output << "*Passing only charges:\n";
inputTest(output, Input().charges(td.c_charges()));
}

void testUsingInputAsBuffer_masses(std::ofstream &output) {
testData td;
output << "*testUsingInputAsBuffer_masses:" << std::endl;
output << "*Passing constant_reference:\n";
try {
workOnMass(Input(td.c_masses(), td.c_charges(), td.c_positions()));
} catch (std::bad_variant_access &) {
output << " should throw\n";
}
output << "*Passing reference:\n";
output << " Initial values\n Masses:";
pv(output, td.masses());
workOnMass(Input(td.masses(), td.charges(), td.positions()));
output << "\n";
// this extra line is to remember that something has changed!!!
output << " Values after call\n Masses:";
pv(output, td.masses());
output << std::endl;
}
void testUsingInputAsBuffer_charges(std::ofstream &output) {
testData td;
output << "*testUsingInputAsBuffer_charges:" << std::endl;
output << "*Passing constant_reference:\n";
try {
workOnCharges(Input(td.c_masses(), td.c_charges(), td.c_positions()));
} catch (std::bad_variant_access &) {
output << " should throw\n";
}
output << "*Passing reference:\n";
output << " Initial values\n Charges:";
pv(output, td.c_charges());
workOnCharges(Input(td.masses(), td.charges(), td.positions()));
output << "\n";
// this extra line is to remember that something has changed!!!
output << " Values after call\n Charges:";
pv(output, td.c_charges());
output << std::endl;
}
void testUsingInputAsBuffer_positions(std::ofstream &output) {
testData td;
output << "*testUsingInputAsBuffer_positions:" << std::endl;
output << "*Passing constant_reference:\n";
try {
workOnPositions(Input(td.c_masses(), td.c_charges(), td.c_positions()));
} catch (std::bad_variant_access &) {
output << " should throw\n";
}
output << "*Passing reference:\n";
output << " Initial values\n Positions:";
pv(output, td.positions());
workOnPositions(Input(td.masses(), td.charges(), td.positions()));
output << "\n";
// this extra line is to remember that something has changed!!!
output << " Values after call\n Positions:";
pv(output, td.positions());
output << std::endl;
}
};

int main() {
std::ofstream output("output");

test t;
t.dotestWithReference(output);
t.dotestWithConstReference(output);
t.testUsingInputAsBuffer_masses(output);
t.testUsingInputAsBuffer_charges(output);
t.testUsingInputAsBuffer_positions(output);
return 0;
}
66 changes: 66 additions & 0 deletions regtest/multicolvar/rt-input/output.reference
Original file line number Diff line number Diff line change
@@ -0,0 +1,66 @@
*dotestWithReference:
*Passing full input:
Masses: (2) 3, 4
Charges: (2) -1, 1
Positions: (2) 0 0 0, 1 1 1

*Passing only positions:
There are no masses
There are no charges
Positions: (2) 0 0 0, 1 1 1

*Passing only masses:
Masses: (2) 3, 4
There are no charges
There are no positions

*Passing only charges:
There are no masses
Charges: (2) -1, 1
There are no positions

*dotestWithConstReference:
*Passing full input:
Masses: (2) 3, 4
Charges: (2) -1, 1
Positions: (2) 0 0 0, 1 1 1

*Passing only positions:
There are no masses
There are no charges
Positions: (2) 0 0 0, 1 1 1

*Passing only masses:
Masses: (2) 3, 4
There are no charges
There are no positions

*Passing only charges:
There are no masses
Charges: (2) -1, 1
There are no positions

*testUsingInputAsBuffer_masses:
*Passing constant_reference:
should throw
*Passing reference:
Initial values
Masses:(2) 3, 4
Values after call
Masses:(2) 6, 8
*testUsingInputAsBuffer_charges:
*Passing constant_reference:
should throw
*Passing reference:
Initial values
Charges:(2) -1, 1
Values after call
Charges:(2) -2, 2
*testUsingInputAsBuffer_positions:
*Passing constant_reference:
should throw
*Passing reference:
Initial values
Positions:(2) 0 0 0, 1 1 1
Values after call
Positions:(2) 0 0 0, 2 2 2
36 changes: 21 additions & 15 deletions src/colvar/Angle.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -100,19 +100,16 @@ Calculate multiple angles.

class Angle : public Colvar {
bool pbc;
std::vector<double> value, masses, charges;
std::vector<double> value;
std::vector<std::vector<Vector> > derivs;
std::vector<Tensor> virial;
public:
MULTICOLVAR_DEFAULT(multiColvars::emptyMode);
explicit Angle(const ActionOptions&);
// active methods:
void calculate() override;
static void registerKeywords( Keywords& keys );
static void parseAtomList( const int& num, std::vector<AtomNumber>& t, ActionAtomistic* aa );
static unsigned getModeAndSetupValues( ActionWithValue* av );
static void calculateCV( const unsigned& mode, const std::vector<double>& masses, const std::vector<double>& charges,
const std::vector<Vector>& pos, std::vector<double>& vals, std::vector<std::vector<Vector> >& derivs,
std::vector<Tensor>& virial, const ActionAtomistic* aa );

};

typedef ColvarShortcut<Angle> AngleShortcut;
Expand All @@ -128,7 +125,7 @@ void Angle::registerKeywords( Keywords& keys ) {
keys.setValueDescription("the ANGLE involving these atoms");
}

void Angle::parseAtomList( const int& num, std::vector<AtomNumber>& atoms, ActionAtomistic* aa ) {
void Angle::parseAtomList( int const num, std::vector<AtomNumber>& atoms, ActionAtomistic* aa ) {
aa->parseAtomList("ATOMS",num,atoms);
if(atoms.size()==3) {
aa->log.printf(" between atoms %d %d %d\n",atoms[0].serial(),atoms[1].serial(),atoms[2].serial());
Expand All @@ -140,8 +137,10 @@ void Angle::parseAtomList( const int& num, std::vector<AtomNumber>& atoms, Actio
} else if( num<0 || atoms.size()>0 ) aa->error("Number of specified atoms should be either 3 or 4");
}

unsigned Angle::getModeAndSetupValues( ActionWithValue* av ) {
av->addValueWithDerivatives(); av->setNotPeriodic(); return 0;
Angle::Modetype Angle::getModeAndSetupValues( ActionWithValue* av ) {
av->addValueWithDerivatives();
av->setNotPeriodic();
return {};
}

Angle::Angle(const ActionOptions&ao):
Expand Down Expand Up @@ -169,19 +168,26 @@ Angle::Angle(const ActionOptions&ao):
void Angle::calculate() {

if(pbc) makeWhole();
calculateCV( 0, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( {}, multiColvars::Input().positions(getPositions()), multiColvars::Ouput(value, derivs, virial), this );
setValue( value[0] );
for(unsigned i=0; i<derivs[0].size(); ++i) setAtomsDerivatives( i, derivs[0][i] );
for(unsigned i=0; i<derivs[0].size(); ++i)
setAtomsDerivatives( i, derivs[0][i] );
setBoxDerivatives( virial[0] );
}

void Angle::calculateCV( const unsigned& mode, const std::vector<double>& masses, const std::vector<double>& charges,
const std::vector<Vector>& pos, std::vector<double>& vals, std::vector<std::vector<Vector> >& derivs,
std::vector<Tensor>& virial, const ActionAtomistic* aa ) {
void Angle::calculateCV( Modetype /*mode*/,
multiColvars::Input const in,
multiColvars::Ouput out,
const ActionAtomistic* aa ) {
auto & vals=out.vals();
auto & derivs=out.derivs();
auto & virial=out.virial();
const auto & pos=in.positions();
Vector dij,dik;
dij=delta(pos[2],pos[3]);
dik=delta(pos[1],pos[0]);
Vector ddij,ddik; PLMD::Angle a;
Vector ddij,ddik;
PLMD::Angle a;
vals[0]=a.compute(dij,dik,ddij,ddik);
derivs[0][0]=ddik; derivs[0][1]=-ddik;
derivs[0][2]=-ddij; derivs[0][3]=ddij;
Expand Down
Loading
Loading