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

Swarm of Trajectories String Method. More...

#include <Swarm.h>

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

Public Member Functions

 Swarm (const MPI_Comm &world, const MPI_Comm &comm, const std::vector< double > &centers, unsigned int maxiterations, const std::vector< double > cvspring, unsigned int frequency, unsigned int InitialSteps, unsigned int HarvestLength, unsigned int NumberTrajectories, unsigned int SwarmLength)
 Constructor. More...
 
void PostIntegration (Snapshot *snapshot, const class CVManager &cvmanager) override
 Method call post integration. More...
 
 ~Swarm ()
 Destructor.
 
- Public Member Functions inherited from SSAGES::StringMethod
 StringMethod (const MPI_Comm &world, const MPI_Comm &comm, const std::vector< double > &centers, unsigned int maxiterations, const std::vector< double > cvspring, unsigned int frequency)
 Constructor. More...
 
void PreSimulation (Snapshot *snapshot, const class CVManager &cvmanager) override
 Method call prior to simulation initiation. More...
 
void PostSimulation (Snapshot *snapshot, const class CVManager &cvmanager) override
 Method call post simulation. More...
 
void SetTolerance (std::vector< double > tol)
 Set the tolerance for quitting method. More...
 
void SetSendRecvNeighbors ()
 Communicate neighbor lists over MPI.
 
virtual ~StringMethod ()
 Destructor.
 
- Public Member Functions inherited from SSAGES::Method
 Method (unsigned int frequency, const MPI_Comm &world, const MPI_Comm &comm)
 Constructor. More...
 
void SetCVMask (const std::vector< unsigned int > &mask)
 Sets the collective variable mask. More...
 
virtual ~Method ()
 Destructor.
 
- Public Member Functions inherited from SSAGES::EventListener
 EventListener (unsigned int frequency)
 Constructor. More...
 
unsigned int GetFrequency () const
 Get frequency of event listener. More...
 
virtual ~EventListener ()
 Destructor.
 

Private Member Functions

void StringUpdate () override
 Updates the positions of the string.
 
bool CVInitialized (const CVList &cvs)
 Helper function check if CVs are initialized correctly. More...
 

Private Attributes

std::vector< double > cv_drift_
 Drift of CVs for one iteration.
 
unsigned int initialize_steps_
 Total number of MD steps for initialization for one iteration.
 
unsigned int harvest_length_
 Length to run before harvesting a trajectory for unrestrained sampling.
 
unsigned int restrained_steps_
 Total number of restrained MD steps for one iteration.
 
unsigned int number_trajectories_
 Number of trajectories per swarm.
 
unsigned int swarm_length_
 Length of unrestrained trajectories.
 
unsigned int unrestrained_steps_
 Total number of unrestrained MD steps for one iteration.
 
unsigned int index_
 For indexing trajectory vectors.
 
std::vector< Labeltraj_atomids_
 Store atom IDs for starting trajecotires.
 
bool sampling_started
 Flag for determing whether to perform initialization or not.
 
bool snapshot_stored
 Flag for whether a snapshot was stored in the umbrella sampling.
 
int initialized
 Flag for whether a system is initialized at a given iteration.
 
bool original_initialized
 Flag for whether a system was initialized before it checked whether other systems were.
 

Additional Inherited Members

- Static Public Member Functions inherited from SSAGES::StringMethod
static StringMethodBuild (const Json::Value &json, const MPI_Comm &world, const MPI_Comm &comm, const std::string &path)
 Build a derived method from JSON node. More...
 
- Static Public Member Functions inherited from SSAGES::Method
static MethodBuildMethod (const Json::Value &json, const MPI_Comm &world, const MPI_Comm &comm, const std::string &path)
 Build a derived method from JSON node. More...
 
- Static Public Member Functions inherited from SSAGES::EventListener
static unsigned int GetWalkerID (const MPI_Comm &world, const MPI_Comm &comm)
 Get walker ID number of specified communicator. More...
 
static unsigned int GetNumWalkers (const MPI_Comm &world, const MPI_Comm &comm)
 Get total number of walkers in the simulation. More...
 
static bool IsMasterRank (const MPI_Comm &comm)
 Check if current processor is master. More...
 
- Protected Member Functions inherited from SSAGES::StringMethod
double distance (const std::vector< double > &x, const std::vector< double > &y) const
 Helper function for calculating distances. More...
 
