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::PairwiseCV Class Reference

Generalized collective variable based on pairwise properties of atoms. More...

#include <PairwiseCV.h>

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

Public Member Functions

 PairwiseCV (const Label &group1, const Label &group2, PairwiseKernel *pk)
 Constructor. More...
 
void Initialize (const Snapshot &snapshot) override
 Initialize necessary variables. More...
 
void Evaluate (const Snapshot &snapshot) override
 Evaluate the CV. More...
 
 ~PairwiseCV ()
 Destructor.
 
- 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...
 
virtual double GetMinimumImage (double) const
 Returns the minimum image of a CV based on the input location. More...
 
virtual double GetPeriodicValue (double location) const
 Apply periodic boundaries to a given value. 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...
 
virtual double GetDifference (double location) const
 

Static Public Member Functions

static PairwiseCVBuild (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 group1_
 IDs of the first group of atoms.
 
Label group2_
 IDs of the second group of atoms.
 
PairwiseKernelpk_
 Pairwise kernel function used for CV.
 

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

Generalized collective variable based on pairwise properties of atoms.

Collective variable on pairwise properties between two groups of atoms. To ensure generality of usage, there are various pairwise kernel functions from which to choose.

Definition at line 40 of file PairwiseCV.h.

Constructor & Destructor Documentation

◆ PairwiseCV()

SSAGES::PairwiseCV::PairwiseCV ( const Label group1,
const Label group2,
PairwiseKernel pk 
)
inline

Constructor.

Parameters
group1IDs of the first group of atoms.
group2IDs of the second group of atoms.
pkpairwise kernel (function) to use.

Construct a PairwiseCV.

Definition at line 57 of file PairwiseCV.h.

57  :
58  group1_(group1), group2_(group2), pk_(pk)
59  {
60  }
Label group1_
IDs of the first group of atoms.
Definition: PairwiseCV.h:43
Label group2_
IDs of the second group of atoms.
Definition: PairwiseCV.h:44
PairwiseKernel * pk_
Pairwise kernel function used for CV.
Definition: PairwiseCV.h:45

Referenced by Build().

Here is the caller graph for this function:

Member Function Documentation

◆ Build()

static PairwiseCV* SSAGES::PairwiseCV::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 200 of file PairwiseCV.h.

201  {
202  Json::ObjectRequirement validator;
203  Json::Value schema;
204  Json::CharReaderBuilder rbuilder;
205  Json::CharReader* reader = rbuilder.newCharReader();
206 
207  reader->parse(JsonSchema::PairwiseCV.c_str(),
208  JsonSchema::PairwiseCV.c_str() + JsonSchema::PairwiseCV.size(),
209  &schema, nullptr);
210  validator.Parse(schema, path);
211 
212  // Validate inputs.
213  validator.Validate(json, path);
214  if(validator.HasErrors())
215  throw BuildException(validator.GetErrors());
216 
217  std::vector<int> group1, group2;
218 
219  for(auto& s : json["group1"])
220  group1.push_back(s.asInt());
221 
222  for(auto& s : json["group2"])
223  group2.push_back(s.asInt());
224 
225  return new PairwiseCV(group1, group2, PairwiseKernel::Build(json["kernel"], path));
226  }
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
PairwiseCV(const Label &group1, const Label &group2, PairwiseKernel *pk)
Constructor.
Definition: PairwiseCV.h:57
static PairwiseKernel * Build(const Json::Value &json, const std::string &path)
Build PairwiseKernel from JSON value.

References SSAGES::PairwiseKernel::Build(), Json::Requirement::GetErrors(), Json::Requirement::HasErrors(), PairwiseCV(), Json::ObjectRequirement::Parse(), 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::PairwiseCV::Evaluate ( const Snapshot snapshot)
inlineoverride

Evaluate the CV.

Parameters
snapshotCurrent simulation snapshot.

Definition at line 112 of file PairwiseCV.h.

113  {
114  // Get local atom indices.
115  std::vector<int> idx1, idx2;
116  snapshot.GetLocalIndices(group1_, &idx1);
117  snapshot.GetLocalIndices(group2_, &idx2);
118 
119  // Get data from snapshot.
120  auto n = snapshot.GetNumAtoms();
121  auto& atomids = snapshot.GetAtomIDs();
122  auto& positions = snapshot.GetPositions();
123  auto& comm = snapshot.GetCommunicator();
124 
125  // Initialize gradient.
126  std::fill(grad_.begin(), grad_.end(), Vector3{0,0,0});
127  grad_.resize(n, Vector3{0,0,0});
128  boxgrad_ = Matrix3::Zero();
129 
130  // Fill eigen matrix with coordinate data.
131  std::vector<double> pos1(3*idx1.size()), pos2(3*idx2.size());
132  std::vector<int> id1(idx1.size()), id2(idx2.size());
133  for(size_t i = 0; i < idx1.size(); ++i)
134  {
135  pos1[3*i+0] = positions[idx1[i]][0];
136  pos1[3*i+1] = positions[idx1[i]][1];
137  pos1[3*i+2] = positions[idx1[i]][2];
138 
139  id1[i] = atomids[idx1[i]];
140  }
141 
142  for(size_t i = 0; i < idx2.size(); ++i)
143  {
144  pos2[3*i+0] = positions[idx2[i]][0];
145  pos2[3*i+1] = positions[idx2[i]][1];
146  pos2[3*i+2] = positions[idx2[i]][2];
147 
148  id2[i] = atomids[idx2[i]];
149  }
150 
151  // Gather across all procs.
152  pos1 = mxx::allgatherv(pos1.data(), pos1.size(), comm);
153  pos2 = mxx::allgatherv(pos2.data(), pos2.size(), comm);
154  id1 = mxx::allgatherv(id1.data(), id1.size(), comm);
155  id2 = mxx::allgatherv(id2.data(), id2.size(), comm);
156 
157  // Create Eigen map for ease of use.
158  using Map = Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, 3, Eigen::RowMajor>>;
159  Map gpos1(pos1.data(), group1_.size(), 3), gpos2(pos2.data(), group2_.size(), 3);
160 
161  val_ = 0;
162  double df = 0.;
163  // Compute gradient and value.
164  for(size_t i = 0; i < group1_.size(); ++i)
165  {
166  auto t1 = id1[i];
167  const auto& pi = gpos1.row(i);
168  for(size_t j = 0; j < group2_.size(); ++j)
169  {
170  auto t2 = id2[j];
171  const auto& pj = gpos2.row(j);
172 
173  // Skip identical pairs.
174  if(t1 == t2)
175  continue;
176 
177  // Compute distance and pairwise function.
178  Vector3 rij = pi - pj;
179  snapshot.ApplyMinimumImage(&rij);
180  auto r = rij.norm();
181  val_ += pk_->Evaluate(r, df);
182 
183  // Get local indices and sum gradients.
184  auto lid1 = snapshot.GetLocalIndex(t1);
185  if(lid1 != -1)
186  grad_[lid1] += df*rij/r;
187 
188  auto lid2 = snapshot.GetLocalIndex(t2);
189  if(lid2 != -1)
190  grad_[lid2] -= df*rij/r;
191 
192  // Only sum box gradient on a single proc.
193  if(comm.rank() == 0)
194  boxgrad_ += df*rij/r*rij.transpose();
195  }
196  }
197  }
Matrix3 boxgrad_
Gradient w.r.t box vectors dCv/dHij.
std::vector< Vector3 > grad_
Gradient vector dCv/dxi.
double val_
Current value of CV.
virtual double Evaluate(double rij, double &df) const =0
Evaluate the pairwise kernel function.
Eigen::Vector3d Vector3
Three-dimensional vector.
Definition: types.h:33

References SSAGES::Snapshot::ApplyMinimumImage(), SSAGES::CollectiveVariable::boxgrad_, SSAGES::PairwiseKernel::Evaluate(), SSAGES::Snapshot::GetAtomIDs(), SSAGES::Snapshot::GetCommunicator(), SSAGES::Snapshot::GetLocalIndex(), SSAGES::Snapshot::GetLocalIndices(), SSAGES::Snapshot::GetNumAtoms(), SSAGES::Snapshot::GetPositions(), SSAGES::CollectiveVariable::grad_, group1_, group2_, pk_, and SSAGES::CollectiveVariable::val_.

Here is the call graph for this function:

◆ Initialize()

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

Initialize necessary variables.

Parameters
snapshotCurrent simulation snapshot.

Definition at line 66 of file PairwiseCV.h.

67  {
68  using std::to_string;
69 
70  auto n1 = group1_.size(), n2 = group2_.size();
71 
72  // Make sure atom ID's are on at least one processor.
73  std::vector<int> found1(n1, 0), found2(n2, 0);
74  for(size_t i = 0; i < n1; ++i)
75  {
76  if(snapshot.GetLocalIndex(group1_[i]) != -1)
77  found1[i] = 1;
78  }
79 
80  for(size_t i = 0; i < n2; ++i)
81  {
82  if(snapshot.GetLocalIndex(group2_[i]) != -1)
83  found2[i] = 1;
84  }
85 
86  MPI_Allreduce(MPI_IN_PLACE, found1.data(), n1, MPI_INT, MPI_SUM, snapshot.GetCommunicator());
87  MPI_Allreduce(MPI_IN_PLACE, found2.data(), n2, MPI_INT, MPI_SUM, snapshot.GetCommunicator());
88 
89  unsigned ntot1 = std::accumulate(found1.begin(), found1.end(), 0, std::plus<int>());
90  unsigned ntot2 = std::accumulate(found2.begin(), found2.end(), 0, std::plus<int>());
91  if(ntot1 != n1)
92  throw BuildException({
93  "Pairwise: Expected to find " +
94  to_string(n1) +
95  " atoms in group 1, but only found " +
96  to_string(ntot1) + "."
97  });
98 
99  if(ntot2 != n2)
100  throw BuildException({
101  "Pairwise: Expected to find " +
102  to_string(n2) +
103  " atoms in group 2, but only found " +
104  to_string(ntot2) + "."
105  });
106  }

References SSAGES::Snapshot::GetCommunicator(), SSAGES::Snapshot::GetLocalIndex(), group1_, and group2_.

Here is the call graph for this function:

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