SSAGES  0.9.3
Software Suite for Advanced General Ensemble Simulations
Public Member Functions | Static Public Member Functions | Private Attributes | List of all members
SSAGES::TorsionalCV Class Reference

Collective variable on the torsion angles. More...

#include <TorsionalCV.h>

Inheritance diagram for SSAGES::TorsionalCV:
Inheritance graph
[legend]

Public Member Functions

 TorsionalCV (int atomid1, int atomid2, int atomid3, int atomid4, bool periodic)
 Constructor. More...
 
void Initialize (const Snapshot &snapshot) override
 Initialize necessary variables. More...
 
void Evaluate (const Snapshot &snapshot) override
 Evaluate the CV. More...
 
double GetPeriodicValue (double Location) const override
 Return value taking periodic boundary conditions into account. More...
 
double GetDifference (const double Location) const override
 Get difference taking periodic boundary conditions into account. More...
 
double GetMinimumImage (const double Location) const override
 Returns the minimum image of a CV based on the input location. More...
 
- Public Member Functions inherited from SSAGES::CollectiveVariable
 CollectiveVariable ()
 Constructor.
 
virtual ~CollectiveVariable ()
 Destructor.
 
virtual void Initialize (const class Snapshot &)
 Initialize CV. More...
 
virtual void Evaluate (const class Snapshot &)=0
 Evaluate CV. More...
 
double GetValue () const
 Get current value of the CV. More...
 
const std::vector< Vector3 > & GetGradient () const
 Get current gradient of the CV. More...
 
const Matrix3GetBoxGradient () const
 Get gradient contribution to box. More...
 
const std::array< double, 2 > & GetBoundaries ()
 Get CV boundaries. More...
 

Static Public Member Functions

static TorsionalCVBuild (const Json::Value &json, const std::string &path)
 Set up collective variable. More...
 
- Static Public Member Functions inherited from SSAGES::CollectiveVariable
static CollectiveVariableBuildCV (const Json::Value &json, const std::string &path)
 Set up collective variable. More...
 

Private Attributes

Label atomids_
 Vector of 4 atom ID's of interest.
 
bool periodic_
 If True, use periodic boundary conditions.
 

Additional Inherited Members

- Protected Attributes inherited from SSAGES::CollectiveVariable
std::vector< Vector3grad_
 Gradient vector dCv/dxi.
 
Matrix3 boxgrad_
 Gradient w.r.t box vectors dCv/dHij.
 
double val_
 Current value of CV.
 
std::array< double, 2 > bounds_
 Bounds on CV.
 

Detailed Description

Collective variable on the torsion angles.

Collective variable on an proper dihedral. This will return the angle between two planes as defined in [4]. Singularities are avoided as described in [1].

Definition at line 42 of file TorsionalCV.h.

Constructor & Destructor Documentation

◆ TorsionalCV()

SSAGES::TorsionalCV::TorsionalCV ( int  atomid1,
int  atomid2,
int  atomid3,
int  atomid4,
bool  periodic 
)
inline

Constructor.

Parameters
atomid1ID of the first atom defining the dihedral angle.
atomid2ID of the second atom defining the dihedral angle.
atomid3ID of the third atom defining the dihedral angle.
atomid4ID of the forth atom defining the dihedral angle.
periodicIf True consider periodic boundary conditions.

Construct an dihedral CV.

Todo:
Bounds needs to be an input and periodic boundary conditions

Definition at line 64 of file TorsionalCV.h.

64  :
65  atomids_({atomid1, atomid2, atomid3, atomid4}), periodic_(periodic)
66  {
67  }
Label atomids_
Vector of 4 atom ID's of interest.
Definition: TorsionalCV.h:46
bool periodic_
If True, use periodic boundary conditions.
Definition: TorsionalCV.h:49

Referenced by Build().

Here is the caller graph for this function:

Member Function Documentation

◆ Build()

static TorsionalCV* SSAGES::TorsionalCV::Build ( const Json::Value &  json,
const std::string &  path 
)
inlinestatic

Set up collective variable.

Parameters
jsonJSON input value.
pathPath for JSON path specification.
Returns
Pointer to the CV built by this function. nullptr in case of unknown error.

Builds a CV from a JSON node. Returns a pointer to the built cv. If an unknown error is encountered, this function will return a nullptr, but generally it will throw a BuildException on failure.

Warning
Object lifetime is the caller's responsibility.

Definition at line 229 of file TorsionalCV.h.