void PrintString (const CVList &CV)
 Prints the CV positions to file. More...
 
void GatherNeighbors (std::vector< double > *lcv0, std::vector< double > *ucv0)
 Gather neighbors over MPI. More...
 
void StringReparam (double alpha_star)
 Reparameterize the string. More...
 
void UpdateWorldString (const CVList &cvs)
 Update the world string over MPI. More...
 
bool TolCheck () const
 Check whether tolerance criteria have been met. More...
 
bool CheckEnd (const CVList &CV)
 Check if method reached one of the exit criteria. More...
 
- Protected Attributes inherited from SSAGES::StringMethod
std::vector< double > centers_
 CV starting location values.
 
std::vector< double > newcenters_
 CV starting location values.
 
std::vector< std::vector< double > > worldstring_
 The world's strings centers for each CV. More...
 
int mpiid_
 The node this belongs to.
 
std::vector< double > tol_
 Tolerance criteria for determining when to stop string (default 0 if no tolerance criteria)
 
int numnodes_
 Number of nodes on a string.
 
unsigned int maxiterator_
 Maximum cap on number of string method iterations performed.
 
std::vector< double > cvspring_
 Vector of spring constants.
 
unsigned int iterator_
 The local method iterator.
 
unsigned int iteration_
 The global method iteration.
 
std::ofstream stringout_
 Output stream for string data.
 
int sendneigh_
 Neighbor to send info to.
 
int recneigh_
 Neighbor to gain info from.
 
std::vector< std::vector< double > > prev_positions_
 Store positions for starting trajectories.
 
std::vector< std::vector< double > > prev_velocities_
 Store velocities for starting trajectories.
 
std::vector< std::vector< int > > prev_IDs_
 Store atom IDs for starting trajectories.
 
- Protected Attributes inherited from SSAGES::Method
mxx::comm world_
 Global MPI communicator.
 
mxx::comm comm_
 Local MPI communicator.
 
std::vector< unsigned int > cvmask_
 Mask which identifies which CVs to act on.
 

Detailed Description

Swarm of Trajectories String Method.

Implementation of the swarm of trajectories string method.

Definition at line 32 of file Swarm.h.

Constructor & Destructor Documentation

◆ Swarm()

SSAGES::Swarm::Swarm ( const MPI_Comm &  world,
const MPI_Comm &  comm,
const std::vector< double > &  centers,
unsigned int  maxiterations,
const std::vector< double >  cvspring,
unsigned int  frequency,
unsigned int  InitialSteps,
unsigned int  HarvestLength,
unsigned int  NumberTrajectories,
unsigned int  SwarmLength 
)
inline

Constructor.

Parameters
worldMPI global communicator.
commMPI local communicator.
centersList of centers.
maxiterationsMaximum number of iterations.
cvspringSpring constants for CVs.
frequencyAplly method with this frequency.
InitialStepsNumber of initial steps.
HarvestLengthLength of trajectory before weighing.
NumberTrajectoriesNumber of trajectories.
SwarmLengthLengt of the swarms.

Constructs an instance of the swarm of trajectories method.

Definition at line 103 of file Swarm.h.

112  :
113  StringMethod(world, comm, centers, maxiterations, cvspring, frequency),
114  cv_drift_(),
115  initialize_steps_(InitialSteps),
116  harvest_length_(HarvestLength),
117  number_trajectories_(NumberTrajectories),
118  swarm_length_(SwarmLength)
119  {
120  cv_drift_.resize(centers_.size(), 0);
124  //Additional initializing
125 
126  index_ = 0;
129  sampling_started = false;
130  snapshot_stored = false;
131 
132  iterator_ = 0; //Override default StringMethod.h initializing
133  }
StringMethod(const MPI_Comm &world, const MPI_Comm &comm, const std::vector< double > &centers, unsigned int maxiterations, const std::vector< double > cvspring, unsigned int frequency)
Constructor.
Definition: StringMethod.h:176
std::vector< std::vector< double > > prev_velocities_
Store velocities for starting trajectories.
Definition: StringMethod.h:88
std::vector< double > centers_
CV starting location values.
Definition: StringMethod.h:43
std::vector< std::vector< double > > prev_positions_
Store positions for starting trajectories.
Definition: StringMethod.h:85
std::vector< std::vector< int > > prev_IDs_
Store atom IDs for starting trajectories.
Definition: StringMethod.h:91
unsigned int iterator_
The local method iterator.
Definition: StringMethod.h:70
unsigned int harvest_length_
Length to run before harvesting a trajectory for unrestrained sampling.
Definition: Swarm.h:43
unsigned int index_
For indexing trajectory vectors.
Definition: Swarm.h:58
unsigned int swarm_length_
Length of unrestrained trajectories.
Definition: Swarm.h:52
bool sampling_started
Flag for determing whether to perform initialization or not.
Definition: Swarm.h:75
unsigned int restrained_steps_
Total number of restrained MD steps for one iteration.
Definition: Swarm.h:46
unsigned int unrestrained_steps_
Total number of unrestrained MD steps for one iteration.
Definition: Swarm.h:55
bool snapshot_stored
Flag for whether a snapshot was stored in the umbrella sampling.
Definition: Swarm.h:78
unsigned int number_trajectories_
Number of trajectories per swarm.
Definition: Swarm.h:49
std::vector< double > cv_drift_
Drift of CVs for one iteration.
Definition: Swarm.h:37
unsigned int initialize_steps_
Total number of MD steps for initialization for one iteration.
Definition: Swarm.h:40

