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

Artificial Neural Network Method. More...

#include <ANN.h>

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

Public Member Functions

 ANN (const MPI_Comm &world, const MPI_Comm &comm, const Eigen::VectorXi &topol, Grid< Eigen::VectorXd > *fgrid, Grid< unsigned int > *hgrid, Grid< double > *ugrid, const std::vector< double > &lowerb, const std::vector< double > &upperb, const std::vector< double > &lowerk, const std::vector< double > &upperk, double temperature, double weight, unsigned int nsweep)
 Constructor. More...
 
void PreSimulation (Snapshot *snapshot, const class CVManager &cvmanager) override
 Pre-simulation hook. More...
 
void PostIntegration (Snapshot *snapshot, const class CVManager &cvmanager) override
 Post-integration hook. More...
 
void PostSimulation (Snapshot *snapshot, const class CVManager &cvmanager) override
 Post-simulation hook. More...
 
void SetPrevWeight (double h)
 Set previous history weight.
 
void SetOutput (const std::string &outfile)
 Set name of output file.
 
void SetOutputOverwrite (bool overwrite)
 Set overwrite flag on output file.
 
void SetConvergeIters (unsigned int citers)
 Set number of iterations after which we turn on full weight.
 
void SetMaxIters (unsigned int iters)
 Set maximum number of training iterations per sweep.
 
void SetMinLoss (double loss)
 Set minimum loss function value (should be zero for production).
 
- 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.
 
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.
 

Static Public Member Functions