230  {
231  Json::ObjectRequirement validator;
232  Json::Value schema;
233  Json::CharReaderBuilder rbuilder;
234  Json::CharReader* reader = rbuilder.newCharReader();
235 
236  reader->parse(JsonSchema::TorsionalCV.c_str(),
237  JsonSchema::TorsionalCV.c_str() + JsonSchema::TorsionalCV.size(),
238  &schema, nullptr);
239  validator.Parse(schema, path);
240 
241  // Validate inputs.
242  validator.Validate(json, path);
243  if(validator.HasErrors())
244  throw BuildException(validator.GetErrors());
245 
246  std::vector<int> atomids;
247  for(auto& s : json["atom_ids"])
248  atomids.push_back(s.asInt());
249 
250  auto periodic = json.get("periodic", true).asBool();
251  return new TorsionalCV(atomids[0], atomids[1], atomids[2], atomids[3], periodic);
252  }
Requirements on an object.
virtual void Parse(Value json, const std::string &path) override
Parse JSON value to generate Requirement(s).
virtual void Validate(const Value &json, const std::string &path) override
Validate JSON value.
std::vector< std::string > GetErrors()
Get list of error messages.
Definition: Requirement.h:92
bool HasErrors()
Check if errors have occured.
Definition: Requirement.h:86
TorsionalCV(int atomid1, int atomid2, int atomid3, int atomid4, bool periodic)
Constructor.
Definition: TorsionalCV.h:64

References Json::Requirement::GetErrors(), Json::Requirement::HasErrors(), Json::ObjectRequirement::Parse(), TorsionalCV(), and Json::ObjectRequirement::Validate().

Referenced by SSAGES::CollectiveVariable::BuildCV().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ Evaluate()

void SSAGES::TorsionalCV::Evaluate ( const Snapshot snapshot)
inlineoverride

Evaluate the CV.

Parameters
snapshotCurrent simulation snapshot.

Definition at line 95 of file TorsionalCV.h.

96  {
97  // Get data from snapshot.
98  auto n = snapshot.GetNumAtoms();
99  const auto& pos = snapshot.GetPositions();
100  auto& comm = snapshot.GetCommunicator();
101 
102  // Initialize gradient.
103  std::fill(grad_.begin(), grad_.end(), Vector3{0,0,0});
104  grad_.resize(n, Vector3{0,0,0});
105 
106  Vector3 xi{0, 0, 0}, xj{0, 0, 0}, xk{0, 0, 0}, xl{0, 0, 0};
107 
108  auto iindex = snapshot.GetLocalIndex(atomids_[0]);
109  auto jindex = snapshot.GetLocalIndex(atomids_[1]);
110  auto kindex = snapshot.GetLocalIndex(atomids_[2]);
111  auto lindex = snapshot.GetLocalIndex(atomids_[3]);
112 
113  if(iindex != -1) xi = pos[iindex];
114  if(jindex != -1) xj = pos[jindex];
115  if(kindex != -1) xk = pos[kindex];
116  if(lindex != -1) xl = pos[lindex];
117 
118  // By performing a reduce, we actually collect all. This can
119  // be converted to a more intelligent allgater on rank then bcast.
120  MPI_Allreduce(MPI_IN_PLACE, xi.data(), 3, MPI_DOUBLE, MPI_SUM, comm);
121  MPI_Allreduce(MPI_IN_PLACE, xj.data(), 3, MPI_DOUBLE, MPI_SUM, comm);
122  MPI_Allreduce(MPI_IN_PLACE, xk.data(), 3, MPI_DOUBLE, MPI_SUM, comm);
123  MPI_Allreduce(MPI_IN_PLACE, xl.data(), 3, MPI_DOUBLE, MPI_SUM, comm);
124 
125  //Calculate pertinent vectors
126  auto F = snapshot.ApplyMinimumImage(xi - xj);
127  auto G = snapshot.ApplyMinimumImage(xj - xk);
128  auto H = snapshot.ApplyMinimumImage(xl - xk);
129  auto A = F.cross(G);
130  auto B = H.cross(G);
131 
132  auto y = B.cross(A).dot(G.normalized());
133  auto x = A.dot(B);
134 
135  val_ = atan2(y, x);
136 
137  auto Zed = F.dot(G.normalized())/A.dot(A);
138  auto Ned = H.dot(G.normalized())/B.dot(B);
139 
140  Vector3 gradi{0,0,0}, gradl{0,0,0};
141  if(iindex != -1) gradi = -G.norm()*A/A.dot(A);
142  if(lindex != -1) gradl = G.norm()*B/B.dot(B);
143  MPI_Allreduce(MPI_IN_PLACE, gradi.data(), 3, MPI_DOUBLE, MPI_SUM, comm);
144  MPI_Allreduce(MPI_IN_PLACE, gradl.data(), 3, MPI_DOUBLE, MPI_SUM, comm);
145  if(iindex != -1) grad_[iindex] = gradi;
146  if(lindex != -1) grad_[lindex] = gradl;
147  if(jindex != -1) grad_[jindex] = Zed*A - Ned*B - gradi;
148  if(kindex != -1) grad_[kindex] = Ned*B - Zed*A - gradl;
149  }
std::vector< Vector3 > grad_
Gradient vector dCv/dxi.
double val_
Current value of CV.
Eigen::Vector3d Vector3
Three-dimensional vector.
Definition: types.h:33