References SSAGES::StringMethod::centers_, cv_drift_, harvest_length_, index_, SSAGES::StringMethod::iterator_, number_trajectories_, SSAGES::StringMethod::prev_IDs_, SSAGES::StringMethod::prev_positions_, SSAGES::StringMethod::prev_velocities_, restrained_steps_, sampling_started, snapshot_stored, swarm_length_, and unrestrained_steps_.

Member Function Documentation

◆ CVInitialized()

bool SSAGES::Swarm::CVInitialized ( const CVList cvs)
private

Helper function check if CVs are initialized correctly.

Parameters
cvsList of CVs to check
Returns
Boolean if CVs are initialized correctly

Definition at line 35 of file Swarm.cpp.

36  {
37  double threshold = 0.2;
38  const double eps = 0.0000000001;
39  double diff;
40 
41  //On the first iteration, check that the CVs are within (threshold*100)% of the center value they're associated to
42  for(size_t i = 0; i < cvs.size(); i++)
43  {
44  if(std::abs(centers_[i]) <= eps)
45  {//e.g. if centers_[i] = 0.0
46  if(std::abs(cvs[i]->GetValue()) <= threshold)
47  {
48  diff = 0; //If current value is approximately zero then go ahead
49  }
50  else
51  {
52  diff = 2*threshold; //Will always be greater than threshold
53  }
54  //diff = std::abs((cvs[i]->GetValue() - (centers_[i]+eps))) / ((cvs[i]->GetValue() + (eps + centers_[i]))/2.0);
55  //std::cout << "MPIID = " << mpiid_ << " and center[" << i << "] = 0.0 and diff = " << diff << std::endl;
56  }
57  else
58  {
59  diff = std::abs((cvs[i]->GetValue() - (centers_[i]))) / ((cvs[i]->GetValue() + (centers_[i]))/2.0);
60  //std::cout << "MPIID = " << mpiid_ << " and center[" << i << "] != 0.0 and diff = " << diff << std::endl;
61  }
62  if(std::abs(diff) >= threshold)
63  {
64  return false; //proceed to initialize again
65  }
66  }
67  return true; //OK to move on to regular sampling
68  }

References SSAGES::StringMethod::centers_.

Referenced by PostIntegration().

Here is the caller graph for this function:

◆ PostIntegration()

void SSAGES::Swarm::PostIntegration ( Snapshot snapshot,
const class CVManager cvmanager 
)
overridevirtual

Method call post integration.

Post-integration hook

Parameters
snapshotPointer to the simulation snapshot.
cvmanagerCollective variable manager.

This function will be called after each integration step.

Implements SSAGES::StringMethod.

Definition at line 70 of file Swarm.cpp.

