Skip to content

Commit

Permalink
Macro for fasterMulticolvar+named modes
Browse files Browse the repository at this point in the history
  • Loading branch information
@Iximiel
Iximiel committed Nov 19, 2024
1 parent c9df402 commit 8f1acce
Show file tree
Hide file tree
Showing 12 changed files with 239 additions and 153 deletions.
24 changes: 13 additions & 11 deletions src/colvar/Angle.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -104,15 +104,13 @@ class Angle : public Colvar {
std::vector<std::vector<Vector> > derivs;
std::vector<Tensor> virial;
public:
MULTICOLVAR_SETTINGS(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 Down Expand Up @@ -140,8 +138,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 +169,21 @@ Angle::Angle(const ActionOptions&ao):
void Angle::calculate() {

if(pbc) makeWhole();
calculateCV( 0, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( {}, masses, charges, getPositions(), 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,
void Angle::calculateCV( Modetype /*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 ) {
std::vector<Tensor>& virial, const ActionAtomistic* /*aa*/ ) {
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
19 changes: 9 additions & 10 deletions src/colvar/DihedralCorrelation.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -68,14 +68,11 @@ class DihedralCorrelation : public Colvar {
std::vector<std::vector<Vector> > derivs;
std::vector<Tensor> virial;
public:
MULTICOLVAR_SETTINGS(multiColvars::emptyMode);
static void registerKeywords( Keywords& keys );
explicit DihedralCorrelation(const ActionOptions&);
static void parseAtomList( const int& num, std::vector<AtomNumber>& t, ActionAtomistic* aa );
static unsigned getModeAndSetupValues( ActionWithValue* av );
void calculate() override;
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<DihedralCorrelation> DihedralCorrelationShortcut;
Expand Down Expand Up @@ -119,22 +116,24 @@ void DihedralCorrelation::parseAtomList( const int& num, std::vector<AtomNumber>
}
}

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

void DihedralCorrelation::calculate() {

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

void DihedralCorrelation::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 DihedralCorrelation::calculateCV(Modetype /*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 ) {
const Vector d10=delta(pos[1],pos[0]);
const Vector d11=delta(pos[2],pos[1]);
const Vector d12=delta(pos[3],pos[2]);
Expand Down
41 changes: 26 additions & 15 deletions src/colvar/Dipole.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -91,12 +91,19 @@ class Dipole : public Colvar {
public:
explicit Dipole(const ActionOptions&);
static void parseAtomList( const int& num, std::vector<AtomNumber>& t, ActionAtomistic* aa );
static unsigned getModeAndSetupValues( ActionWithValue* av );
void calculate() override;
static void registerKeywords(Keywords& keys);
static void calculateCV( const unsigned& mode, const std::vector<double>& masses, 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 );
MULTICOLVAR_SETTINGS_MODE(multiColvars::components);
MULTICOLVAR_SETTINGS_SETUPF();
//Declaring this playinly becasue we are using modifying the charges
static void calculateCV( Modetype mode,
const std::vector<double>& masses,
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<Dipole> DipoleShortcut;
Expand All @@ -123,8 +130,9 @@ Dipole::Dipole(const ActionOptions&ao):
derivs(1),
virial(1)
{
parseAtomList(-1,ga_lista,this); charges.resize(ga_lista.size());
components=(getModeAndSetupValues(this)==1);
parseAtomList(-1,ga_lista,this);
charges.resize(ga_lista.size());
components=(getModeAndSetupValues(this)==Modetype::withCompontents);
if( components ) {
value.resize(3); derivs.resize(3); virial.resize(3);
valuex=getPntrToComponent("x");
Expand Down Expand Up @@ -152,15 +160,17 @@ void Dipole::parseAtomList( const int& num, std::vector<AtomNumber>& t, ActionAt
}
}

unsigned Dipole::getModeAndSetupValues( ActionWithValue* av ) {
bool c; av->parseFlag("COMPONENTS",c);
Dipole::Modetype Dipole::getModeAndSetupValues( ActionWithValue* av ) {
bool c;
av->parseFlag("COMPONENTS",c);
if( c ) {
av->addComponentWithDerivatives("x"); av->componentIsNotPeriodic("x");
av->addComponentWithDerivatives("y"); av->componentIsNotPeriodic("y");
av->addComponentWithDerivatives("z"); av->componentIsNotPeriodic("z");
return 1;
return Modetype::withCompontents;
}
av->addValueWithDerivatives(); av->setNotPeriodic(); return 0;
av->addValueWithDerivatives(); av->setNotPeriodic();
return Modetype::noComponents;
}

// calculator
Expand All @@ -171,12 +181,12 @@ void Dipole::calculate()
for(unsigned i=0; i<N; ++i) charges[i]=getCharge(i);

if(!components) {
calculateCV( 0, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( Modetype::noComponents, masses, charges, getPositions(), value, derivs, virial, this );
for(unsigned i=0; i<N; i++) setAtomsDerivatives(i,derivs[0][i]);
setBoxDerivatives(virial[0]);
setValue(value[0]);
} else {
calculateCV( 1, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( Modetype::withCompontents, masses, charges, getPositions(), value, derivs, virial, this );
for(unsigned i=0; i<N; i++) {
setAtomsDerivatives(valuex,i,derivs[0][i]);
setAtomsDerivatives(valuey,i,derivs[1][i]);
Expand All @@ -191,7 +201,7 @@ void Dipole::calculate()
}
}

void Dipole::calculateCV( const unsigned& mode, const std::vector<double>& masses, std::vector<double>& charges,
void Dipole::calculateCV( Modetype mode, const std::vector<double>& masses, 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 ) {
unsigned N=pos.size(); double ctot=0.;
Expand All @@ -200,10 +210,11 @@ void Dipole::calculateCV( const unsigned& mode, const std::vector<double>& masse

Vector dipje;
for(unsigned i=0; i<N; ++i) {
charges[i]-=ctot; dipje += charges[i]*pos[i];
charges[i]-=ctot;
dipje += charges[i]*pos[i];
}

if( mode==1 ) {
if( mode==Modetype::withCompontents ) {
for(unsigned i=0; i<N; i++) {
derivs[0][i]=charges[i]*Vector(1.0,0.0,0.0);
derivs[1][i]=charges[i]*Vector(0.0,1.0,0.0);
Expand Down
51 changes: 29 additions & 22 deletions src/colvar/Distance.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -135,12 +135,9 @@ class Distance : public Colvar {
static void registerKeywords( Keywords& keys );
explicit Distance(const ActionOptions&);
static void parseAtomList( const int& num, std::vector<AtomNumber>& t, ActionAtomistic* aa );
static unsigned getModeAndSetupValues( ActionWithValue* av );
// active methods:
void calculate() override;
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 );
MULTICOLVAR_SETTINGS(multiColvars::scaledComponents);
};

typedef ColvarShortcut<Distance> DistanceShortcut;
Expand Down Expand Up @@ -186,11 +183,18 @@ Distance::Distance(const ActionOptions&ao):
if(pbc) log.printf(" using periodic boundary conditions\n");
else log.printf(" without periodic boundary conditions\n");

unsigned mode = getModeAndSetupValues( this );
if(mode==1) components=true; else if(mode==2) scaled_components=true;
Modetype mode = getModeAndSetupValues( this );
if(mode==Modetype::withCompontents) {
components=true;
} else if(mode==Modetype::scaledComponents) {
scaled_components=true;
}
if( components || scaled_components ) {
value.resize(3); derivs.resize(3); virial.resize(3);
for(unsigned i=0; i<3; ++i) derivs[i].resize(2);
value.resize(3);
derivs.resize(3);
virial.resize(3);
for(unsigned i=0; i<3; ++i)
derivs[i].resize(2);
}
requestAtoms(atoms);
}
Expand All @@ -200,25 +204,28 @@ void Distance::parseAtomList( const int& num, std::vector<AtomNumber>& t, Action
if( t.size()==2 ) aa->log.printf(" between atoms %d %d\n",t[0].serial(),t[1].serial());
}

unsigned Distance::getModeAndSetupValues( ActionWithValue* av ) {
bool c; av->parseFlag("COMPONENTS",c);
bool sc; av->parseFlag("SCALED_COMPONENTS",sc);
if( c && sc ) av->error("COMPONENTS and SCALED_COMPONENTS are not compatible");

Distance::Modetype Distance::getModeAndSetupValues( ActionWithValue* av ) {
bool c;
av->parseFlag("COMPONENTS",c);
bool sc;
av->parseFlag("SCALED_COMPONENTS",sc);
if( c && sc ) {
av->error("COMPONENTS and SCALED_COMPONENTS are not compatible");
}
if(c) {
av->addComponentWithDerivatives("x"); av->componentIsNotPeriodic("x");
av->addComponentWithDerivatives("y"); av->componentIsNotPeriodic("y");
av->addComponentWithDerivatives("z"); av->componentIsNotPeriodic("z");
av->log<<" WARNING: components will not have the proper periodicity - see manual\n";
return 1;
return Modetype::withCompontents;
} else if(sc) {
av->addComponentWithDerivatives("a"); av->componentIsPeriodic("a","-0.5","+0.5");
av->addComponentWithDerivatives("b"); av->componentIsPeriodic("b","-0.5","+0.5");
av->addComponentWithDerivatives("c"); av->componentIsPeriodic("c","-0.5","+0.5");
return 2;
return Modetype::scaledComponents;
}
av->addValueWithDerivatives(); av->setNotPeriodic();
return 0;
return Modetype::noComponents;
}

// calculator
Expand All @@ -227,7 +234,7 @@ void Distance::calculate() {
if(pbc) makeWhole();

if( components ) {
calculateCV( 1, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( Modetype::withCompontents, masses, charges, getPositions(), value, derivs, virial, this );
Value* valuex=getPntrToComponent("x");
Value* valuey=getPntrToComponent("y");
Value* valuez=getPntrToComponent("z");
Expand All @@ -244,7 +251,7 @@ void Distance::calculate() {
setBoxDerivatives(valuez,virial[2]);
valuez->set(value[2]);
} else if( scaled_components ) {
calculateCV( 2, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( Modetype::scaledComponents, masses, charges, getPositions(), value, derivs, virial, this );

Value* valuea=getPntrToComponent("a");
Value* valueb=getPntrToComponent("b");
Expand All @@ -256,21 +263,21 @@ void Distance::calculate() {
for(unsigned i=0; i<2; ++i) setAtomsDerivatives(valuec,i,derivs[2][i] );
valuec->set(value[2]);
} else {
calculateCV( 0, masses, charges, getPositions(), value, derivs, virial, this );
calculateCV( Modetype::noComponents, masses, charges, getPositions(), value, derivs, virial, this );
for(unsigned i=0; i<2; ++i) setAtomsDerivatives(i,derivs[0][i] );
setBoxDerivatives(virial[0]);
setValue (value[0]);
}
}

void Distance::calculateCV( const unsigned& mode, const std::vector<double>& masses, const std::vector<double>& charges,
void Distance::calculateCV( Modetype 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 ) {
Vector distance=delta(pos[0],pos[1]);
const double value=distance.modulo();
const double invvalue=1.0/value;

if(mode==1) {
if(mode==Modetype::withCompontents) {
derivs[0][0] = Vector(-1,0,0);
derivs[0][1] = Vector(+1,0,0);
vals[0] = distance[0];
Expand All @@ -283,7 +290,7 @@ void Distance::calculateCV( const unsigned& mode, const std::vector<double>& mas
derivs[2][1] = Vector(0,0,+1);
vals[2] = distance[2];
setBoxDerivativesNoPbc( pos, derivs, virial );
} else if(mode==2) {
} else if(mode==Modetype::scaledComponents) {
Vector d=aa->getPbc().realToScaled(distance);
derivs[0][0] = matmul(aa->getPbc().getInvBox(),Vector(-1,0,0));
derivs[0][1] = matmul(aa->getPbc().getInvBox(),Vector(+1,0,0));
Expand Down
Loading

0 comments on commit 8f1acce

Please sign in to comment.