Difference for src/colvar.C from version 1.34 to 1.35

version 1.34version 1.35
Line 1
Line 1
  
 // -*- c++ -*- // -*- c++ -*-
  
 // This file is part of the Collective Variables module (Colvars). // This file is part of the Collective Variables module (Colvars).
 // The original version of Colvars and its updates are located at: // The original version of Colvars and its updates are located at:
 // https://github.com/colvars/colvars // https://github.com/colvars/colvars
Line 7
Line 7
 // If you wish to distribute your changes, please submit them to the // If you wish to distribute your changes, please submit them to the
 // Colvars repository at GitHub. // Colvars repository at GitHub.
  
  
 #include "colvarmodule.h" #include "colvarmodule.h"
 #include "colvarvalue.h" #include "colvarvalue.h"
 #include "colvarparse.h" #include "colvarparse.h"
 #include "colvar.h" #include "colvar.h"
 #include "colvarcomp.h" #include "colvarcomp.h"
 #include "colvarscript.h" #include "colvarscript.h"
  
  // used in build_atom_list()
 #include <algorithm> #include <algorithm>
  
  
Line 25
Line 26
  
  
 colvar::colvar() colvar::colvar()
    : prev_timestep(-1)
 { {
   // Initialize static array once and for all   // Initialize static array once and for all
   init_cv_requires();   init_cv_requires();
Line 223
Line 225
  
   // at this point, the colvar's type is defined   // at this point, the colvar's type is defined
   f.type(value());   f.type(value());
   f_accumulated.type(value()); 
  
   x_old.type(value());   x_old.type(value());
   v_fdiff.type(value());   v_fdiff.type(value());
Line 239
Line 240
  
   reset_bias_force();   reset_bias_force();
  
    get_keyval(conf, "timeStepFactor", time_step_factor, 1);
    if (time_step_factor < 0) {
      cvm::error("Error: timeStepFactor must be positive.\n");
      return COLVARS_ERROR;
    }
    if (time_step_factor != 1) {
      enable(f_cv_multiple_ts);
    }
  
   // TODO use here information from the CVCs' own natural boundaries   // TODO use here information from the CVCs' own natural boundaries
   error_code |= init_grid_parameters(conf);   error_code |= init_grid_parameters(conf);
  
   get_keyval(conf, "timeStepFactor", time_step_factor, 1); 
  
   error_code |= init_extended_Lagrangian(conf);   error_code |= init_extended_Lagrangian(conf);
   error_code |= init_output_flags(conf);   error_code |= init_output_flags(conf);
  
   // Start in active state by default   // Now that the children are defined we can solve dependencies
   enable(f_cv_active);   enable(f_cv_active);
   // Make sure dependency side-effects are correct 
   refresh_deps(); 
  
   if (cvm::b_analysis)   if (cvm::b_analysis)
     parse_analysis(conf);     parse_analysis(conf);
Line 326
Line 332
       std::string const walls_conf("\n\       std::string const walls_conf("\n\
 harmonicWalls {\n\ harmonicWalls {\n\
     name "+this->name+"w\n\     name "+this->name+"w\n\
     colvars "+this->name+"\n"+lw_conf+uw_conf+     colvars "+this->name+"\n"+lw_conf+uw_conf+"\
      timeStepFactor "+cvm::to_str(time_step_factor)+"\n"+
                              "}\n");                              "}\n");
       cv->append_new_config(walls_conf);       cv->append_new_config(walls_conf);
     }     }
Line 511
Line 518
     if ( (cvcp->period != 0.0) || (cvcp->wrap_center != 0.0) ) {     if ( (cvcp->period != 0.0) || (cvcp->wrap_center != 0.0) ) {
       if ( (cvcp->function_type != std::string("distance_z")) &&       if ( (cvcp->function_type != std::string("distance_z")) &&
            (cvcp->function_type != std::string("dihedral")) &&            (cvcp->function_type != std::string("dihedral")) &&
             (cvcp->function_type != std::string("polar_phi")) &&
            (cvcp->function_type != std::string("spin_angle")) ) {            (cvcp->function_type != std::string("spin_angle")) ) {
         cvm::error("Error: invalid use of period and/or "         cvm::error("Error: invalid use of period and/or "
                    "wrapAround in a \""+                    "wrapAround in a \""+
Line 563
Line 571
     "on an axis", "distanceZ");     "on an axis", "distanceZ");
   error_code |= init_components_type<distance_xy>(conf, "distance projection "   error_code |= init_components_type<distance_xy>(conf, "distance projection "
     "on a plane", "distanceXY");     "on a plane", "distanceXY");
    error_code |= init_components_type<polar_theta>(conf, "spherical polar angle theta",
      "polarTheta");
    error_code |= init_components_type<polar_phi>(conf, "spherical azimuthal angle phi",
      "polarPhi");
   error_code |= init_components_type<distance_inv>(conf, "average distance "   error_code |= init_components_type<distance_inv>(conf, "average distance "
     "weighted by inverse power", "distanceInv");     "weighted by inverse power", "distanceInv");
   error_code |= init_components_type<distance_pairs>(conf, "N1xN2-long vector "   error_code |= init_components_type<distance_pairs>(conf, "N1xN2-long vector "
Line 615
Line 627
 } }
  
  
 int colvar::refresh_deps() void colvar::do_feature_side_effects(int id)
 { {
   // If enabled features are changed upstream, the features below should be refreshed   // If enabled features are changed upstream, the features below should be refreshed
   if (is_enabled(f_cv_total_force_calc)) {   switch (id) {
      case f_cv_total_force_calc:
     cvm::request_total_force();     cvm::request_total_force();
   }       break;
   if (is_enabled(f_cv_collect_gradient) && atom_ids.size() == 0) {     case f_cv_collect_gradient:
        if (atom_ids.size() == 0) {
     build_atom_list();     build_atom_list();
   }   }
   return COLVARS_OK;       break;
    }
 } }
  
  
Line 765
Line 780
  
 colvar::~colvar() colvar::~colvar()
 { {
    // There is no need to call free_children_deps() here
    // because the children are cvcs and will be deleted
    // just below
  
 //   Clear references to this colvar's cvcs as children //   Clear references to this colvar's cvcs as children
 //   for dependency purposes //   for dependency purposes
   remove_all_children();   remove_all_children();
Line 981
Line 1000
       (cvcs[i])->calc_gradients();       (cvcs[i])->calc_gradients();
       // if requested, propagate (via chain rule) the gradients above       // if requested, propagate (via chain rule) the gradients above
       // to the atoms used to define the roto-translation       // to the atoms used to define the roto-translation
       // This could be integrated in the CVC base class      (cvcs[i])->calc_fit_gradients();
       for (size_t ig = 0; ig < cvcs[i]->atom_groups.size(); ig++) {       if ((cvcs[i])->is_enabled(f_cvc_debug_gradient))
         if (cvcs[i]->atom_groups[ig]->b_fit_gradients)         (cvcs[i])->debug_gradients();
           cvcs[i]->atom_groups[ig]->calc_fit_gradients(); 
  
         if (cvcs[i]->is_enabled(f_cvc_debug_gradient)) { 
           cvm::log("Debugging gradients for " + cvcs[i]->description); 
           cvcs[i]->debug_gradients(cvcs[i]->atom_groups[ig]); 
         } 
       } 
     }     }
  
     cvm::decrease_depth();     cvm::decrease_depth();
Line 1008
Line 1020
   size_t i;   size_t i;
  
   if (is_enabled(f_cv_collect_gradient)) {   if (is_enabled(f_cv_collect_gradient)) {
  
     if (is_enabled(f_cv_scripted)) { 
       cvm::error("Collecting atomic gradients is not implemented for " 
                  "scripted colvars.", COLVARS_NOT_IMPLEMENTED); 
       return COLVARS_NOT_IMPLEMENTED; 
     } 
  
     // Collect the atomic gradients inside colvar object     // Collect the atomic gradients inside colvar object
     for (unsigned int a = 0; a < atomic_gradients.size(); a++) {     for (unsigned int a = 0; a < atomic_gradients.size(); a++) {
       atomic_gradients[a].reset();       atomic_gradients[a].reset();
Line 1212
Line 1217
   f.type(value());   f.type(value());
   f.reset();   f.reset();
  
    // If we are not active at this timestep, that's all we have to do
    // return with energy == zero
    if (!is_enabled(f_cv_active)) return 0.;
  
   // add the biases' force, which at this point should already have   // add the biases' force, which at this point should already have
   // been summed over each bias using this colvar   // been summed over each bias using this colvar
   f += fb;   f += fb;
Line 1234
Line 1243
     }     }
  
     cvm::real dt = cvm::dt();     cvm::real dt = cvm::dt();
  
      if (prev_timestep > -1) {
        // Keep track of slow timestep to integrate MTS colvars
        // the colvar checks the interval after waking up twice
        int n_timesteps = cvm::step_relative() - prev_timestep;
        if (n_timesteps != 0 && n_timesteps != time_step_factor) {
          cvm::error("Error: extended-Lagrangian " + description + " has timeStepFactor " +
            cvm::to_str(time_step_factor) + ", but was activated after " + cvm::to_str(n_timesteps) +
            " steps at timestep " + cvm::to_str(cvm::step_absolute()) + " (relative step: " +
            cvm::to_str(cvm::step_relative()) + ").\n" +
            "Make sure that this colvar is requested by biases at multiples of timeStepFactor.\n");
          return 0.;
        }
      }
      prev_timestep = cvm::step_relative();
  
      // Integrate with slow timestep
      dt *= cvm::real(time_step_factor);
  
     colvarvalue f_ext(fr.type()); // force acting on the extended variable     colvarvalue f_ext(fr.type()); // force acting on the extended variable
     f_ext.reset();     f_ext.reset();
  
Line 1280
Line 1308
   // bypass the extended Lagrangian mass)   // bypass the extended Lagrangian mass)
   f += fb_actual;   f += fb_actual;
  
   // Store force to be applied, possibly summed over several timesteps 
   f_accumulated += f; 
  
   if (is_enabled(f_cv_fdiff_velocity)) {   if (is_enabled(f_cv_fdiff_velocity)) {
     // set it for the next step     // set it for the next step
     x_old = x;     x_old = x;
Line 1303
Line 1328
   size_t i;   size_t i;
   if (cvm::debug()) {   if (cvm::debug()) {
     cvm::log("Communicating forces from colvar \""+this->name+"\".\n");     cvm::log("Communicating forces from colvar \""+this->name+"\".\n");
     cvm::log("Force to be applied: " + cvm::to_str(f_accumulated) + "\n");     cvm::log("Force to be applied: " + cvm::to_str(f) + "\n");
   }   }
  
   if (is_enabled(f_cv_scripted)) {   if (is_enabled(f_cv_scripted)) {
Line 1330
Line 1355
       if (!cvcs[i]->is_enabled()) continue;       if (!cvcs[i]->is_enabled()) continue;
       // cvc force is colvar force times colvar/cvc Jacobian       // cvc force is colvar force times colvar/cvc Jacobian
       // (vector-matrix product)       // (vector-matrix product)
       (cvcs[i])->apply_force(colvarvalue(f_accumulated.as_vector() * func_grads[grad_index++],       (cvcs[i])->apply_force(colvarvalue(f.as_vector() * func_grads[grad_index++],
                              cvcs[i]->value().type()));                              cvcs[i]->value().type()));
     }     }
   } else if (x.type() == colvarvalue::type_scalar) {   } else if (x.type() == colvarvalue::type_scalar) {
  
     for (i = 0; i < cvcs.size(); i++) {     for (i = 0; i < cvcs.size(); i++) {
       if (!cvcs[i]->is_enabled()) continue;       if (!cvcs[i]->is_enabled()) continue;
       (cvcs[i])->apply_force(f_accumulated * (cvcs[i])->sup_coeff *       (cvcs[i])->apply_force(f * (cvcs[i])->sup_coeff *
                               cvm::real((cvcs[i])->sup_np) *                               cvm::real((cvcs[i])->sup_np) *
                               (std::pow((cvcs[i])->value().real_value,                               (std::pow((cvcs[i])->value().real_value,
                                       (cvcs[i])->sup_np-1)) );                                       (cvcs[i])->sup_np-1)) );
Line 1347
Line 1372
  
     for (i = 0; i < cvcs.size(); i++) {     for (i = 0; i < cvcs.size(); i++) {
       if (!cvcs[i]->is_enabled()) continue;       if (!cvcs[i]->is_enabled()) continue;
       (cvcs[i])->apply_force(f_accumulated * (cvcs[i])->sup_coeff);       (cvcs[i])->apply_force(f * (cvcs[i])->sup_coeff);
     }     }
   }   }
  
   // Accumulated forces have been applied, impulse-style 
   // Reset to start accumulating again 
   f_accumulated.reset(); 
  
   if (cvm::debug())   if (cvm::debug())
     cvm::log("Done communicating forces from colvar \""+this->name+"\".\n");     cvm::log("Done communicating forces from colvar \""+this->name+"\".\n");
 } }
Line 1391
Line 1412
       cvm::error("ERROR: All CVCs are disabled for colvar " + this->name +"\n");       cvm::error("ERROR: All CVCs are disabled for colvar " + this->name +"\n");
       return COLVARS_ERROR;       return COLVARS_ERROR;
     }     }
     cvc_flags.resize(0);     cvc_flags.clear();
   }   }
  
   return COLVARS_OK;   return COLVARS_OK;

Legend:
Removed in v.1.34 
changed lines
 Added in v.1.35


Made by using version 1.53 of cvs2html