71  {
72  auto cvs = cvmanager.GetCVs(cvmask_);
73  auto& forces = snapshot->GetForces();
74  auto& positions = snapshot->GetPositions();
75  auto& velocities = snapshot->GetVelocities();
76 
77  //First check if a previous snapshot was stored; the system is definitely initialized if so
78  if(snapshot_stored)
79  {
80  initialized = true;
81  }
82  else
83  {
84  initialized = CVInitialized(cvs); //If no snapshot stored, evalute whether the system finish initializing
85  }
86 
87  original_initialized = initialized; //Will remember whether this system was initialized before the MPI call overrides it
88 
89  //If system finished initializing and hadn't entered the sampling blocks, check if a snapshot was stored. If not, store one and don't store any future snapshots.
91  {
92  if(!snapshot_stored)
93  {
94  prev_positions_[index_] = snapshot->SerializePositions();
95  prev_velocities_[index_] = snapshot->SerializeVelocities();
96  prev_IDs_[index_] = snapshot->SerializeIDs();
97  snapshot_stored = true;
98  }
99  }
100  MPI_Allreduce(MPI_IN_PLACE, &initialized, 1, MPI_INT, MPI_LAND, world_); //Check if all systems are initialized; initialize is only true if everything node finished initializing
101  //If any system wasn't initialized and the sampling block were not entered, begin restraining each node
103  {
104  //Do restrained sampling; NO trajectory harvesting
106  {
107  //If the system was not initialized originally, restrain it
108  for(size_t i = 0; i < cvs.size(); i++)
109  {
110  //Get current CV and gradient
111  auto& cv = cvs[i];
112  auto& grad = cv->GetGradient();
113 
114  //Compute dV/dCV
115  auto D = cvspring_[i]*(cv->GetDifference(centers_[i]));
116 
117  //Update forces
118  for(size_t j = 0; j < forces.size(); j++)
119  {
120  for(int k = 0; k < forces[j].size(); k++)
121  {
122  forces[j][k] -= (double)D*grad[j][k];
123  }
124  }
125  }
126  }
127  else
128  {
129  //If the system was already initialized, simply reset the positions and velocities; this keeps them in their initialized states
130  index_ = 0;
131  for(auto& force: forces)
132  force.setZero();
133 
134  for(size_t i = 0; i < prev_IDs_[index_].size(); i++)
135  {
136  auto localindex = snapshot->GetLocalIndex(prev_IDs_[index_][i]);
137  if(localindex!= -1)
138  {
139  positions[localindex][0] = prev_positions_[index_][i*3];
140  positions[localindex][1] = prev_positions_[index_][i*3 + 1];
141  positions[localindex][2] = prev_positions_[index_][i*3 + 2];
142 
143  velocities[localindex][0] = prev_velocities_[index_][i*3];
144  velocities[localindex][1] = prev_velocities_[index_][i*3 + 1];
145  velocities[localindex][2] = prev_velocities_[index_][i*3 + 2];
146 
147  }
148  }
149  }
150  }
151  else
152  {
153  //Enter this block is every node was done initializing, or the sampling had already begun
154  if(!sampling_started)
155  {
156  sampling_started = true; //If the node got here right after initializing, set the flag that sampling has begun (preventing unneeded initializing in the future)
157  }
159  {
160  //Do restrained sampling, and do harvest trajectories
161  for(size_t i = 0; i < cvs.size(); i++)
162  {
163  if(iterator_ == 0)
164  {
165  index_ = 0; //Reset index when starting
166  }
167  //Get current CV and gradient
168  auto& cv = cvs[i];
169  auto& grad = cv->GetGradient();
170 
171  //Compute dV/dCV
172  auto D = cvspring_[i]*(cv->GetDifference(centers_[i]));
173 
174  //Update forces
175  for(size_t j = 0; j < forces.size(); j++)
176  {
177  for(int k = 0; k < forces[j].size(); k++)
178  {
179  forces[j][k] -= (double)D*grad[j][k];
180  }
181  }
182  }
184  {
185  //Harvest a trajectory every harvest_length_ steps
186  if(iterator_ % harvest_length_ == 0)
187  {
188  prev_positions_[index_] = snapshot->SerializePositions();
189  prev_velocities_[index_] = snapshot->SerializeVelocities();
190  prev_IDs_[index_] = snapshot->SerializeIDs();
191  index_++;
192  }
193  }
195  {
196  //Reset positions and forces before first call to unrestrained sampling
197  index_ = 0;
198  for(auto& force: forces)
199  force.setZero();
200 
201  for(size_t i = 0; i < prev_IDs_[index_].size(); i++)
202  {
203  auto localindex = snapshot->GetLocalIndex(prev_IDs_[index_][i]);
204  if(localindex!= -1)
205  {
206  positions[localindex][0] = prev_positions_[index_][i*3];
207  positions[localindex][1] = prev_positions_[index_][i*3 + 1];
208  positions[localindex][2] = prev_positions_[index_][i*3 + 2];
209 
210  velocities[localindex][0] = prev_velocities_[index_][i*3];
211  velocities[localindex][1] = prev_velocities_[index_][i*3 + 1];
212  velocities[localindex][2] = prev_velocities_[index_][i*3 + 2];
213 
214  }
215  }
216  }
217  iterator_++;
218  }
220  {
221  //Launch unrestrained trajectories
223  {
224  //End of trajectory, harvest drift
225  for(size_t i = 0; i < cv_drift_.size(); i++)
226  {
227  cv_drift_[i] = (cv_drift_[i]*index_ + cvs[i]->GetMinimumImage(centers_[i]) - centers_[i]) / (index_+1); //Calculate running average of drifts
228  }
229  //Set up for next trajectory
230  index_++;
231 
233  {
234  for(auto& force: forces)
235  force.setZero();
236 
237  for(size_t i = 0; i < prev_IDs_[index_].size(); i++)
238  {
239  auto localindex = snapshot->GetLocalIndex(prev_IDs_[index_][i]);
240  if(localindex!= -1)
241  {
242  positions[localindex][0] = prev_positions_[index_][i*3];
243  positions[localindex][1] = prev_positions_[index_][i*3 + 1];
244  positions[localindex][2] = prev_positions_[index_][i*3 + 2];
245 
246  velocities[localindex][0] = prev_velocities_[index_][i*3];
247  velocities[localindex][1] = prev_velocities_[index_][i*3 + 1];
248  velocities[localindex][2] = prev_velocities_[index_][i*3 + 2];
249 
250  }
251  }
252  }
253  }
254  iterator_++;
255  }
256  else
257  {
258  //Evolve CVs, reparametrize, and reset vectors
259  iteration_++;
260 
261  MPI_Barrier(world_);
262  StringUpdate();
263  CheckEnd(cvs);
264  MPI_Barrier(world_);
265  UpdateWorldString(cvs);
266  PrintString(cvs);
267 
268  iterator_ = 0;
269  index_ = 0;
270  snapshot_stored = false;
271 
272  for(size_t i = 0; i < cv_drift_.size(); i++)
273  {
274  cv_drift_[i] = 0;
275  }
276  }
277  }
278  }
std::vector< unsigned int > cvmask_
Mask which identifies which CVs to act on.
Definition: Method.h:50
mxx::comm world_
Global MPI communicator.
Definition: Method.h:46
std::vector< double > cvspring_
Vector of spring constants.
Definition: StringMethod.h:67
void UpdateWorldString(const CVList &cvs)
Update the world string over MPI.
bool CheckEnd(const CVList &CV)
Check if method reached one of the exit criteria.
unsigned int iteration_
The global method iteration.
Definition: StringMethod.h:73
void PrintString(const CVList &CV)
Prints the CV positions to file.
void StringUpdate() override
Updates the positions of the string.
Definition: Swarm.cpp:280
int initialized
Flag for whether a system is initialized at a given iteration.
Definition: Swarm.h:81
bool original_initialized
Flag for whether a system was initialized before it checked whether other systems were.
Definition: Swarm.h:84
bool CVInitialized(const CVList &cvs)
Helper function check if CVs are initialized correctly.
Definition: Swarm.cpp:35

References SSAGES::StringMethod::centers_, SSAGES::StringMethod::CheckEnd(), cv_drift_, CVInitialized(), SSAGES::Method::cvmask_, SSAGES::StringMethod::cvspring_, SSAGES::CVManager::GetCVs(), SSAGES::Snapshot::GetForces(), SSAGES::Snapshot::GetLocalIndex(), SSAGES::Snapshot::GetPositions(), SSAGES::Snapshot::GetVelocities(), harvest_length_, index_, initialize_steps_, initialized, SSAGES::StringMethod::iteration_, SSAGES::StringMethod::iterator_, number_trajectories_, original_initialized, SSAGES::StringMethod::prev_IDs_, SSAGES::StringMethod::prev_positions_, SSAGES::StringMethod::prev_velocities_, SSAGES::StringMethod::PrintString(), restrained_steps_, sampling_started, SSAGES::Snapshot::SerializeIDs(), SSAGES::Snapshot::SerializePositions(), SSAGES::Snapshot::SerializeVelocities(), snapshot_stored, StringUpdate(), swarm_length_, unrestrained_steps_, SSAGES::StringMethod::UpdateWorldString(), and SSAGES::Method::world_.

Here is the call graph for this function:

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