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

Collective variable on the distance between two particles' centers of mass. More...

#include <ParticleSeparationCV.h>

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

Public Member Functions

 ParticleSeparationCV (const Label &group1, const Label &group2)
 Constructor. More...
 
 ParticleSeparationCV (const Label &group1, const Label &group2, bool dimx, bool dimy, bool dimz)
 Constructor. More...
 
void Initialize (const Snapshot &snapshot) override
 Initialize necessary variables. More...
 
void Evaluate (const Snapshot &snapshot) override
 Evaluate the CV. 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...
 
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 ParticleSeparationCVBuild (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_
 Atom ID's of group 1.
 
Label group2_
 Atom ID's of group 2.
 
Bool3 dim_
 Each dimension determines if it is applied by the 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

Collective variable on the distance between two particles' centers of mass.

Collective variable on two particle positions. This CV will return the distance between two specific atom groups of the simulation.

Definition at line 39 of file ParticleSeparationCV.h.

Constructor & Destructor Documentation

◆ ParticleSeparationCV() [1/2]

SSAGES::ParticleSeparationCV::ParticleSeparationCV ( const Label group1,
const Label group2 
)
inline

Constructor.

Parameters
group1Vector of atom ID's of the first particle.
group2Vector of atom ID's of the second particle.

Construct a paarticle separation CV.

Todo:
bounds needs to be an input.

Definition at line 58 of file ParticleSeparationCV.h.

58  :
59  group1_(group1), group2_(group2), dim_{true, true, true}
60  {}
Bool3 dim_
Each dimension determines if it is applied by the CV.
Label group1_
Atom ID's of group 1.
Label group2_
Atom ID's of group 2.

Referenced by Build().

Here is the caller graph for this function:

◆ ParticleSeparationCV() [2/2]

SSAGES::ParticleSeparationCV::ParticleSeparationCV ( const Label group1,
const Label group2,
bool  dimx,
bool  dimy,
bool  dimz 
)
inline

Constructor.

Parameters
group1Vector of atom ID's of the first particle.
group2Vector of atom ID's of the second particle.
dimxIf True, include x dimension.
dimyIf True, include y dimension.
dimzIf True, include z dimension.

Construct a particle separation CV.

Todo:
Bounds needs to be an input.

Definition at line 74 of file ParticleSeparationCV.h.

74  :
75  group1_(group1), group2_(group2), dim_{dimx, dimy, dimz}
76  {}

Member Function Documentation

◆ Build()

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

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, nullptr);
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  }
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
ParticleSeparationCV(const Label &group1, const Label &group2)
Constructor.

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

Evaluate the CV.

Parameters
snapshotCurrent simulation snapshot.

Definition at line 128 of file ParticleSeparationCV.h.

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);
148  Vector3 com2 = snapshot.CenterOfMass(idx2);
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  }
Matrix3 boxgrad_
Gradient w.r.t box vectors dCv/dHij.
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(), SSAGES::CollectiveVariable::boxgrad_, SSAGES::Snapshot::CenterOfMass(), dim_, SSAGES::Snapshot::GetLocalIndices(), SSAGES::Snapshot::GetMasses(), SSAGES::Snapshot::GetNumAtoms(), SSAGES::CollectiveVariable::grad_, group1_, group2_, SSAGES::Snapshot::TotalMass(), and SSAGES::CollectiveVariable::val_.

Here is the call graph for this function:

◆ Initialize()

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

Initialize necessary variables.

Parameters
snapshotCurrent simulation snapshot.

Definition at line 82 of file ParticleSeparationCV.h.

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  }

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: