22 #ifndef __PLUMED_secondarystructure_SecondaryStructureRMSD_h
23 #define __PLUMED_secondarystructure_SecondaryStructureRMSD_h
25 #include "core/ActionAtomistic.h"
26 #include "core/ActionWithValue.h"
27 #include "vesselbase/ActionWithVessel.h"
32 class SingleDomainRMSD;
36 namespace secondarystructure {
60 std::vector< std::vector<Vector> >
der;
61 std::vector<Tensor>
vir;
72 std::vector<Vector>
pos;
78 void readBackboneAtoms(
const std::vector<std::string>& backnames, std::vector<unsigned>& chain_lengths );
80 void addColvar(
const std::vector<unsigned>& newatoms );
bool isPeriodic()
Are the base quantities periodic.
unsigned current
The numerical index of the task we are curently performing.
void setSecondaryStructure(std::vector< Vector > &structure, double bondlength, double units)
Set a reference configuration.
std::vector< double > forcesToApply
void apply()
Apply an Action.
bool pbcon
Are we using pbc.
A class for storing a list that changes which members are active as a function of time...
std::vector< Vector > pos
Tempory variables for getting positions of atoms and applying forces.
int updateFreq
Everything for controlling the updating of neighbor lists.
static void registerKeywords(Keywords &keys)
void readBackboneAtoms(const std::vector< std::string > &backnames, std::vector< unsigned > &chain_lengths)
Get the atoms in the backbone.
std::string alignType
The type of rmsd we are calculating.
void clearDerivativesAfterTask(const unsigned &)
std::vector< std::vector< Vector > > der
Stuff for derivatives.
void addColvar(const std::vector< unsigned > &newatoms)
Add a set of atoms to calculat ethe rmsd from.
unsigned getAtomIndex(const unsigned &iatom)
Get the index of an atom.
Used to create a PLMD::Action that has some scalar or vectorial output that may or may not have some ...
DynamicList< AtomNumber > all_atoms
List of all the atoms we require.
std::vector< Tensor > vir
This class holds the keywords and their documentation.
unsigned getNumberOfFunctionsInAction()
This class is used to bring the relevant information to the Action constructor.
Action used to create objects that access the positions of the atoms from the MD code.
unsigned getNumberOfDerivatives()
void calculate()
Calculate an Action.
SecondaryStructureRMSD(const ActionOptions &)
bool align_strands
Variables for strands cutoff.
void performTask(const unsigned &j)
Calculate one of the functions in the distribution.
void mergeDerivatives(const unsigned &, const double &)
unsigned getNumberOfAtoms() const
Get number of available atoms.
void setAtomsFromStrands(const unsigned &atom1, const unsigned &atom2)
Setup a pair of atoms to use for strands cutoff.
std::vector< RMSD * > secondary_rmsd
The reference configurations.
unsigned closest
Tempory integer to say which refernce configuration is the closest.
std::vector< std::vector< unsigned > > colvar_atoms
The atoms involved in each of the secondary structure segments.
std::vector< DRMSD * > secondary_drmsd
void prepare()
Prepare an Action for calculation This can be used by Action if they need some special preparation be...
This is used to create PLMD::Action objects that are computed by calculating the same function multip...
virtual ~SecondaryStructureRMSD()
Base action for calculating things like AlphRMSD, AntibetaRMSD, etc.