SSAGES  0.8.5
Software Suite for Advanced General Ensemble Simulations
Snapshot.h
1 
21 #pragma once
22 #include <numeric>
23 #include <vector>
24 #include <memory>
25 #include <mxx/comm.hpp>
26 #include "types.h"
27 
28 namespace SSAGES
29 {
31 
35  inline double roundf(double x)
36  {
37  return ( x >= 0 ) ? floor( x + 0.5 ) : ceil( x - 0.5 );
38  }
39 
41 
47  class Snapshot
48  {
49  private:
51  mxx::comm comm_;
52 
53  unsigned int wid_;
54  unsigned int nlocal_;
55 
56  std::string ID_;
57 
60 
62 
64 
66 
67  std::vector<Vector3> positions_;
68  std::vector<Integer3> images_;
69  std::vector<Vector3> velocities_;
70  std::vector<Vector3> forces_;
71  std::vector<double> masses_;
72  std::vector<double> charges_;
75  std::vector<std::vector<double> > sigma_;
76 
77  size_t iteration_;
78  size_t targetiter_;
79  double temperature_;
80  double energy_;
81  double kb_;
82 
83  double dielectric_;
84  double qqrd2e_;
85  bool changed_;
86 
87  public:
89 
96  Snapshot(const MPI_Comm& comm, unsigned int wid) :
97  comm_(comm), wid_(wid), H_(), Hinv_(), virial_(Matrix3::Zero()),
98  origin_({0,0,0}), isperiodic_({true, true, true}), positions_(0),
99  images_(0), velocities_(0), forces_(0), masses_(0), atomids_(0),
100  types_(0), iteration_(0), temperature_(0), energy_(0), kb_(0)
101  {}
102 
104 
107  size_t GetIteration() const {return iteration_; }
108 
110 
113  size_t GetTargetIterations() const {return targetiter_; }
114 
116 
119  double GetTemperature() const {return temperature_; }
120 
122 
125  double GetEnergy() const { return energy_; }
126 
128 
131  const Matrix3& GetHMatrix() const { return H_; }
132 
134 
137  const Matrix3& GetVirial() const { return virial_; }
138 
140 
143  Matrix3& GetVirial() { return virial_; }
144 
146 
149  const Vector3& GetOrigin() const { return origin_; }
150 
152 
155  const Bool3& IsPeriodic() const {return isperiodic_; }
156 
158 
161  double GetVolume() const { return H_.determinant(); }
162 
164 
167  double GetKb() const { return kb_; }
168 
170 
173  double GetDielectric() const { return dielectric_; }
174 
176 
179  double Getqqrd2e() const { return qqrd2e_; }
180 
182 
188  const mxx::comm& GetCommunicator() const
189  {
190  return comm_;
191  }
192 
194 
197  unsigned GetWalkerID() const { return wid_; }
198 
199 
201 
204  unsigned GetNumAtoms() const { return nlocal_; }
205 
207 
210  void SetIteration(size_t iteration)
211  {
212  iteration_ = iteration;
213  changed_ = true;
214  }
215 
217 
220  void SetTargetIterations(int target)
221  {
222  targetiter_ = target;
223  changed_ = true;
224  }
225 
227 
230  void SetTemperature(double temperature)
231  {
232  temperature_ = temperature;
233  changed_ = true;
234  }
235 
237 
240  void SetEnergy(double energy)
241  {
242  energy_ = energy;
243  changed_ = true;
244  }
245 
247 
250  void SetHMatrix(const Matrix3& hmat)
251  {
252  H_ = hmat;
253  Hinv_ = hmat.inverse();
254  changed_ = true;
255  }
256 
258 
261  void SetVirial(const Matrix3& virial)
262  {
263  virial_ = virial;
264  changed_ = true;
265  }
266 
268 
271  void SetOrigin(const Vector3& origin)
272  {
273  origin_ = origin;
274  changed_ = true;
275  }
276 
278 
281  void SetPeriodicity(const Bool3& isperiodic)
282  {
283  isperiodic_ = isperiodic;
284  changed_ = true;
285  }
286 
288 
291  void SetKb(double kb)
292  {
293  kb_ = kb;
294  changed_ = true;
295  }
296 
298 
301  void SetDielectric(double dielectric)
302  {
303  dielectric_ = dielectric;
304  changed_ = true;
305  }
306 
308 
311  void Setqqrd2e(double qqrd2e)
312  {
313  qqrd2e_ = qqrd2e;
314  changed_ = true;
315  }
316 
318 
321  void SetNumAtoms(unsigned int natoms) { nlocal_ = natoms; }
322 
324 
327  const std::vector<Vector3>& GetPositions() const { return positions_; }
328 
330  std::vector<Vector3>& GetPositions()
331  {
332  changed_ = true;
333  return positions_;
334  }
335 
337 
340  const std::vector<Integer3>& GetImageFlags() const { return images_; }
341 
343  std::vector<Integer3>& GetImageFlags()
344  {
345  changed_ = true;
346  return images_;
347  }
348 
350 
353  const std::vector<Vector3>& GetVelocities() const { return velocities_; }
354 
356  std::vector<Vector3>& GetVelocities()
357  {
358  changed_ = true;
359  return velocities_;
360  }
361 
363 
366  const std::vector<Vector3>& GetForces() const { return forces_; }
367 
369  std::vector<Vector3>& GetForces()
370  {
371  changed_ = true;
372  return forces_;
373  }
374 
376 
382  const std::vector<double>& GetMasses() const { return masses_; }
383 
385  std::vector<double>& GetMasses()
386  {
387  changed_ = true;
388  return masses_;
389  }
390 
392 
396  Vector3 ScaleVector(const Vector3& v) const
397  {
398  return Hinv_*(v-origin_);
399  }
400 
402 
413  Vector3 UnwrapVector(const Vector3& v, const Integer3& image) const
414  {
415  return H_*image.cast<double>()+v;
416  }
417 
419 
422  void ApplyMinimumImage(Vector3* v) const
423  {
424  Vector3 scaled = Hinv_*(*v);
425 
426  for(int i = 0; i < 3; ++i)
427  scaled[i] -= isperiodic_[i]*roundf(scaled[i]);
428 
429  *v = H_*scaled;
430  }
431 
433 
438  {
439  Vector3 scaled = Hinv_*v;
440 
441  for(int i = 0; i < 3; ++i)
442  scaled[i] -= isperiodic_[i]*roundf(scaled[i]);
443 
444  return H_*scaled;
445  }
446 
448 
455  double TotalMass(const Label& indices) const
456  {
457  auto mtot = 0.;
458  for(auto& i : indices)
459  mtot += masses_[i];
460  MPI_Allreduce(MPI_IN_PLACE, &mtot, 1, MPI_DOUBLE, MPI_SUM, comm_);
461  return mtot;
462  }
463 
466 
470  Vector3 CenterOfMass(const Label& indices) const
471  {
472  // Get total mass.
473  auto mtot = TotalMass(indices);
474 
475  return CenterOfMass(indices, mtot);
476  }
477 
480 
487  Vector3 CenterOfMass(const Label& indices, double mtot) const
488  {
489  // Store coorinates and masses in vectors to gather.
490  std::vector<double> pos, mass, gpos, gmass;
491  std::vector<int> pcounts(comm_.size(), 0), mcounts(comm_.size(), 0);
492  std::vector<int> pdispls(comm_.size()+1, 0), mdispls(comm_.size()+1, 0);
493 
494  pcounts[comm_.rank()] = 3*indices.size();
495  mcounts[comm_.rank()] = indices.size();
496 
497  // Reduce counts.
498  MPI_Allreduce(MPI_IN_PLACE, pcounts.data(), pcounts.size(), MPI_INT, MPI_SUM, comm_);
499  MPI_Allreduce(MPI_IN_PLACE, mcounts.data(), mcounts.size(), MPI_INT, MPI_SUM, comm_);
500 
501  // Compute displacements.
502  std::partial_sum(pcounts.begin(), pcounts.end(), pdispls.begin() + 1);
503  std::partial_sum(mcounts.begin(), mcounts.end(), mdispls.begin() + 1);
504 
505  // Fill up mass and position vectors.
506  for(auto& idx : indices)
507  {
508  auto& p = positions_[idx];
509  pos.push_back(p[0]);
510  pos.push_back(p[1]);
511  pos.push_back(p[2]);
512  mass.push_back(masses_[idx]);
513  }
514 
515  // Re-size receiving vectors.
516  gpos.resize(pdispls.back(), 0);
517  gmass.resize(mdispls.back(), 0);
518 
519  // All-gather data.
520  MPI_Allgatherv(pos.data(), pos.size(), MPI_DOUBLE, gpos.data(), pcounts.data(), pdispls.data(), MPI_DOUBLE, comm_);
521  MPI_Allgatherv(mass.data(), mass.size(), MPI_DOUBLE, gmass.data(), mcounts.data(), mdispls.data(), MPI_DOUBLE, comm_);
522 
523  // Loop through atoms and compute mass weighted sum.
524  // We march linearly through list and find nearest image
525  // to each successive particle to properly unwrap object.
526  Vector3 ppos = {gpos[0], gpos[1], gpos[2]}; // Previous unwrapped position.
527  Vector3 cpos = ppos;
528  Vector3 xcm = gmass[0]*cpos;
529 
530  for(size_t i = 1, j = 3; i < gmass.size(); ++i, j += 3)
531  {
532  cpos = {gpos[j], gpos[j+1], gpos[j+2]};
533  cpos = ApplyMinimumImage(cpos - ppos) + ppos;
534  xcm += gmass[i]*cpos;
535  ppos = cpos;
536  }
537 
538  return xcm/mtot;
539  }
540 
542 
549  std::vector<Vector3> UnwrapPositions(const Label& indices) const
550  {
551  // Store coordinates in vectors to gather.
552  std::vector<Vector3> upos;
553  std::vector<double> pos, gpos;
554  std::vector<int> pcounts(comm_.size(), 0);
555  std::vector<int> pdispls(comm_.size()+1, 0);
556 
557  pcounts[comm_.rank()] = 3*indices.size();
558 
559  // Reduce counts.
560  MPI_Allreduce(MPI_IN_PLACE, pcounts.data(), pcounts.size(), MPI_INT, MPI_SUM, comm_);
561 
562  // Compute displacements.
563  std::partial_sum(pcounts.begin(), pcounts.end(), pdispls.begin() + 1);
564 
565  // Fill up position vectors.
566  for(auto& idx : indices)
567  {
568  auto& p = positions_[idx];
569  pos.push_back(p[0]);
570  pos.push_back(p[1]);
571  pos.push_back(p[2]);
572  }
573 
574  // Re-size receiving vectors.
575  gpos.resize(pdispls.back(), 0);
576 
577  // All-gather data.
578  MPI_Allgatherv(pos.data(), pos.size(), MPI_DOUBLE, gpos.data(), pcounts.data(), pdispls.data(), MPI_DOUBLE, comm_);
579 
580  // Loop through atoms and unwrap positions.
581  // We march linearly through list and find nearest image
582  // to each successive particle to properly unwrap object.
583  Vector3 ppos = {gpos[0], gpos[1], gpos[2]}; // Previous unwrapped position.
584  Vector3 cpos = ppos;
585  upos.push_back(cpos); // Used as reference position.
586 
587  for(size_t i = 1; i < gpos.size(); ++i)
588  {
589  cpos = {gpos[3*i], gpos[3*i+1], gpos[3*i+2]};
590  cpos = ApplyMinimumImage(cpos - ppos) + ppos;
591  upos.push_back(cpos);
592  ppos = cpos;
593  }
594 
595  return upos;
596  }
597 
599 
602  const Label& GetAtomIDs() const { return atomids_; }
603 
606  {
607  changed_ = true;
608  return atomids_;
609  }
610 
612 
616  int GetLocalIndex(int id) const
617  {
618 
619  auto s = std::find(atomids_.begin(), atomids_.end(), id);
620  if(s == atomids_.end())
621  return -1;
622  else
623  return s - atomids_.begin();
624  }
625 
628 
634  void GetLocalIndices(const Label& ids, Label* indices) const
635  {
636  for(auto& id : ids)
637  {
638  auto idx = GetLocalIndex(id);
639  if(idx != -1)
640  indices->push_back(idx);
641  }
642  }
643 
645 
648  const std::vector<double>& GetCharges() const { return charges_; }
649 
651  std::vector<double>& GetCharges()
652  {
653  changed_ = true;
654  return charges_;
655  }
656 
658 
661  const Label& GetAtomTypes() const { return types_; }
662 
665  {
666  changed_ = true;
667  return types_;
668  }
669 
671 
674  const std::vector<std::vector<double>>& GetSigmas() const { return sigma_; }
675 
677  std::vector<std::vector<double>>& GetSigmas()
678  {
679  changed_ = true;
680  return sigma_;
681  }
682 
684 
687  const std::string& GetSnapshotID() const { return ID_; }
688 
690  std::string& GetSnapshotID()
691  {
692  changed_ = true;
693  return ID_;
694  }
695 
697 
700  bool HasChanged() const { return changed_; }
701 
703 
706  void Changed(bool state) { changed_ = state; }
707 
709 
712  std::vector<double> SerializePositions()
713  {
714 
715  std::vector<int> pcounts(comm_.size(), 0);
716  std::vector<int> pdispls(comm_.size()+1, 0);
717 
718  pcounts[comm_.rank()] = 3*nlocal_;
719 
720  // Reduce counts.
721  MPI_Allreduce(MPI_IN_PLACE, pcounts.data(), pcounts.size(), MPI_INT, MPI_SUM, comm_);
722 
723  // Compute displacements.
724  std::partial_sum(pcounts.begin(), pcounts.end(), pdispls.begin() + 1);
725 
726  // Re-size receiving vectors.
727  std::vector<double> positions;
728  positions.resize(pdispls.back(), 0);
729 
730  std::vector<double> ptemp;
731 
732  for(auto& p : positions_)
733  {
734  ptemp.push_back(p[0]);
735  ptemp.push_back(p[1]);
736  ptemp.push_back(p[2]);
737  }
738 
739  // All-gather data.
740  MPI_Allgatherv(ptemp.data(), ptemp.size(), MPI_DOUBLE, positions.data(), pcounts.data(), pdispls.data(), MPI_DOUBLE, comm_);
741  return positions;
742  }
743 
745 
748  std::vector<double> SerializeVelocities()
749  {
750  std::vector<int> vcounts(comm_.size(), 0);
751  std::vector<int> vdispls(comm_.size()+1, 0);
752 
753  vcounts[comm_.rank()] = 3*nlocal_;
754 
755  // Reduce counts.
756  MPI_Allreduce(MPI_IN_PLACE, vcounts.data(), vcounts.size(), MPI_INT, MPI_SUM, comm_);
757 
758  // Compute displacements.
759  std::partial_sum(vcounts.begin(), vcounts.end(), vdispls.begin() + 1);
760 
761  // Re-size receiving vectors.
762  std::vector<double> velocities;
763  velocities.resize(vdispls.back(), 0);
764 
765  std::vector<double> vtemp;
766 
767  for(auto& v : velocities_)
768  {
769  vtemp.push_back(v[0]);
770  vtemp.push_back(v[1]);
771  vtemp.push_back(v[2]);
772  }
773 
774  // All-gather data.
775  MPI_Allgatherv(vtemp.data(), vtemp.size(), MPI_DOUBLE, velocities.data(), vcounts.data(), vdispls.data(), MPI_DOUBLE, comm_);
776  return velocities;
777  }
778 
780 
783  std::vector<int> SerializeIDs()
784  {
785  std::vector<int> mcounts(comm_.size(), 0);
786  std::vector<int> mdispls(comm_.size()+1, 0);
787 
788  mcounts[comm_.rank()] = nlocal_;
789 
790  // Reduce counts.
791  MPI_Allreduce(MPI_IN_PLACE, mcounts.data(), mcounts.size(), MPI_INT, MPI_SUM, comm_);
792 
793  // Compute displacements.
794  std::partial_sum(mcounts.begin(), mcounts.end(), mdispls.begin() + 1);
795 
796  // Re-size receiving vectors.
797  std::vector<int> IDs;
798  IDs.resize(mdispls.back(), 0);
799 
800  // All-gather data.
801  MPI_Allgatherv(atomids_.data(), atomids_.size(), MPI_INT, IDs.data(), mcounts.data(), mdispls.data(), MPI_INT, comm_);
802  return IDs;
803  }
804 
807  };
808 }
809 
bool changed_
TRUE is Simulation state changed
Definition: Snapshot.h:85
std::vector< Integer3 > & GetImageFlags()
Access the particles image flags.
Definition: Snapshot.h:343
void SetOrigin(const Vector3 &origin)
Change the box origin.
Definition: Snapshot.h:271
unsigned GetWalkerID() const
Get walker ID.
Definition: Snapshot.h:197
std::vector< int > Label
List of integers.
Definition: types.h:48
Eigen::Matrix3d Matrix3
3x3 matrix.
Definition: types.h:42
const std::vector< double > & GetCharges() const
Access the atom charges.
Definition: Snapshot.h:648
Eigen::Vector3i Integer3
Three-dimensional integer vector.
Definition: types.h:39
std::vector< Integer3 > images_
Unwrapped positions.
Definition: Snapshot.h:68
unsigned GetNumAtoms() const
Get number of atoms in this snapshot.
Definition: Snapshot.h:204
const Matrix3 & GetVirial() const
Get box virial.
Definition: Snapshot.h:137
void Changed(bool state)
Set the "changed" flag of the Snapshot.
Definition: Snapshot.h:706
std::vector< Vector3 > forces_
Forces.
Definition: Snapshot.h:70
void SetNumAtoms(unsigned int natoms)
Set number of atoms in this snapshot.
Definition: Snapshot.h:321
const std::vector< Integer3 > & GetImageFlags() const
Access the particles image flags.
Definition: Snapshot.h:340
size_t iteration_
Iteration of Simulation.
Definition: Snapshot.h:77
Matrix3 H_
Parrinello-Rahman box H-matrix.
Definition: Snapshot.h:58
double kb_
Kb from the MD driver.
Definition: Snapshot.h:81
Label & GetAtomIDs()
Access the atom IDs.
Definition: Snapshot.h:605
double TotalMass(const Label &indices) const
Compute the total mass of a group of particles based on index.
Definition: Snapshot.h:455
Class containing a snapshot of the current simulation in time.
Definition: Snapshot.h:47
std::vector< Vector3 > positions_
Positions.
Definition: Snapshot.h:67
const Vector3 & GetOrigin() const
Get origin of the system.
Definition: Snapshot.h:149
std::vector< Vector3 > & GetForces()
Access the per-particle forces.
Definition: Snapshot.h:369
Vector3 ScaleVector(const Vector3 &v) const
Scale a vector into fractional coordinates.
Definition: Snapshot.h:396
Matrix3 Hinv_
Parinello-Rahman box inverse.
Definition: Snapshot.h:59
void SetVirial(const Matrix3 &virial)
Change the box virial.
Definition: Snapshot.h:261
const Bool3 & IsPeriodic() const
Get periodicity of three dimensions.
Definition: Snapshot.h:155
void SetTargetIterations(int target)
Set target iterations.
Definition: Snapshot.h:220
double energy_
Average per-particle energy.
Definition: Snapshot.h:80
std::vector< double > SerializePositions()
Return the serialized positions across all local cores.
Definition: Snapshot.h:712
std::vector< double > charges_
Charges.
Definition: Snapshot.h:72
bool HasChanged() const
Query if Snapshot was modified.
Definition: Snapshot.h:700
double GetDielectric() const
Get dielectric constant.
Definition: Snapshot.h:173
double temperature_
Current temperature.
Definition: Snapshot.h:79
void SetEnergy(double energy)
Change the energy.
Definition: Snapshot.h:240
std::vector< int > SerializeIDs()
Return the serialized IDs across all local cores.
Definition: Snapshot.h:783
std::vector< double > & GetCharges()
Access the atom charges.
Definition: Snapshot.h:651
int GetLocalIndex(int id) const
Gets the local atom index corresponding to an atom ID.
Definition: Snapshot.h:616
Matrix3 virial_
Virial tensor.
Definition: Snapshot.h:61
std::vector< std::vector< double > > & GetSigmas()
Access the atom sigmas.
Definition: Snapshot.h:677
Bool3 isperiodic_
Periodicity of box.
Definition: Snapshot.h:65
const mxx::comm & GetCommunicator() const
Get communicator for walker.
Definition: Snapshot.h:188
const std::vector< std::vector< double > > & GetSigmas() const
Access the atom sigmas.
Definition: Snapshot.h:674
std::vector< std::vector< double > > sigma_
Sigma.
Definition: Snapshot.h:75
Snapshot(const MPI_Comm &comm, unsigned int wid)
Constructor.
Definition: Snapshot.h:96
void GetLocalIndices(const Label &ids, Label *indices) const
Definition: Snapshot.h:634
double Getqqrd2e() const
Get qqrd2e value.
Definition: Snapshot.h:179
void SetPeriodicity(const Bool3 &isperiodic)
Change the periodicity of the system.
Definition: Snapshot.h:281
mxx::comm comm_
Local snapshot (walker) communicator.
Definition: Snapshot.h:51
const Label & GetAtomIDs() const
Access the atom IDs.
Definition: Snapshot.h:602
std::vector< double > masses_
Masses.
Definition: Snapshot.h:71
Label types_
List of Atom types.
Definition: Snapshot.h:74
double GetTemperature() const
Get current temperature.
Definition: Snapshot.h:119
size_t GetTargetIterations() const
Get target iterations.
Definition: Snapshot.h:113
const Label & GetAtomTypes() const
Access the atom types.
Definition: Snapshot.h:661
double GetVolume() const
Get system volume.
Definition: Snapshot.h:161
Vector3 CenterOfMass(const Label &indices) const
Definition: Snapshot.h:470
size_t targetiter_
Iteration target of simulation.
Definition: Snapshot.h:78
~Snapshot()
Destructor.
Definition: Snapshot.h:806
double dielectric_
Dielectric.
Definition: Snapshot.h:83
const std::vector< double > & GetMasses() const
Const access to the particle masses.
Definition: Snapshot.h:382
Label & GetAtomTypes()
Access the atom types.
Definition: Snapshot.h:664
double GetEnergy() const
Get per-particle energy.
Definition: Snapshot.h:125
std::vector< Vector3 > velocities_
Velocities.
Definition: Snapshot.h:69
void SetDielectric(double dielectric)
Set the dielectric constant.
Definition: Snapshot.h:301
void SetKb(double kb)
Change the kb.
Definition: Snapshot.h:291
void SetTemperature(double temperature)
Change the temperature.
Definition: Snapshot.h:230
double roundf(double x)
Quick helper function to round a double.
Definition: Snapshot.h:35
unsigned int wid_
Walker ID.
Definition: Snapshot.h:53
std::vector< Vector3 > & GetVelocities()
Access the particle velocities.
Definition: Snapshot.h:356
Vector3 ApplyMinimumImage(const Vector3 &v) const
Apply minimum image to a vector.
Definition: Snapshot.h:437
Vector3 CenterOfMass(const Label &indices, double mtot) const
Definition: Snapshot.h:487
Vector3 origin_
Box origin.
Definition: Snapshot.h:63
double GetKb() const
Get system Kb.
Definition: Snapshot.h:167
std::string ID_
ID string.
Definition: Snapshot.h:56
Eigen::Vector3d Vector3
Three-dimensional vector.
Definition: types.h:33
const std::string & GetSnapshotID() const
Access the snapshot ID.
Definition: Snapshot.h:687
void ApplyMinimumImage(Vector3 *v) const
Apply minimum image to a vector.
Definition: Snapshot.h:422
void SetHMatrix(const Matrix3 &hmat)
Change the Box H-matrix.
Definition: Snapshot.h:250
std::string & GetSnapshotID()
Access the snapshot ID.
Definition: Snapshot.h:690
std::vector< double > SerializeVelocities()
Return the serialized velocities across all local cores.
Definition: Snapshot.h:748
Label atomids_
List of Atom IDs.
Definition: Snapshot.h:73
const std::vector< Vector3 > & GetVelocities() const
Access the particle velocities.
Definition: Snapshot.h:353
double qqrd2e_
qqrd2e
Definition: Snapshot.h:84
void Setqqrd2e(double qqrd2e)
Set the value for qqrd2e.
Definition: Snapshot.h:311
Vector3 UnwrapVector(const Vector3 &v, const Integer3 &image) const
Unwrap a vector&#39;s real coordinates according to its image replica count.
Definition: Snapshot.h:413
void SetIteration(size_t iteration)
Set the iteration.
Definition: Snapshot.h:210
const std::vector< Vector3 > & GetForces() const
Access the per-particle forces.
Definition: Snapshot.h:366
const std::vector< Vector3 > & GetPositions() const
Access the particle positions.
Definition: Snapshot.h:327
std::vector< Vector3 > UnwrapPositions(const Label &indices) const
Unwrap the particle positions, based on first position.
Definition: Snapshot.h:549
std::vector< double > & GetMasses()
Const access to the particle masses.
Definition: Snapshot.h:385
const Matrix3 & GetHMatrix() const
Get system H-matrix.
Definition: Snapshot.h:131
unsigned int nlocal_
Number of atoms located on this snapshot.
Definition: Snapshot.h:54
std::vector< Vector3 > & GetPositions()
Access the particle positions.
Definition: Snapshot.h:330
size_t GetIteration() const
Get the current iteration.
Definition: Snapshot.h:107
Matrix3 & GetVirial()
Get box virial.
Definition: Snapshot.h:143
Eigen::Matrix< bool, 3, 1 > Bool3
Three-dimensional boolean.
Definition: types.h:36