Update Colvars to version 2018-08-29 47/4547/1
authorGiacomo Fiorin <giacomo.fiorin@gmail.com>
Thu, 6 Sep 2018 22:44:07 +0000 (18:44 -0400)
committerGiacomo Fiorin <giacomo.fiorin@gmail.com>
Thu, 6 Sep 2018 22:44:07 +0000 (18:44 -0400)
This update includes several minor fixes: relevant commits in
https://github.com/Colvars/colvars are listed below.

6539dd9 2018-08-28 Add cv resetindexgroups command [Giacomo Fiorin]
a9c786e 2018-08-28 Invoke new-style colvarscript commands when available [Giacomo Fiorin]
44a4b34 2018-08-28 Allow multiple instances of indexFile [Giacomo Fiorin]
6f630a4 2018-08-28 Store the module's config string with comments, add command to access it [Giacomo Fiorin]
8272aab 2018-08-17 Raise error when the same keyword is being used repeatedly [Giacomo Fiorin]
10b2840 2018-07-30 More explicit documentation for periodic custom cv [Jérôme Hénin]
788cb54 2018-07-25 Check for matching braces before beginning to parse [Giacomo Fiorin]
ac2ff20 2018-07-24 Initialize correctly OpenMP lock state data structure [Giacomo Fiorin]
acdf322 2018-07-20 NAMD/VMD: accept PDB flags in load_coords [Jérôme Hénin]
f3f7033 2018-07-20 Doc: warn about fittingGroup in PBC [Jérôme Hénin]
61b5f6f 2018-07-07 Forbid atoms with zero mass [Jérôme Hénin]
127375e 2018-06-21 Revert ext lagrangian timestep if necessary [Jérôme Hénin]
fb67424 2018-06-06 Implement modifycvcs subcommand [Giacomo Fiorin]
8a4d86c 2018-06-06 Disable echoing of errno message until sorting out what changes errno [Giacomo Fiorin]
bc09f23 2018-06-05 Prevent renaming CVCs after initialization [Giacomo Fiorin]
33794f3 2018-06-05 Fix handling of previous/current step total forces with SMP [Giacomo Fiorin]

Change-Id: I6d73d3695e9c3f069871f6f39fa06ca6a6e34de4

20 files changed:
colvars/src/colvar.cpp
colvars/src/colvar.h
colvars/src/colvaratoms.cpp
colvars/src/colvarcomp.cpp
colvars/src/colvarcomp.h
colvars/src/colvarcomp_distances.cpp
colvars/src/colvarmodule.cpp
colvars/src/colvarmodule.h
colvars/src/colvarparse.cpp
colvars/src/colvarparse.h
colvars/src/colvarproxy.cpp
colvars/src/colvarproxy.h
colvars/src/colvars_version.h
colvars/src/colvarscript.cpp
colvars/src/colvarscript.h
src/colvarproxy_namd.C
src/colvarproxy_namd.h
src/colvarproxy_namd_version.h
ug/ug_colvars-cv.tex
ug/ug_colvars.tex