static ANNBuild (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...
 

Private Member Functions

void TrainNetwork ()
 Trains the neural network.
 
void WriteBias ()
 Writes out the bias to file.
 

Private Attributes

Eigen::VectorXi topol_
 Neural network topology.
 
unsigned int citers_
 Number of iterations after which we turn on full weight.
 
nnet::neural_net net_
 Neural network.
 
Grid< Eigen::VectorXd > * fgrid_
 Force grid.
 
Grid< unsigned int > * hgrid_
 Histogram grid.
 
Grid< double > * ugrid_
 Unbiased histogram grid.
 
std::string outfile_
 Output filename.
 
bool overwrite_
 Overwrite outputs?
 
unsigned int sweep_
 
unsigned int nsweep_
 
double pweight_
 
double weight_
 
double temp_
 
double kbt_
 
Eigen::MatrixXd hist_
 
Eigen::MatrixXd bias_
 
std::vector< double > lowerb_
 
std::vector< double > upperb_
 
std::vector< double > lowerk_
 
std::vector< double > upperk_
 

Additional Inherited Members

- 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

Artificial Neural Network Method.

Implementation of the Artificial Neural Network Method based on [3]

Definition at line 35 of file ANN.h.

Constructor & Destructor Documentation

◆ ANN()

SSAGES::ANN::ANN ( const MPI_Comm &  world,
const MPI_Comm &  comm,
const Eigen::VectorXi &  topol,
Grid< Eigen::VectorXd > *  fgrid,
Grid< unsigned int > *  hgrid,
Grid< double > *  ugrid,
const std::vector< double > &  lowerb,
const std::vector< double > &  upperb,
const std::vector< double > &  lowerk,
const std::vector< double > &  upperk,
double  temperature,
double  weight,
unsigned int  nsweep 
)

Constructor.

Parameters
worldMPI global communicator.
commMPI local communicator.
topolTopology of network.
fgridGrid containing biasing forces.
hgridGrid containing histogram.
ugridGrid containing unbiased histogram.
lowerbLower bounds for CVs.
upperbUpper bounds for CVs.
lowerkLower bound restraints for CVs.
upperkUpper bound restraints for CVs.
temperatureTemperature of the simulation.
weightRelative weight of the statistics in sweep.
nsweepNumber of iterations in the sweep.

Constructs an instance of Artificial Neural Network method.

Definition at line 34 of file ANN.cpp.

References SSAGES::Grid< T >::begin(), SSAGES::Grid< T >::end(), SSAGES::GridBase< T >::GetDimension(), hgrid_, hist_, net_, and SSAGES::GridBase< T >::size().

Referenced by Build().

46  :
47  Method(1, world, comm), topol_(topol), sweep_(0), nsweep_(nsweep),
48  citers_(0), net_(topol), pweight_(1.), weight_(weight), temp_(temperature),
49  kbt_(0), fgrid_(fgrid), hgrid_(hgrid), ugrid_(ugrid), hist_(), bias_(),
50  lowerb_(lowerb), upperb_(upperb), lowerk_(lowerk), upperk_(upperk),
51  outfile_("ann.out"), overwrite_(true)
52  {
53  // Create histogram grid matrix.
54  hist_.resize(hgrid_->size(), hgrid_->GetDimension());
55 
56  // Fill it up.
57  size_t i = 0;
58  for(auto it = hgrid_->begin(); it != hgrid_->end(); ++it)
59  {
60  auto coord = it.coordinates();
61  for(size_t j = 0; j < coord.size(); ++j)
62  hist_(i, j) = coord[j];
63  ++i;
64  }
65 
66  // Initialize FES vector.
67  bias_.resize(hgrid_->size(), 1);
68  net_.forward_pass(hist_);
69  bias_.array() = net_.get_activation().col(0).array();
70  }
size_t size() const
Get the size of the internal storage vector.
Definition: GridBase.h:321
std::vector< double > lowerk_
Definition: ANN.h:83
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:190
unsigned int citers_
Number of iterations after which we turn on full weight.
Definition: ANN.h:47
Grid< Eigen::VectorXd > * fgrid_
Force grid.
Definition: ANN.h:63
unsigned int sweep_
Definition: ANN.h:43
iterator begin()
Return iterator at first grid point.
Definition: Grid.h:527
double temp_
Definition: ANN.h:59
Grid< unsigned int > * hgrid_
Histogram grid.
Definition: ANN.h:66
Eigen::VectorXi topol_
Neural network topology.
Definition: ANN.h:39
bool overwrite_
Overwrite outputs?
Definition: ANN.h:90
std::string outfile_
Output filename.
Definition: ANN.h:87
nnet::neural_net net_
Neural network.
Definition: ANN.h:50
std::vector< double > lowerb_
Definition: ANN.h:78
Eigen::MatrixXd hist_
Definition: ANN.h:73
Grid< double > * ugrid_
Unbiased histogram grid.
Definition: ANN.h:69
Method(unsigned int frequency, const MPI_Comm &world, const MPI_Comm &comm)
Constructor.
Definition: Method.h:61
double pweight_
Definition: ANN.h:54
iterator end()
Return iterator after last valid grid point.
Definition: Grid.h:540
Here is the call graph for this function:
Here is the caller graph for this function:

Member Function Documentation

◆ Build()

ANN * SSAGES::ANN::Build ( const Json::Value &  json,
const MPI_Comm &  world,
const MPI_Comm &  comm,
const std::string &  path 
)
static

Build a derived method from JSON node.

Parameters
jsonJSON Value containing all input information.
worldMPI global communicator.
commMPI local communicator.
pathPath for JSON path specification.
Returns
Pointer to the Method built. nullptr if an unknown error occurred.

This function builds a registered method from a JSON node. The difference between this function and "Build" is that this automatically determines the appropriate derived type based on the JSON node information.

Note
Object lifetime is the caller's responsibility.

Definition at line 234 of file ANN.cpp.

References ANN(), SSAGES::Grid< T >::BuildGrid(), Json::Requirement::GetErrors(), Json::Requirement::HasErrors(), Json::ObjectRequirement::Parse(), and Json::ObjectRequirement::Validate().

Referenced by SetMinLoss().

239  {
240  ObjectRequirement validator;
241  Value schema;
242  CharReaderBuilder rbuilder;
243  CharReader* reader = rbuilder.newCharReader();
244 
245  reader->parse(JsonSchema::ANNMethod.c_str(),
246  JsonSchema::ANNMethod.c_str() + JsonSchema::ANNMethod.size(),
247  &schema, NULL);
248  validator.Parse(schema, path);
249 
250  // Validate inputs.
251  validator.Validate(json, path);
252  if(validator.HasErrors())
253  throw BuildException(validator.GetErrors());
254 
255  // Grid.
256  auto* fgrid = Grid<VectorXd>::BuildGrid(json.get("grid", Json::Value()));
257  auto* hgrid = Grid<unsigned int>::BuildGrid(json.get("grid", Json::Value()));
258  auto* ugrid = Grid<double>::BuildGrid(json.get("grid", Json::Value()));
259 
260  // Topology.
261  auto nlayers = json["topology"].size() + 2;
262  VectorXi topol(nlayers);
263  topol[0] = fgrid->GetDimension();
264  topol[nlayers-1] = 1;
265  for(int i = 0; i < static_cast<int>(json["topology"].size()); ++i)
266  topol[i+1] = json["topology"][i].asInt();
267 
268  auto weight = json.get("weight", 1.).asDouble();
269  auto temp = json["temperature"].asDouble();
270  auto nsweep = json["nsweep"].asUInt();
271 
272  // Assume all vectors are the same size.
273  std::vector<double> lowerb, upperb, lowerk, upperk;
274  for(int i = 0; i < static_cast<int>(json["lower_bound_restraints"].size()); ++i)
275  {
276  lowerk.push_back(json["lower_bound_restraints"][i].asDouble());
277  upperk.push_back(json["upper_bound_restraints"][i].asDouble());
278  lowerb.push_back(json["lower_bounds"][i].asDouble());
279  upperb.push_back(json["upper_bounds"][i].asDouble());
280  }
281 
282  auto* m = new ANN(world, comm, topol, fgrid, hgrid, ugrid, lowerb, upperb, lowerk, upperk, temp, weight, nsweep);
283 
284  // Set optional params.
285  m->SetPrevWeight(json.get("prev_weight", 1).asDouble());
286  m->SetOutput(json.get("output_file", "ann.out").asString());
287  m->SetOutputOverwrite( json.get("overwrite_output", true).asBool());
288  m->SetConvergeIters(json.get("converge_iters", 0).asUInt());
289  m->SetMaxIters(json.get("max_iters", 1000).asUInt());
290  m->SetMinLoss(json.get("min_loss", 0).asDouble());
291 
292  return m;
293  }
bool HasErrors()
Check if errors have occured.
Definition: Requirement.h:86
ANN(const MPI_Comm &world, const MPI_Comm &comm, const Eigen::VectorXi &topol, Grid< Eigen::VectorXd > *fgrid, Grid< unsigned int > *hgrid, Grid< double > *ugrid, const std::vector< double > &lowerb, const std::vector< double > &upperb, const std::vector< double > &lowerk, const std::vector< double > &upperk, double temperature, double weight, unsigned int nsweep)
Constructor.
Definition: ANN.cpp:34
virtual void Parse(Value json, const std::string &path) override
Parse JSON value to generate Requirement(s).
std::vector< std::string > GetErrors()
Get list of error messages.
Definition: Requirement.h:92
Requirements on an object.
static Grid< T > * BuildGrid(const Json::Value &json)
Set up the grid.
Definition: Grid.h:127
virtual void Validate(const Value &json, const std::string &path) override
Validate JSON value.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ PostIntegration()

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

Post-integration hook.

Parameters
snapshotCurrent simulation snapshot.
cvmanagerCollective variable manager.

Implements SSAGES::Method.

Definition at line 84 of file ANN.cpp.

References SSAGES::GridBase< T >::at(), citers_, SSAGES::Method::comm_, SSAGES::Method::cvmask_, SSAGES::CVManager::GetCVs(), SSAGES::Snapshot::GetForces(), SSAGES::Snapshot::GetIteration(), SSAGES::GridBase< T >::GetLower(), SSAGES::GridBase< T >::GetUpper(), SSAGES::Snapshot::GetVirial(), hgrid_, lowerb_, lowerk_, net_, pweight_, TrainNetwork(), SSAGES::Method::world_, and WriteBias().

85  {
86  if(snapshot->GetIteration() && snapshot->GetIteration() % nsweep_ == 0)
87  {
88  // Switch to full blast.
89  if(citers_ && snapshot->GetIteration() > citers_)
90  pweight_ = 1.0;
91 
92  TrainNetwork();
93  if(world_.rank() == 0)
94  WriteBias();
95  }
96 
97  // Get CV vals.
98  auto cvs = cvmanager.GetCVs(cvmask_);
99  auto n = cvs.size();
100 
101  // Determine if we are in bounds.
102  RowVectorXd vec(n);
103  std::vector<double> val(n);
104  bool inbounds = true;
105  for(size_t i = 0; i < n; ++i)
106  {
107  val[i] = cvs[i]->GetValue();
108  vec[i] = cvs[i]->GetValue();
109  if(val[i] < hgrid_->GetLower(i) || val[i] > hgrid_->GetUpper(i))
110  inbounds = false;
111  }
112 
113  // If in bounds, bias.
114  VectorXd derivatives = VectorXd::Zero(n);
115  if(inbounds)
116  {
117  // Record histogram hit and get gradient.
118  // Only record hits on master processes since we will
119  // reduce later.
120  if(comm_.rank() == 0)
121  hgrid_->at(val) += 1;
122  //derivatives = (*fgrid_)[val];
123  net_.forward_pass(vec);
124  derivatives = net_.get_gradient(0);
125  }
126  else
127  {
128  if(comm_.rank() == 0)
129  {
130  std::cerr << "ANN (" << snapshot->GetIteration() << "): out of bounds ( ";
131  for(auto& v : val)
132  std::cerr << v << " ";
133  std::cerr << ")" << std::endl;
134  }
135  }
136 
137  // Restraints.
138  for(size_t i = 0; i < n; ++i)
139  {
140  auto cval = cvs[i]->GetValue();
141  if(cval < lowerb_[i])
142  derivatives[i] += lowerk_[i]*cvs[i]->GetDifference(cval - lowerb_[i]);
143  else if(cval > upperb_[i])
144  derivatives[i] += upperk_[i]*cvs[i]->GetDifference(cval - upperb_[i]);
145  }
146 
147  // Apply bias to atoms.
148  auto& forces = snapshot->GetForces();
149  auto& virial = snapshot->GetVirial();
150 
151  for(size_t i = 0; i < cvs.size(); ++i)
152  {
153  auto& grad = cvs[i]->GetGradient();
154  auto& boxgrad = cvs[i]->GetBoxGradient();
155 
156  // Update the forces in snapshot by adding in the force bias from each
157  // CV to each atom based on the gradient of the CV.
158  for (size_t j = 0; j < forces.size(); ++j)
159  forces[j] -= derivatives[i]*grad[j];
160 
161  virial += derivatives[i]*boxgrad;
162  }
163  }
const std::vector< double > GetLower() const
Return the lower edges of the Grid.
Definition: GridBase.h:226
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:541
std::vector< double > lowerk_
Definition: ANN.h:83
mxx::comm comm_
Local MPI communicator.
Definition: Method.h:47
mxx::comm world_
Global MPI communicator.
Definition: Method.h:46
unsigned int citers_
Number of iterations after which we turn on full weight.
Definition: ANN.h:47
void WriteBias()
Writes out the bias to file.
Definition: ANN.cpp:215
Grid< unsigned int > * hgrid_
Histogram grid.
Definition: ANN.h:66
const std::vector< double > GetUpper() const
Return the upper edges of the Grid.
Definition: GridBase.h:257
void TrainNetwork()
Trains the neural network.
Definition: ANN.cpp:169
std::vector< unsigned int > cvmask_
Mask which identifies which CVs to act on.
Definition: Method.h:50
nnet::neural_net net_
Neural network.
Definition: ANN.h:50
std::vector< double > lowerb_
Definition: ANN.h:78
double pweight_
Definition: ANN.h:54
Here is the call graph for this function:

◆ PostSimulation()

void SSAGES::ANN::PostSimulation ( Snapshot snapshot,
const class CVManager cvmanager 
)
overridevirtual

Post-simulation hook.

Parameters
snapshotCurrent simulation snapshot.
cvmanagerCollective variable manager.

Implements SSAGES::Method.

Definition at line 165 of file ANN.cpp.

166  {
167  }

◆ PreSimulation()

void SSAGES::ANN::PreSimulation ( Snapshot snapshot,
const class CVManager cvmanager 
)
overridevirtual

Pre-simulation hook.

Parameters
snapshotCurrent simulation snapshot.
cvmanagerCollective variable manager.

Implements SSAGES::Method.

Definition at line 72 of file ANN.cpp.

References SSAGES::Grid< T >::begin(), SSAGES::Grid< T >::end(), fgrid_, SSAGES::GridBase< T >::GetDimension(), SSAGES::Snapshot::GetKb(), hgrid_, temp_, and ugrid_.

73  {
74  auto ndim = hgrid_->GetDimension();
75  kbt_ = snapshot->GetKb()*temp_;
76 
77  // Zero out forces and histogram.
78  VectorXd vec = VectorXd::Zero(ndim);
79  std::fill(hgrid_->begin(), hgrid_->end(), 0);
80  std::fill(ugrid_->begin(), ugrid_->end(), 1.0);
81  std::fill(fgrid_->begin(), fgrid_->end(), vec);
82  }
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:190
Grid< Eigen::VectorXd > * fgrid_
Force grid.
Definition: ANN.h:63
iterator begin()
Return iterator at first grid point.
Definition: Grid.h:527
double temp_
Definition: ANN.h:59
Grid< unsigned int > * hgrid_
Histogram grid.
Definition: ANN.h:66
Grid< double > * ugrid_
Unbiased histogram grid.
Definition: ANN.h:69
iterator end()
Return iterator after last valid grid point.
Definition: Grid.h:540
Here is the call graph for this function:

Member Data Documentation

◆ hist_

Eigen::MatrixXd SSAGES::ANN::hist_
private

Eigen matrices of grids.

Definition at line 73 of file ANN.h.

Referenced by ANN(), TrainNetwork(), and WriteBias().

◆ lowerb_

std::vector<double> SSAGES::ANN::lowerb_
private

Bounds

Definition at line 78 of file ANN.h.

Referenced by PostIntegration().

◆ lowerk_

std::vector<double> SSAGES::ANN::lowerk_
private

Bound restraints.

Definition at line 83 of file ANN.h.

Referenced by PostIntegration().

◆ pweight_

double SSAGES::ANN::pweight_
private

Previous and current histogram weight.

Definition at line 54 of file ANN.h.

Referenced by PostIntegration(), and TrainNetwork().

◆ sweep_

unsigned int SSAGES::ANN::sweep_
private

Number of iterations per sweep.

Definition at line 43 of file ANN.h.

Referenced by TrainNetwork(), and WriteBias().

◆ temp_

double SSAGES::ANN::temp_
private

System temperature and energy units.

Definition at line 59 of file ANN.h.

Referenced by PreSimulation().


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