IMP logo
IMP Reference Guide  2.20.2
The Integrative Modeling Platform
MolecularDynamicsWithWte.h
Go to the documentation of this file.
1 /**
2  * \file IMP/spb/MolecularDynamicsWithWte.h
3  * \brief Simple molecular dynamics optimizer.
4  *
5  * Copyright 2007-2022 IMP Inventors. All rights reserved.
6  *
7  */
8 
9 #ifndef IMPSPB_MOLECULAR_DYNAMICS_WITH_WTE_H
10 #define IMPSPB_MOLECULAR_DYNAMICS_WITH_WTE_H
11 
12 #include <IMP/Optimizer.h>
13 #include <IMP/Particle.h>
14 #include <boost/scoped_array.hpp>
15 #include "IMP/atom/Simulator.h"
16 #include "IMP/atom/atom_macros.h"
17 #include <IMP/spb/spb_config.h>
18 
19 IMPSPB_BEGIN_NAMESPACE
20 
21 //! Simple molecular dynamics optimizer.
22 /** The particles to be optimized must have optimizable x,y,z attributes
23  and a non-optimizable mass attribute; this optimizer assumes the score
24  to be energy in kcal/mol, the xyz coordinates to be in angstroms, and
25  the mass to be in AMU (g/mol).
26 
27  Particles without optimized x,y,z and nonoptimized mass are skipped.
28  \see VelocityScalingOptimizerState
29  \see LangevinThermostatOptimizerState
30  \see BerendsenThermostatOptimizerState
31  \see RemoveRigidMotionOptimizerState
32  */
33 class IMPSPBEXPORT MolecularDynamicsWithWte : public atom::Simulator {
34  private:
35  double min_, max_, sigma_, gamma_, dx_, w0_, currentscore_;
36  boost::scoped_array<double> bias_;
37  int nbin_;
38  void update_bias(double score);
39  double get_derivative(double score) const;
40  double deriv_to_acceleration_;
41 
42  public:
43  /** Score based on the provided model */
44  MolecularDynamicsWithWte(Model *m, double emin, double emax, double sigma,
45  double gamma, double w0);
46 
47  double get_bias(double score) const;
48 
49  Floats get_bias_buffer() const {
50  Floats buffer(bias_.get(), bias_.get() + 2 * nbin_);
51  return buffer;
52  }
53 
54  int get_nbin() const { return nbin_; }
55 
56  void set_w0(double w0) { w0_ = w0; }
57 
58  void set_bias(const Floats &bias) {
59  IMP_USAGE_CHECK(static_cast<int>(bias.size()) == 2 * nbin_, "Don't match");
60  // Getting over a warning about comparing unsigned to signed int.
61  // It won't work if x >= INT_MIN (max integer size)
62  std::copy(bias.begin(), bias.end(), bias_.get());
63  }
64 
65  //! \return the current kinetic energy of the system, in kcal/mol
66  virtual Float get_kinetic_energy() const;
67 
68  //! \return the current kinetic temperature of the system
69  /** \param[in] ekinetic kinetic energy, e.g. from get_kinetic_energy()
70  */
71  Float get_kinetic_temperature(Float ekinetic) const;
72 
73  //! Set maximum velocity in A/fs
74  /** At each dynamics time step, the absolute value of each velocity
75  component is capped at this value. This prevents spurious strong forces
76  (occasionally encountered with frustrated conformations) from causing
77  large oscillations in the system.
78  By default, velocities are not capped.
79 
80  \note The actual velocities that are capped are the half-step velocities
81  in the velocity Verlet algorithm.
82  */
83  void set_velocity_cap(Float velocity_cap) { velocity_cap_ = velocity_cap; }
84 
85  //! Assign velocities representative of the given temperature
86  virtual void assign_velocities(Float temperature);
87 
88  //! Rescale velocities
89  void rescale_velocities(Float rescale);
90 
91  // IMP_SIMULATOR(MolecularDynamicsWithWte);
92  virtual void setup(const ParticleIndexes &ps) override;
93  virtual double do_step(const ParticleIndexes &sc, double dt) override;
94  virtual bool get_is_simulation_particle(ParticleIndex p) const override;
95 
96  protected:
97  void initialize();
98 
99  virtual void setup_degrees_of_freedom(const ParticleIndexes &ps);
100 
101  //! First part of velocity Verlet (update coordinates and half-step velocity)
102  virtual void propagate_coordinates(const ParticleIndexes &ps,
103  double step_size);
104 
105  //! Second part of velocity Verlet (update velocity)
106  virtual void propagate_velocities(const ParticleIndexes &ps,
107  double step_size);
108 
109  //! Cap a velocity component to the maximum value.
110  inline void cap_velocity_component(Float &vel) {
111  if (vel >= 0.0) {
112  vel = std::min(vel, velocity_cap_);
113  } else {
114  vel = std::max(vel, -velocity_cap_);
115  }
116  }
117 
118  //! Keys of the xyz velocities
119  FloatKey vs_[3];
120 
121  //! Number of degrees of freedom in the system
123 
124  //! Maximum absolute value of a single velocity component
126 };
127 
128 IMPSPB_END_NAMESPACE
129 
130 #endif /* IMPSPB_MOLECULAR_DYNAMICS_WITH_WTE_H */
int degrees_of_freedom_
Number of degrees of freedom in the system.
The base class for simulators.
Definition: Simulator.h:34
Simple molecular dynamics optimizer.
virtual bool get_is_simulation_particle(ParticleIndex p) const =0
Return true if the passed particle is appropriate for the simulation.
void set_velocity_cap(Float velocity_cap)
Set maximum velocity in A/fs.
Base class for all optimizers.
virtual double do_step(const ParticleIndexes &sc, double dt)=0
Perform a single time step.
Class for storing model, its restraints, constraints, and particles.
Definition: Model.h:86
virtual void setup(const ParticleIndexes &)
Definition: Simulator.h:178
void cap_velocity_component(Float &vel)
Cap a velocity component to the maximum value.
Base class for "simulators", such as molecular dynamics.
Classes to handle individual model particles. (Note that implementation of inline functions is in int...
Macros for maintaining molecular hierarchies.
double Float
Basic floating-point value (could be float, double...)
Definition: types.h:19
#define IMP_USAGE_CHECK(expr, message)
A runtime test for incorrect usage of a class or method.
Definition: check_macros.h:168
Float velocity_cap_
Maximum absolute value of a single velocity component.