index c5b6d90..43deef3 100644 (file)
@@ -677,6 +677,7 @@ template<typename def_class_name> int colvar::init_components_type(std::string c
     if (cvcp != NULL) {
       cvcs.push_back(cvcp);
       cvcp->check_keywords(def_conf, def_config_key);
+      cvcp->config_key = def_config_key;
       if (cvm::get_error()) {
         cvm::error("Error: in setting up component \""+
                    std::string(def_config_key)+"\".\n", INPUT_ERROR);
@@ -1022,13 +1023,13 @@ int colvar::calc()
 
 int colvar::calc_cvcs(int first_cvc, size_t num_cvcs)
 {
-  colvarproxy *proxy = cvm::main()->proxy;
-
-  int error_code = COLVARS_OK;
   if (cvm::debug())
     cvm::log("Calculating colvar \""+this->name+"\", components "+
              cvm::to_str(first_cvc)+" through "+cvm::to_str(first_cvc+num_cvcs)+".\n");
 
+  colvarproxy *proxy = cvm::main()->proxy;
+  int error_code = COLVARS_OK;
+
   error_code |= check_cvc_range(first_cvc, num_cvcs);
   if (error_code != COLVARS_OK) {
     return error_code;
@@ -1059,9 +1060,10 @@ int colvar::collect_cvc_data()
   if (cvm::debug())
     cvm::log("Calculating colvar \""+this->name+"\"'s properties.\n");
 
+  colvarproxy *proxy = cvm::main()->proxy;
   int error_code = COLVARS_OK;
 
-  if (cvm::step_relative() > 0) {
+  if ((cvm::step_relative() > 0) && (!proxy->total_forces_same_step())){
     // Total force depends on Jacobian derivative from previous timestep
     // collect_cvc_total_forces() uses the previous value of jd
     error_code |= collect_cvc_total_forces();
@@ -1069,6 +1071,10 @@ int colvar::collect_cvc_data()
   error_code |= collect_cvc_values();
   error_code |= collect_cvc_gradients();
   error_code |= collect_cvc_Jacobians();
+  if (proxy->total_forces_same_step()){
+    // Use Jacobian derivative from this timestep
+    error_code |= collect_cvc_total_forces();
+  }
   error_code |= calc_colvar_properties();
 
   if (cvm::debug())
@@ -1394,6 +1400,13 @@ int colvar::calc_colvar_properties()
       vr.reset(); // (already 0; added for clarity)
     }
 
+    // Special case of a repeated timestep (eg. multiple NAMD "run" statements)
+    // revert values of the extended coordinate and velocity prior to latest integration
+    if (cvm::step_relative() == prev_timestep) {
+      xr = prev_xr;
+      vr = prev_vr;
+    }
+
     // report the restraint center as "value"
     x_reported = xr;
     v_reported = vr;
@@ -1450,7 +1463,6 @@ cvm::real colvar::update_forces_energy()
   // extended variable if there is one
 
   if (is_enabled(f_cv_extended_Lagrangian)) {
-
     if (cvm::debug()) {
       cvm::log("Updating extended-Lagrangian degree of freedom.\n");
     }
@@ -1502,6 +1514,11 @@ cvm::real colvar::update_forces_energy()
       ft_reported = f_ext;
     }
 
+    // backup in case we need to revert this integration timestep
+    // if the same MD timestep is re-run
+    prev_xr = xr;
+    prev_vr = vr;
+
     // leapfrog: starting from x_i, f_i, v_(i-1/2)
     vr  += (0.5 * dt) * f_ext / ext_mass;
     // Because of leapfrog, kinetic energy at time i is approximate
@@ -1638,18 +1655,26 @@ int colvar::set_cvc_flags(std::vector<bool> const &flags)
 }
 
 
+void colvar::update_active_cvc_square_norm()
+{
+  active_cvc_square_norm = 0.0;
+  for (size_t i = 0; i < cvcs.size(); i++) {
+    if (cvcs[i]->is_enabled()) {
+      active_cvc_square_norm += cvcs[i]->sup_coeff * cvcs[i]->sup_coeff;
+    }
+  }
+}
+
+
 int colvar::update_cvc_flags()
 {
   // Update the enabled/disabled status of cvcs if necessary
   if (cvc_flags.size()) {
     n_active_cvcs = 0;
-    active_cvc_square_norm = 0.;
-
     for (size_t i = 0; i < cvcs.size(); i++) {
       cvcs[i]->set_enabled(f_cvc_active, cvc_flags[i]);
       if (cvcs[i]->is_enabled()) {
         n_active_cvcs++;
-        active_cvc_square_norm += cvcs[i]->sup_coeff * cvcs[i]->sup_coeff;
       }
     }
     if (!n_active_cvcs) {
@@ -1657,12 +1682,49 @@ int colvar::update_cvc_flags()
       return COLVARS_ERROR;
     }
     cvc_flags.clear();
+
+    update_active_cvc_square_norm();
   }
 
   return COLVARS_OK;
 }
 
 
+int colvar::update_cvc_config(std::vector<std::string> const &confs)
+{
+  cvm::log("Updating configuration for colvar \""+name+"\"");
+
+  if (confs.size() != cvcs.size()) {
+    return cvm::error("Error: Wrong number of CVC config strings.  "
+                      "For those CVCs that are not being changed, try passing "
+                      "an empty string.", INPUT_ERROR);
+  }
+
+  int error_code = COLVARS_OK;
+  int num_changes = 0;
+  for (size_t i = 0; i < cvcs.size(); i++) {
+    if (confs[i].size()) {
+      std::string conf(confs[i]);
+      cvm::increase_depth();
+      error_code |= cvcs[i]->colvar::cvc::init(conf);
+      error_code |= cvcs[i]->check_keywords(conf,
+                                            cvcs[i]->config_key.c_str());
+      cvm::decrease_depth();
+      num_changes++;
+    }
+  }
+
+  if (num_changes == 0) {
+    cvm::log("Warning: no changes were applied through modifycvcs; "
+             "please check that its argument is a list of strings.\n");
+  }
+
+  update_active_cvc_square_norm();
+
+  return error_code;
+}
+
+
 // ******************** METRIC FUNCTIONS ********************
 // Use the metrics defined by \link cvc \endlink objects
 
index 271231f..989d551 100644 (file)
@@ -171,8 +171,12 @@ protected:
   // Options for extended_lagrangian
   /// Restraint center
   colvarvalue xr;
+  /// Previous value of the restraint center;
+  colvarvalue prev_xr;
   /// Velocity of the restraint center
   colvarvalue vr;
+  /// Previous velocity of the restraint center
+  colvarvalue prev_vr;
   /// Mass of the restraint center
   cvm::real ext_mass;
   /// Restraint force constant
@@ -352,6 +356,9 @@ public:
   /// \brief Updates the flags in the CVC objects, and their number
   int update_cvc_flags();
 
+  /// \brief Modify the configuration of CVCs (currently, only base class data)
+  int update_cvc_config(std::vector<std::string> const &confs);
+
 protected:
   /// \brief Number of CVC objects with an active flag
   size_t n_active_cvcs;
@@ -359,10 +366,17 @@ protected:
   /// Sum of square coefficients for active cvcs
   cvm::real active_cvc_square_norm;
 
+  /// Update the sum of square coefficients for active cvcs
+  void update_active_cvc_square_norm();
+
   /// \brief Absolute timestep number when this colvar was last updated
   int prev_timestep;
 
 public:
+
+  /// \brief Return the number of CVC objects defined
+  inline size_t num_cvcs() const { return cvcs.size(); }
+
   /// \brief Return the number of CVC objects with an active flag (as set by update_cvc_flags)
   inline size_t num_active_cvcs() const { return n_active_cvcs; }
 
@@ -371,21 +385,21 @@ public:
   ///
   /// Handles correctly symmetries and periodic boundary conditions
   cvm::real dist2(colvarvalue const &x1,
-                   colvarvalue const &x2) const;
+                  colvarvalue const &x2) const;
 
   /// \brief Use the internal metrics (as from \link cvc
   /// \endlink objects) to calculate square distances and gradients
   ///
   /// Handles correctly symmetries and periodic boundary conditions
   colvarvalue dist2_lgrad(colvarvalue const &x1,
-                           colvarvalue const &x2) const;
+                          colvarvalue const &x2) const;
 
   /// \brief Use the internal metrics (as from \link cvc
   /// \endlink objects) to calculate square distances and gradients
   ///
   /// Handles correctly symmetries and periodic boundary conditions
   colvarvalue dist2_rgrad(colvarvalue const &x1,
-                           colvarvalue const &x2) const;
+                          colvarvalue const &x2) const;
 
   /// \brief Use the internal metrics (as from \link cvc
   /// \endlink objects) to wrap a value into a standard interval
index c7b6d3d..3315007 100644 (file)
@@ -397,7 +397,7 @@ int cvm::atom_group::parse(std::string const &group_conf)
       }
 
       // NOTE: calls to add_atom() and/or add_atom_id() are in the proxy-implemented function
-      cvm::load_atoms(atoms_file_name.c_str(), *this, atoms_col, atoms_col_value);
+      parse_error |= cvm::load_atoms(atoms_file_name.c_str(), *this, atoms_col, atoms_col_value);
     }
   }
 