References SSAGES::Snapshot::ApplyMinimumImage(), atomids_, SSAGES::Snapshot::GetCommunicator(), SSAGES::Snapshot::GetLocalIndex(), SSAGES::Snapshot::GetNumAtoms(), SSAGES::Snapshot::GetPositions(), SSAGES::CollectiveVariable::grad_, and SSAGES::CollectiveVariable::val_.

Here is the call graph for this function:

◆ GetDifference()

double SSAGES::TorsionalCV::GetDifference ( const double  Location) const
inlineoverridevirtual

Get difference taking periodic boundary conditions into account.

Parameters
LocationInput value
Returns
Wrapped or unwrapped difference depending on whether periodic boundaries are used.

If periodic boundaries are used, this function calculates the difference and wraps the result into the range (-pi, pi). Otherwise, the simple difference is returned.

Reimplemented from SSAGES::CollectiveVariable.

Definition at line 188 of file TorsionalCV.h.

189  {
190  double PeriodicDiff = val_ - Location;
191 
192  if(!periodic_)
193  return PeriodicDiff;
194 
195  PeriodicDiff = GetPeriodicValue(PeriodicDiff);
196 
197  if(PeriodicDiff > M_PI)
198  PeriodicDiff -= 2.0*M_PI;
199  else if(PeriodicDiff < -M_PI)
200  PeriodicDiff += 2.0*M_PI;
201 
202  return PeriodicDiff;
203  }
double GetPeriodicValue(double Location) const override
Return value taking periodic boundary conditions into account.
Definition: TorsionalCV.h:161

References GetPeriodicValue(), periodic_, and SSAGES::CollectiveVariable::val_.

Here is the call graph for this function:

◆ GetMinimumImage()

double SSAGES::TorsionalCV::GetMinimumImage ( const double  Location) const
inlineoverridevirtual

Returns the minimum image of a CV based on the input location.

Parameters
LocationValue against which the minimum image is calculated.
Returns
Minimum image of the CV

Takes the input location and applies the periodic boundary conditions to return a minimum image of the CV.

Reimplemented from SSAGES::CollectiveVariable.

Definition at line 213 of file TorsionalCV.h.

214  {
215  double PeriodicDiff = val_ - Location;
216 
217  if(!periodic_)
218  return val_;
219 
220  if(PeriodicDiff > M_PI)
221  return (val_ - 2.0*M_PI);
222  else if(PeriodicDiff < -M_PI)
223  return (val_ + 2.0*M_PI);
224 
225  return val_;
226  }

References periodic_, and SSAGES::CollectiveVariable::val_.

◆ GetPeriodicValue()

double SSAGES::TorsionalCV::GetPeriodicValue ( double  Location) const
inlineoverridevirtual

Return value taking periodic boundary conditions into account.

Parameters
LocationInput value.
Returns
Wrapped or unwrapped input value depending on whether periodic boundaries are used.

If periodic boundaries are used, this function wraps the input value into the range (-pi, pi). Otherwise the input value is returned unmodified.

Reimplemented from SSAGES::CollectiveVariable.

Definition at line 161 of file TorsionalCV.h.

162  {
163  if(!periodic_)
164  return Location;
165 
166  int n = (int)(Location/(2.0*M_PI));
167  double PeriodicLocation;
168 
169  PeriodicLocation = Location - n*M_PI;
170  if(PeriodicLocation < -M_PI)
171  PeriodicLocation += 2.0*M_PI;
172  else if (PeriodicLocation > M_PI)
173  PeriodicLocation -= 2.0*M_PI;
174 
175  return PeriodicLocation;
176  }

References periodic_.

Referenced by GetDifference().

Here is the caller graph for this function:

◆ Initialize()

void SSAGES::TorsionalCV::Initialize ( const Snapshot snapshot)
inlineoverride

Initialize necessary variables.

Parameters
snapshotCurrent simulation snapshot.

Definition at line 73 of file TorsionalCV.h.

74  {
75  using std::to_string;
76 
77  std::vector<int> found;
78  snapshot.GetLocalIndices(atomids_, &found);
79  int nfound = found.size();
80  MPI_Allreduce(MPI_IN_PLACE, &nfound, 1, MPI_INT, MPI_SUM, snapshot.GetCommunicator());
81 
82  if(nfound != 4)
83  throw BuildException({
84  "TorsionalCV: Expected to find " +
85  to_string(4) +
86  " atoms, but only found " +
87  to_string(nfound) + "."
88  });
89  }

References atomids_, SSAGES::Snapshot::GetCommunicator(), and SSAGES::Snapshot::GetLocalIndices().

Here is the call graph for this function:

The documentation for this class was generated from the following file: