SSAGES  0.8.3
Software Suite for Advanced General Ensemble Simulations
ParticleSeparationCV.h
1 
21 #pragma once
22 
23 #include "CollectiveVariable.h"
24 #include "Validator/ObjectRequirement.h"
25 #include "Drivers/DriverException.h"
26 #include "Snapshot.h"
27 #include "schema.h"
28 #include <numeric>
29 
30 namespace SSAGES
31 {
33 
40  {
41  private:
44 
47 
48  public:
50 
58  ParticleSeparationCV(const Label& group1, const Label& group2) :
59  group1_(group1), group2_(group2), dim_{true, true, true}
60  {}
61 
63 
74  ParticleSeparationCV(const Label& group1, const Label& group2, bool dimx, bool dimy, bool dimz) :
75  group1_(group1), group2_(group2), dim_{dimx, dimy, dimz}
76  {}
77 
79 
82  void Initialize(const Snapshot& snapshot) override
83  {
84  using std::to_string;
85 
86  auto n1 = group1_.size(), n2 = group2_.size();
87 
88  // Make sure atom ID's are on at least one processor.
89  std::vector<int> found1(n1, 0), found2(n2, 0);
90  for(size_t i = 0; i < n1; ++i)
91  {
92  if(snapshot.GetLocalIndex(group1_[i]) != -1)
93  found1[i] = 1;
94  }
95 
96  for(size_t i = 0; i < n2; ++i)
97  {
98  if(snapshot.GetLocalIndex(group2_[i]) != -1)
99  found2[i] = 1;
100  }
101 
102  MPI_Allreduce(MPI_IN_PLACE, found1.data(), n1, MPI_INT, MPI_SUM, snapshot.GetCommunicator());
103  MPI_Allreduce(MPI_IN_PLACE, found2.data(), n2, MPI_INT, MPI_SUM, snapshot.GetCommunicator());
104 
105  unsigned ntot1 = std::accumulate(found1.begin(), found1.end(), 0, std::plus<int>());
106  unsigned ntot2 = std::accumulate(found2.begin(), found2.end(), 0, std::plus<int>());
107  if(ntot1 != n1)
108  throw BuildException({
109  "ParticleSeparationCV: Expected to find " +
110  to_string(n1) +
111  " atoms in particle 1, but only found " +
112  to_string(ntot1) + "."
113  });
114 
115  if(ntot2 != n2)
116  throw BuildException({
117  "ParticleSeparationCV: Expected to find " +
118  to_string(n2) +
119  " atoms in particle 2, but only found " +
120  to_string(ntot2) + "."
121  });
122  }
123 
125 
128  void Evaluate(const Snapshot& snapshot) override
129  {
130  // Get local atom indices.
131  std::vector<int> idx1, idx2;
132  snapshot.GetLocalIndices(group1_, &idx1);
133  snapshot.GetLocalIndices(group2_, &idx2);
134 
135  // Get data from snapshot.
136  auto n = snapshot.GetNumAtoms();
137  const auto& masses = snapshot.GetMasses();
138 
139  // Initialize gradient.
140  std::fill(grad_.begin(), grad_.end(), Vector3{0,0,0});
141  grad_.resize(n, Vector3{0,0,0});
142  boxgrad_ = Matrix3::Zero();
143 
144  // Get centers of mass.
145  auto mtot1 = snapshot.TotalMass(idx1);
146  auto mtot2 = snapshot.TotalMass(idx2);
147  Vector3 com1 = snapshot.CenterOfMass(idx1, mtot1);
148  Vector3 com2 = snapshot.CenterOfMass(idx2, mtot2);
149 
150  // Account for pbc.
151  Vector3 rij = snapshot.ApplyMinimumImage(com1 - com2).cwiseProduct(dim_.cast<double>());
152 
153  // Compute gradient.
154  val_ = rij.norm();
155 
156  for(auto& id : idx1)
157  {
158  grad_[id] = rij/val_*masses[id]/mtot1;
159  boxgrad_ += grad_[id]*rij.transpose();
160  }
161 
162  for(auto& id : idx2)
163  grad_[id] = -rij/val_*masses[id]/mtot2;
164  }
165 
167  static ParticleSeparationCV* Build(const Json::Value& json, const std::string& path)
168  {
169  Json::ObjectRequirement validator;
170  Json::Value schema;
171  Json::CharReaderBuilder rbuilder;
172  Json::CharReader* reader = rbuilder.newCharReader();
173 
174  reader->parse(JsonSchema::ParticleSeparationCV.c_str(),
175  JsonSchema::ParticleSeparationCV.c_str() + JsonSchema::ParticleSeparationCV.size(),
176  &schema, NULL);
177  validator.Parse(schema, path);
178 
179  // Validate inputs.
180  validator.Validate(json, path);
181  if(validator.HasErrors())
182  throw BuildException(validator.GetErrors());
183 
184  std::vector<int> group1, group2;
185 
186  for(auto& s : json["group1"])
187  group1.push_back(s.asInt());
188 
189  for(auto& s : json["group2"])
190  group2.push_back(s.asInt());
191 
193  if(json.isMember("dimension"))
194  {
195  auto fixx = json["dimension"][0].asBool();
196  auto fixy = json["dimension"][1].asBool();
197  auto fixz = json["dimension"][2].asBool();
198  c = new ParticleSeparationCV(group1, group2, fixx, fixy, fixz);
199  }
200  else
201  c = new ParticleSeparationCV(group1, group2);
202 
203  return c;
204  }
205  };
206 }
std::vector< int > Label
List of integers.
Definition: types.h:48
bool HasErrors()
Check if errors have occured.
Definition: Requirement.h:86
static ParticleSeparationCV * Build(const Json::Value &json, const std::string &path)
Set up collective variable.
unsigned GetNumAtoms() const
Get number of atoms in this snapshot.
Definition: Snapshot.h:200
void ApplyMinimumImage(Vector3 *v) const
Apply minimum image to a vector.
Definition: Snapshot.h:418
const mxx::comm & GetCommunicator() const
Get communicator for walker.
Definition: Snapshot.h:184
ParticleSeparationCV(const Label &group1, const Label &group2, bool dimx, bool dimy, bool dimz)
Constructor.
Class containing a snapshot of the current simulation in time.
Definition: Snapshot.h:43
int GetLocalIndex(int id) const
Gets the local atom index corresponding to an atom ID.
Definition: Snapshot.h:555
double TotalMass(const Label &indices) const
Compute the total mass of a group of particles based on index.
Definition: Snapshot.h:451
void Initialize(const Snapshot &snapshot) override
Initialize necessary variables.
Label group2_
Atom ID&#39;s of group 2.
virtual void Parse(Value json, const std::string &path) override
Parse JSON value to generate Requirement(s).
Collective variable on the distance between two particles&#39; centers of mass.
std::vector< Vector3 > grad_
Gradient vector dCv/dxi.
Exception to be thrown when building the Driver fails.
std::vector< std::string > GetErrors()
Get list of error messages.
Definition: Requirement.h:92
ParticleSeparationCV(const Label &group1, const Label &group2)
Constructor.
Requirements on an object.
Vector3 CenterOfMass(const Label &indices) const
Definition: Snapshot.h:466
const std::vector< double > & GetMasses() const
Const access to the particle masses.
Definition: Snapshot.h:378
Eigen::Vector3d Vector3
Three-dimensional vector.
Definition: types.h:33
Bool3 dim_
Each dimension determines if it is applied by the CV.
double val_
Current value of CV.
Abstract class for a collective variable.
Matrix3 boxgrad_
Gradient w.r.t box vectors dCv/dHij.
void Evaluate(const Snapshot &snapshot) override
Evaluate the CV.
Label group1_
Atom ID&#39;s of group 1.
void GetLocalIndices(const Label &ids, Label *indices) const
Definition: Snapshot.h:573
virtual void Validate(const Value &json, const std::string &path) override
Validate JSON value.
Eigen::Matrix< bool, 3, 1 > Bool3
Three-dimensional boolean.
Definition: types.h:36