index 9ecd650..cb272ee 100644 (file)
@@ -46,12 +46,24 @@ int colvar::cvc::init(std::string const &conf)
   if (cvm::debug())
     cvm::log("Initializing cvc base object.\n");
 
-  get_keyval(conf, "name", this->name, this->name);
+  std::string const old_name(name);
+
   if (name.size() > 0) {
-    // Temporary description until child object is initialized
-    description = "cvc " + name;
-  } else {
-    description = "uninitialized cvc";
+    cvm::log("Updating configuration for component \""+name+"\"");
+  }
+
+  if (get_keyval(conf, "name", name, name)) {
+    if (name.size() > 0) {
+      description = "cvc \"" + name + "\" of type " + function_type;
+    } else {
+      description = "unnamed cvc";
+    }
+    if ((name != old_name) && (old_name.size() > 0)) {
+      cvm::error("Error: cannot rename component \""+old_name+
+                 "\" after initialization (new name = \""+name+"\")",
+                 INPUT_ERROR);
+      name = old_name;
+    }
   }
 
   get_keyval(conf, "componentCoeff", sup_coeff, sup_coeff);
index 23ecdee..1a6df37 100644 (file)
@@ -82,6 +82,9 @@ public:
   /// this variable definition should be set within the constructor.
   std::string function_type;
 
+  /// Keyword used in the input to denote this CVC
+  std::string config_key;
+
   /// \brief Coefficient in the polynomial combination (default: 1.0)
   cvm::real sup_coeff;
   /// \brief Exponent in the polynomial combination (default: 1)
index 57b2a9a..a2b1a5c 100644 (file)
@@ -939,7 +939,7 @@ colvar::rmsd::rmsd(std::string const &conf)
 
   bool b_Jacobian_derivative = true;
   if (atoms->fitting_group != NULL && b_Jacobian_derivative) {
-    cvm::log("The option \"refPositionsGroup\" (alternative group for fitting) was enabled: "
+    cvm::log("The option \"fittingGroup\" (alternative group for fitting) was enabled: "
               "Jacobian derivatives of the RMSD will not be calculated.\n");
     b_Jacobian_derivative = false;
   }
index dc6242c..87b08b1 100644 (file)
@@ -139,7 +139,7 @@ int colvarmodule::read_config_file(char const  *config_filename)
   // read the config file into a string
   std::string conf = "";
   std::string line;
-  while (colvarparse::getline_nocomments(config_s, line)) {
+  while (parse->read_config_line(config_s, line)) {
     // Delete lines that contain only white space after removing comments
     if (line.find_first_not_of(colvarparse::white_space) != std::string::npos)
       conf.append(line+"\n");
@@ -159,11 +159,12 @@ int colvarmodule::read_config_string(std::string const &config_str)
   // strip the comments away
   std::string conf = "";
   std::string line;
-  while (colvarparse::getline_nocomments(config_s, line)) {
+  while (parse->read_config_line(config_s, line)) {
     // Delete lines that contain only white space after removing comments
     if (line.find_first_not_of(colvarparse::white_space) != std::string::npos)
       conf.append(line+"\n");
   }
+
   return parse_config(conf);
 }
 
@@ -191,6 +192,12 @@ int colvarmodule::parse_config(std::string &conf)
 {
   extra_conf.clear();
 
+  // Check that the input has matching braces
+  if (colvarparse::check_braces(conf, 0) != COLVARS_OK) {
+    return cvm::error("Error: unmatched curly braces in configuration.\n",
+                      INPUT_ERROR);
+  }
+
   // Parse global options
   if (catch_input_errors(parse_global_params(conf))) {
     return get_error();
@@ -235,6 +242,12 @@ int colvarmodule::parse_config(std::string &conf)
 }
 
 
+std::string const & colvarmodule::get_config() const
+{
+  return parse->get_config();
+}
+
+
 int colvarmodule::append_new_config(std::string const &new_conf)
 {
   extra_conf += new_conf;
@@ -246,9 +259,13 @@ int colvarmodule::parse_global_params(std::string const &conf)
 {
   colvarmodule *cvm = cvm::main();
 
-  std::string index_file_name;
-  if (parse->get_keyval(conf, "indexFile", index_file_name)) {
-    cvm->read_index_file(index_file_name.c_str());
+  {
+    std::string index_file_name;
+    size_t pos = 0;
+    while (parse->key_lookup(conf, "indexFile", &index_file_name, &pos)) {
+      cvm->read_index_file(index_file_name.c_str());
+      index_file_name.clear();
+    }
   }
 
   if (parse->get_keyval(conf, "smp", proxy->b_smp_active, proxy->b_smp_active)) {
@@ -1073,10 +1090,10 @@ colvarmodule::~colvarmodule()
 
 int colvarmodule::reset()
 {
-  parse->init();
-
   cvm::log("Resetting the Collective Variables module.\n");
 
+  parse->init();
+
   // Iterate backwards because we are deleting the elements as we go
   for (std::vector<colvarbias *>::reverse_iterator bi = biases.rbegin();
        bi != biases.rend();
index f3d99a2..3d93798 100644 (file)
@@ -131,7 +131,7 @@ public:
 
   /// Module-wide error state
   /// see constants at the top of this file
-protected:
+private:
 
   static int errorCode;
 
@@ -274,6 +274,9 @@ public:
   /// \brief Parse a "clean" config string (no comments)
   int parse_config(std::string &conf);
 
+  /// Get the configuration string read so far (includes comments)
+  std::string const & get_config() const;
+
   // Parse functions (setup internal data based on a string)
 
   /// Allow reading from Windows text files using using std::getline
@@ -296,6 +299,9 @@ public:
 
 private:
 
+  /// Configuration string read so far by the module (includes comments)
+  std::string config_string;
+
   /// Auto-generated configuration during parsing (e.g. to implement
   /// back-compatibility)
   std::string extra_conf;
index 43c8c69..d8b3a35 100644 (file)
@@ -43,9 +43,10 @@ template<typename TYPE> bool colvarparse::_get_keyval_scalar_(std::string const
     }
   } while (b_found);
 
-  if (found_count > 1)
-    cvm::log("Warning: found more than one instance of \""+
-             std::string(key)+"\".\n");
+  if (found_count > 1) {
+    cvm::error("Error: found more than one instance of \""+
+               std::string(key)+"\".\n", INPUT_ERROR);
+  }
 
   if (data.size()) {
     std::istringstream is(data);
@@ -98,9 +99,10 @@ bool colvarparse::_get_keyval_scalar_string_(std::string const &conf,
     }
   } while (b_found);
 
-  if (found_count > 1)
-    cvm::log("Warning: found more than one instance of \""+
-             std::string(key)+"\".\n");
+  if (found_count > 1) {
+    cvm::error("Error: found more than one instance of \""+
+               std::string(key)+"\".\n", INPUT_ERROR);
+  }
 
   if (data.size()) {
     std::istringstream is(data);
@@ -162,9 +164,10 @@ template<typename TYPE> bool colvarparse::_get_keyval_vector_(std::string const
     }
   } while (b_found);
 
-  if (found_count > 1)
-    cvm::log("Warning: found more than one instance of \""+
-             std::string(key)+"\".\n");
+  if (found_count > 1) {
+    cvm::error("Error: found more than one instance of \""+
+               std::string(key)+"\".\n", INPUT_ERROR);
+  }
 
   if (data.size()) {
     std::istringstream is(data);
@@ -319,9 +322,10 @@ bool colvarparse::get_keyval(std::string const &conf,
     }
   } while (b_found);
 
-  if (found_count > 1)
-    cvm::log("Warning: found more than one instance of \""+
-             std::string(key)+"\".\n");
+  if (found_count > 1) {
+    cvm::error("Error: found more than one instance of \""+
+               std::string(key)+"\".\n", INPUT_ERROR);
+  }
 
   if (data.size()) {
     if ( (data == std::string("on")) ||
@@ -535,6 +539,19 @@ int colvarparse::check_keywords(std::string &conf, char const *key)
 }
 
 
+std::istream & colvarparse::read_config_line(std::istream &is,
+                                             std::string &line)
+{
+  cvm::getline(is, line);
+  config_string += line+'\n';
+  size_t const comment = line.find('#');
+  if (comment != std::string::npos) {
+    line.erase(comment);
+  }
+  return is;
+}
+
+
 std::istream & colvarparse::getline_nocomments(std::istream &is,
                                                std::string &line)
 {
@@ -607,7 +624,7 @@ bool colvarparse::key_lookup(std::string const &conf,
     }
 
     // check that there are matching braces between here and the end of conf
-    bool const b_not_within_block = brace_check(conf, pos);
+    bool const b_not_within_block = (check_braces(conf, pos) == COLVARS_OK);
 
     bool const b_isolated = (b_isolated_left && b_isolated_right &&
                              b_not_within_block);
@@ -781,19 +798,15 @@ std::istream & operator>> (std::istream &is, colvarparse::read_block const &rb)
 }
 
 
-bool colvarparse::brace_check(std::string const &conf,
+int colvarparse::check_braces(std::string const &conf,
                               size_t const start_pos)
 {
-  size_t brace_count = 0;
+  int brace_count = 0;
   size_t brace = start_pos;
-  while ( (brace = conf.find_first_of("{}", brace)) != std::string::npos) {
+  while ((brace = conf.find_first_of("{}", brace)) != std::string::npos) {
     if (conf[brace] == '{') brace_count++;
     if (conf[brace] == '}') brace_count--;
     brace++;
   }
-
-  if (brace_count != 0)
-    return false;
-  else
-    return true;
+  return (brace_count != 0) ? INPUT_ERROR : COLVARS_OK;
 }
index 9389bc4..28ad3c0 100644 (file)
@@ -24,7 +24,7 @@
 /// need to parse input inherit from this
 class colvarparse {
 
-private:
+protected:
 
   /// \brief List of legal keywords for this object: this is updated
   /// by each call to colvarparse::get_keyval() or
@@ -47,7 +47,7 @@ private:
   /// \brief Remove all the values from the config string
   void strip_values(std::string &conf);
 
-  /// \brief Configuration string of the object
+  /// \brief Configuration string of the object (includes comments)
   std::string config_string;
 
 public:
@@ -72,7 +72,7 @@ public:
   }
 
   /// Set a new config string for this object
-  inline void init(const std::string& conf)
+  inline void init(std::string const &conf)
   {
     if (! config_string.size()) {
       init();
@@ -80,7 +80,8 @@ public:
     }
   }
 
-  inline const std::string& get_config()
+  /// Get the configuration string (includes comments)
+  inline std::string const & get_config() const
   {
     return config_string;
   }
@@ -284,14 +285,19 @@ public:
                   std::string *data = NULL,
                   size_t *save_pos = NULL);
 
+  /// \brief Reads a configuration line, adds it to config_string, and returns
+  /// the stream \param is Input stream \param s String that will hold the
+  /// configuration line, with comments stripped
+  std::istream & read_config_line(std::istream &is, std::string &line);
+
   /// \brief Works as std::getline() but also removes everything
   /// between a comment character and the following newline
-  static std::istream & getline_nocomments(std::istream &is,
-                                           std::string &s);
+  static std::istream & getline_nocomments(std::istream &is, std::string &s);
 
-  /// Check if the content of the file has matching braces
-  bool brace_check(std::string const &conf,
-                   size_t const start_pos = 0);
+  /// \brief Check if the content of a config string has matching braces
+  /// \param conf The configuration string \param start_pos Start the count
+  /// from this position
+  static int check_braces(std::string const &conf, size_t const start_pos);
 
 };
 
index 86338df..da9257e 100644 (file)
@@ -289,13 +289,23 @@ colvarproxy_smp::colvarproxy_smp()
   omp_lock_state = NULL;
 #if defined(_OPENMP)
   if (smp_thread_id() == 0) {
+    omp_lock_state = reinterpret_cast<void *>(new omp_lock_t);
     omp_init_lock(reinterpret_cast<omp_lock_t *>(omp_lock_state));
   }
 #endif
 }
 
 
-colvarproxy_smp::~colvarproxy_smp() {}
+colvarproxy_smp::~colvarproxy_smp()
+{
+#if defined(_OPENMP)
+  if (smp_thread_id() == 0) {
+    if (omp_lock_state) {
+      delete reinterpret_cast<omp_lock_t *>(omp_lock_state);
+    }
+  }
+#endif
+}
 
 
 int colvarproxy_smp::smp_enabled()
@@ -499,6 +509,14 @@ char const *colvarproxy_script::script_obj_to_str(unsigned char *obj)
 }
 
 
+std::vector<std::string> colvarproxy_script::script_obj_to_str_vector(unsigned char *obj)
+{
+  cvm::error("Error: trying to print a script object without a scripting "
+             "language interface.\n", BUG_ERROR);
+  return std::vector<std::string>();
+}
+
+
 int colvarproxy_script::run_force_callback()
 {
   return COLVARS_NOT_IMPLEMENTED;
index 845f93c..9f636d8 100644 (file)
@@ -184,7 +184,8 @@ public:
   /// Get the mass of the given atom
   inline cvm::real get_atom_mass(int index) const
   {
-    return atoms_masses[index];
+    cvm::real m = atoms_masses[index];
+    return ( m > 0. ? m : 1.0); // Avoid atoms with zero mass
   }
 
   /// Get the charge of the given atom
@@ -461,6 +462,9 @@ public:
   /// Convert a script object (Tcl or Python call argument) to a C string
   virtual char const *script_obj_to_str(unsigned char *obj);
 
+  /// Convert a script object (Tcl or Python call argument) to a vector of strings
+  virtual std::vector<std::string> script_obj_to_str_vector(unsigned char *obj);
+
   /// Pointer to the scripting interface object
   /// (does not need to be allocated in a new interface)
   colvarscript *script;
index 817bb2a..a0f083a 100644 (file)
@@ -1,5 +1,5 @@
 #ifndef COLVARS_VERSION
-#define COLVARS_VERSION "2018-05-29"
+#define COLVARS_VERSION "2018-08-29"
 // This file is part of the Collective Variables module (Colvars).
 // The original version of Colvars and its updates are located at:
 // https://github.com/colvars/colvars
index a55e4c6..366758a 100644 (file)
@@ -82,6 +82,15 @@ int colvarscript::run(int objc, unsigned char *const objv[])
 
   int error_code = COLVARS_OK;
 
+  // If command is found in map, execute it
+  std::string const cmd_key("cv_"+cmd);
+  if (comm_str_map.count(cmd_key) > 0) {
+  cvm::log("cmd_key = "+cmd_key);
+    error_code |= (*(comm_fns[comm_str_map[cmd_key]]))(
+                      reinterpret_cast<void *>(this), objc, objv);
+    return error_code;
+  }
+
   if (cmd == "colvar") {
     if (objc < 3) {
       result = "Missing parameters\n" + help_string();
@@ -295,11 +304,12 @@ int colvarscript::proc_colvar(colvar *cv, int objc, unsigned char *const objv[])
   }
 
   if (subcmd == "delete") {
-    size_t i;
-    for (i = 0; i < cv->biases.size(); i++) {
+    while (cv->biases.size() > 0) {
+      size_t i = cv->biases.size()-1;
+      cvm::log("Warning: before deleting colvar " + cv->name
+        + ", deleting parent bias " + cv->biases[i]->name);
       delete cv->biases[i];
     }
-    cv->biases.clear();
     // colvar destructor is tasked with the cleanup
     delete cv;
     // TODO this could be done by the destructors
@@ -373,6 +383,23 @@ int colvarscript::proc_colvar(colvar *cv, int objc, unsigned char *const objv[])
     return COLVARS_OK;
   }
 
+  if (subcmd == "modifycvcs") {
+    if (objc < 4) {
+      result = "cvcflags: missing parameter: vector of strings";
+      return COLVARSCRIPT_ERROR;
+    }
+    std::vector<std::string> const confs(proxy->script_obj_to_str_vector(objv[3]));
+    cvm::increase_depth();
+    int res = cv->update_cvc_config(confs);
+    cvm::decrease_depth();
+    if (res != COLVARS_OK) {
+      result = "Error setting CVC flags";
+      return COLVARSCRIPT_ERROR;
+    }
+    result = "0";
+    return COLVARS_OK;
+  }
+
   if ((subcmd == "get") || (subcmd == "set") || (subcmd == "state")) {
     return proc_features(cv, objc, objv);
   }
@@ -547,6 +574,8 @@ std::string colvarscript::help_string() const
 Managing the Colvars module:\n\
   configfile <file name>      -- read configuration from a file\n\
   config <string>             -- read configuration from the given string\n\
+  getconfig                   -- get the module's configuration string\n\
+  resetindexgroups            -- clear the index groups loaded so far\n\
   reset                       -- delete all internal configuration\n\
   delete                      -- delete this Colvars module instance\n\
   version                     -- return version of Colvars code\n\
@@ -579,6 +608,7 @@ Accessing collective variables:\n\
   colvar <name> gettotalforce -- return total force of colvar <name>\n\
   colvar <name> getconfig     -- return config string of colvar <name>\n\
   colvar <name> cvcflags <fl> -- enable or disable cvcs according to 0/1 flags\n\
+  colvar <name> modifycvcs <str> -- pass new config strings to each CVC\n\
   colvar <name> get <f>       -- get the value of the colvar feature <f>\n\
   colvar <name> set <f> <val> -- set the value of the colvar feature <f>\n\
 \n\
index 39cd089..313dbd6 100644 (file)
@@ -71,8 +71,10 @@ public:
     cv_help,
     cv_version,
     cv_config,
+    cv_getconfig,
     cv_configfile,
     cv_reset,
+    cv_resetindexgroups,
     cv_delete,
     cv_list,
     cv_list_biases,
@@ -254,6 +256,23 @@ extern "C" {
            return COLVARSCRIPT_ERROR;
            )
 
+  CVSCRIPT(cv_getconfig,
+           "Get the module's configuration string read so far",
+           0, 0,
+           { },
+           script->set_str_result(cvm::main()->get_config());
+           return COLVARS_OK;
+           )
+
+  CVSCRIPT(cv_resetindexgroups,
+           "Clear the index groups loaded so far, allowing to replace them",
+           0, 0,
+           { },
+           cvm::main()->index_group_names.clear();
+           cvm::main()->index_groups.clear();
+           return COLVARS_OK;
+           )
+
   CVSCRIPT(cv_addenergy,
            "Add an energy to the MD engine",
            1, 1,
index 3cdc5e8..22e1f6e 100644 (file)
@@ -548,7 +548,7 @@ void colvarproxy_namd::fatal_error(std::string const &message)
 {
   log(message);
   if (errno) {
-    log(strerror(errno));
+    // log(strerror(errno));
     NAMD_err("Error in the collective variables module");
   } else {
     NAMD_die("Error in the collective variables module: exiting.\n");
@@ -827,13 +827,14 @@ int colvarproxy_namd::load_coords(char const *pdb_filename,
         break;
     }
 
-    if ((ipos < pos.size()) || (current_index != indices.end()))
-      cvm::error("Error: the number of records in the PDB file \""+
-                 std::string(pdb_filename)+
-                 "\" does not appear to match either the total number of atoms,"+
-                 " or the number of coordinates requested at this point("+
-                 cvm::to_str(pos.size())+").\n", BUG_ERROR);
-
+    if (ipos < pos.size() || (!use_pdb_field && current_index != indices.end())) {
+      size_t n_requested = use_pdb_field ? pos.size() : indices.size();
+      cvm::error("Error: number of matching records in the PDB file \""+
+                 std::string(pdb_filename)+"\" ("+cvm::to_str(ipos)+
+                 ") does not match the number of requested coordinates ("+
+                 cvm::to_str(n_requested)+").\n", INPUT_ERROR);
+      return COLVARS_ERROR;
+    }
   } else {
 
     // when the PDB contains exactly the number of atoms of the array,
@@ -992,6 +993,33 @@ char const *colvarproxy_namd::script_obj_to_str(unsigned char *obj)
 #endif
 }
 
+
+std::vector<std::string> colvarproxy_namd::script_obj_to_str_vector(unsigned char *obj)
+{
+  if (cvm::debug()) {
+    cvm::log("Called colvarproxy_namd::script_obj_to_str_vector().\n");
+  }
+  std::vector<std::string> result;
+#ifdef NAMD_TCL
+  Tcl_Interp *interp = reinterpret_cast<Tcl_Interp *>(_tcl_interp);
+  Tcl_Obj *tcl_obj = reinterpret_cast<Tcl_Obj *>(obj);
+  Tcl_Obj **tcl_list_elems = NULL;
+  int count = 0;
+  if (Tcl_ListObjGetElements(interp, tcl_obj, &count, &tcl_list_elems) ==
+      TCL_OK) {
+    result.reserve(count);
+    for (int i = 0; i < count; i++) {
+      result.push_back(Tcl_GetString(tcl_list_elems[i]));
+    }
+  } else {
+    Tcl_SetResult(interp,
+                  const_cast<char *>("Cannot parse Tcl list."), TCL_STATIC);
+  }
+#endif
+  return result;
+}
+
+
 int colvarproxy_namd::init_atom_group(std::vector<int> const &atoms_ids)
 {
   if (cvm::debug())
index 15279ee..36e7bbb 100644 (file)
@@ -283,6 +283,7 @@ public:
   int backup_file(char const *filename);
 
   char const *script_obj_to_str(unsigned char *obj);
+  std::vector<std::string> script_obj_to_str_vector(unsigned char *obj);
 };
 
 
index 0d45121..edd6778 100644 (file)
@@ -1,5 +1,5 @@
 #ifndef COLVARPROXY_VERSION
-#define COLVARPROXY_VERSION "2018-04-29"
+#define COLVARPROXY_VERSION "2018-08-29"
 // This file is part of the Collective Variables module (Colvars).
 // The original version of Colvars and its updates are located at:
 // https://github.com/colvars/colvars
index 82b9223..136d2bb 100644 (file)
@@ -91,8 +91,11 @@ The following sections explain the syntax of the various sub-commands of \texttt
 \item \texttt{colvar} \emph{$<$name$>$} \texttt{gettotalforce}: get the total force acting on colvar \emph{$<$name$>$} at the previous step (see \ref{sec:cvc_sys_forces});
 }
 \item \texttt{colvar} \emph{$<$name$>$} \texttt{getconfig}: return config string of colvar \emph{$<$name$>$}.
-\item \texttt{colvar} \emph{$<$name$>$} \texttt{cvcflags} \emph{$<$flags$>$}: for a colvar with several cvcs (numbered according to their name
-string order), set which cvcs are enabled or disabled in subsequent evaluations according to a list of 0/1 flags (one per cvc).
+\item \texttt{colvar} \emph{$<$name$>$} \texttt{cvcflags} \emph{$<$flags$>$}: for a colvar with several CVCs (numbered according to their name string order), set which CVCs are enabled or disabled in subsequent evaluations according to a list of 0/1 flags (one per CVC).
+\item \texttt{colvar} \emph{$<$name$>$} \texttt{modifycvcs} \emph{$<$strings$>$}: for a colvar with one or more CVCs (numbered according to their name string order), pass a list of new configuration strings to modify each CVC without needing to delete the colvar.
+This option is currently limited to changing the values of \refkey{componentCoeff}{sec:cvc_superp} and \refkey{componentExp}{sec:cvc_superp} (e.g.{} to update the polynomial superposition parameters on the fly), of \refkey{period}{sec:cvc_periodic} and \refkey{wrapAround}{sec:cvc_periodic} (e.g.{} to update the period of variables such as \texttt{distanceZ}), and of the \texttt{forceNoPBC} option for those CVCs that support it.
+Changes in configuration done by this command are not saved to state files, and are lost when restarting a simulation, deleting the parent colvar, or resetting the module with \cvvmdonly{\texttt{cv delete} or} \texttt{cv reset}.
+
 \end{itemize}
 
 \cvparagraph{Accessing biases}{sec:cv_command_bias}
index 7eeb4b9..e84a568 100644 (file)
@@ -1072,9 +1072,15 @@ Other components are instead affected by global rotations or translations: howev
 Finally, a few components are defined by convention using a roto-translated frame (e.g. the minimal RMSD): for these components, \texttt{centerReference} and \texttt{rotateReference} are enabled by default.
 In typical applications, the default settings result in the expected behavior.
 
+\paragraph*{Warning on rotating frames of reference and periodic boundary conditions.}
+\texttt{rotateReference} affects coordinates that depend on minimum-image distances in periodic boundary conditions (PBC).
+After rotation of the coordinates, the periodic cell vectors become irrelevant: the rotated system is effectively non-periodic.
+A safe way to handle this is to ensure that the relevant inter-group distance vectors remain smaller than the half-size of the periodic cell.
+If this is not desirable, one should avoid the rotating frame of reference, and apply orientational restraints to the reference group instead, in order to keep the orientation of the reference group consistent with the orientation of the periodic cell.
+
+\paragraph*{Warning on rotating frames of reference and ABF.}
 Note that \texttt{centerReference} and \texttt{rotateReference} may affect the Jacobian derivative of colvar components in a way that is not taken into account by default.
 Be careful when using these options in ABF simulations or when using total force values.
-At this point, only the RMSD component will correctly adjust the Jacobian derivative computation to all situations.
 
 \begin{itemize}
 
@@ -1704,8 +1710,9 @@ spherical coordinates, for the center of mass of a group of atoms
 described by the block \texttt{atoms}.  It returns an angle
 (in degrees) within the interval $[0:180]$.
 To obtain spherical coordinates in a frame of reference tied to
-another group of atoms, use the \texttt{fittingGroup}\ref{sec:colvar_atom_groups_ref_frame} option
+another group of atoms, use the \texttt{fittingGroup} (\ref{sec:colvar_atom_groups_ref_frame}) option
 within the \texttt{atoms} block.
+An example is provided in file \texttt{examples/11\_polar\_angles.in} of the Colvars public repository.
 
 \begin{cvcoptions}
 \item %
@@ -1728,6 +1735,11 @@ described by the block \texttt{atoms}. It returns an angle
 calculates all the distances between two angles taking into account
 periodicity.  For instance, reference values for restraints or range
 boundaries can be defined by using any real number of choice.
+To obtain spherical coordinates in a frame of reference tied to
+another group of atoms, use the \texttt{fittingGroup} (\ref{sec:colvar_atom_groups_ref_frame}) option
+within the \texttt{atoms} block.
+An example is provided in file \texttt{examples/11\_polar\_angles.in} of the Colvars public repository.
+
 
 \begin{cvcoptions}
 \item %
@@ -2065,6 +2077,13 @@ reflects the deviation from reference coordinates in a separate, moving
 reference frame.
 \end{enumerate}
 
+\cvscriptonly{
+\cvsubsubsec{Path collective variables}{sec:pathcv}
+
+An application of the \texttt{rmsd} component is "path collective variables",\cite{Branduardi2007}
+which are implemented as Tcl-scripted combinations or RMSDs.
+The implementation is available as file \texttt{colvartools/pathCV.tcl}, and
+an example is provided in file \texttt{examples/10\_pathCV.namd} of the Colvars public repository.}
 
 \cvsubsubsec{\texttt{eigenvector}: projection of the atomic  coordinates on a vector.}{sec:cvc_eigenvector}
 The block \texttt{eigenvector~\{...\}} defines the projection of the coordinates
@@ -2676,14 +2695,20 @@ In certain conditions, \texttt{distanceZ} can also be periodic, namely
 when periodic boundary conditions (PBCs) are defined in the simulation
 and \texttt{distanceZ}'s axis is parallel to a unit cell vector.
 
-The following keywords can be used within periodic components (and are
-illegal elsewhere):
+In addition, a custom \cvscriptonly{or scripted} scalar colvar may be periodic
+depending on its user-defined expression. It will only be treated as such by
+the Colvars module if the period is specified using the \texttt{period} keyword,
+while \texttt{wrapAround} is optional.
+
+The following keywords can be used within periodic components, and within
+custom \cvscriptonly{or scripted} colvars (\ref{sec:colvar_custom_function}
+\cvscriptonly{, \ref{sec::colvar_scripted}}).
 
 \begin{itemize}
 \item %
   \keydef
     {period}{%
-    \texttt{distanceZ}}{%
+    \texttt{distanceZ}, custom colvars}{%
     Period of the component}{%
     positive decimal}{%
     0.0}{%
@@ -2696,24 +2721,20 @@ illegal elsewhere):
 \item %
   \keydef
     {wrapAround}{%
-    \texttt{distanceZ}, \texttt{dihedral} or \texttt{spinAngle}}{%
-    Center
-    of the wrapping interval for periodic variables}{%
+    \texttt{distanceZ}, \texttt{dihedral}, \texttt{spinAngle}, custom colvars}{%
+    Center of the wrapping interval for periodic variables}{%
     decimal}{%
     0.0}{%
     By default, values of the periodic components are centered around zero, ranging from $-P/2$ to $P/2$, where $P$ is the period.
     Setting this number centers the interval around this value.
     This can be useful for convenience of output, or to set the walls for a \texttt{harmonicWalls} in an order that would not otherwise be allowed.}
 \end{itemize}
+
 Internally, all differences between two values of a periodic colvar
 follow the minimum image convention: they are calculated based on
 the two periodic images that are closest to each other.
 
-\emph{Note: linear or polynomial combinations of periodic components
-  may become meaningless when components cross the periodic boundary.
-  Use such combinations carefully: estimate the range of possible values
-  of each component in a given simulation, and make use of
-  \texttt{wrapAround} to limit this problem whenever possible.}
+\emph{Note: linear or polynomial combinations of periodic components (see \ref{sec:cvc_superp}) may become meaningless when components cross the periodic boundary.  Use such combinations carefully: estimate the range of possible values of each component in a given simulation, and make use of \texttt{wrapAround} to limit this problem whenever possible.}
 
 
 \cvsubsubsec{Non-scalar components.}{sec:cvc_non_scalar}
@@ -2890,8 +2911,8 @@ the keyword \texttt{period}, and the optional keyword \texttt{wrapAround}, with
 same meaning as in periodic components (see \ref{sec:cvc_periodic} for details).
 
 An example of elaborate scripted colvar is given in example 10, in the
-form of path-based collective variables as defined by Branduardi et al\cite{Branduardi2007}.
-The required Tcl procedures are provided in the \texttt{colvartools} directory.
+form of path-based collective variables as defined by Branduardi et al\cite{Branduardi2007}
+(\ref{sec:pathcv}).
 
 \begin{itemize}
  \item %
@@ -3539,7 +3560,7 @@ Calculate UI estimator of the free energy?}
 When enabled, it triggers calculation of the free energy following the UI estimator.}
 \end{itemize}
 
-\paragraph{Usage for multiple--replica eABF.}
+\paragraph*{Usage for multiple--replica eABF.}
 The eABF algorithm can be associated with a multiple--walker strategy \cite{Minoukadeh2010,Comer2014c} (\ref{sec:mw-ABF}).
 To run a multiple--replica eABF simulation, start a multiple-replica
 NAMD run (option \texttt{+replicas}) and set {\tt shared on} in the Colvars config file to enable