Difference for src/colvar.C from version 1.29 to 1.30

version 1.29version 1.30
Line 212
Line 212
   get_keyval(conf, "width", width, 1.0);   get_keyval(conf, "width", width, 1.0);
   if (width <= 0.0) {   if (width <= 0.0) {
     cvm::error("Error: \"width\" must be positive.\n", INPUT_ERROR);     cvm::error("Error: \"width\" must be positive.\n", INPUT_ERROR);
      return;
   }   }
  
   // NOTE: not porting wall stuff to new deps, as this will change to a separate bias   // NOTE: not porting wall stuff to new deps, as this will change to a separate bias
Line 761
Line 762
     return error_code;     return error_code;
   }   }
  
    if (cvm::step_relative() > 0) {
      // Total force depends on Jacobian derivative from previous timestep
      error_code |= calc_cvc_total_force(first_cvc, num_cvcs);
    }
    // atom coordinates are updated by the next line
   error_code |= calc_cvc_values(first_cvc, num_cvcs);   error_code |= calc_cvc_values(first_cvc, num_cvcs);
   error_code |= calc_cvc_gradients(first_cvc, num_cvcs);   error_code |= calc_cvc_gradients(first_cvc, num_cvcs);
   error_code |= calc_cvc_total_force(first_cvc, num_cvcs); 
   error_code |= calc_cvc_Jacobians(first_cvc, num_cvcs);   error_code |= calc_cvc_Jacobians(first_cvc, num_cvcs);
  
   if (cvm::debug())   if (cvm::debug())
Line 780
Line 785
  
   int error_code = COLVARS_OK;   int error_code = COLVARS_OK;
  
    if (cvm::step_relative() > 0) {
      // Total force depends on Jacobian derivative from previous timestep
      error_code |= collect_cvc_total_forces();
    }
   error_code |= collect_cvc_values();   error_code |= collect_cvc_values();
   error_code |= collect_cvc_gradients();   error_code |= collect_cvc_gradients();
   error_code |= collect_cvc_total_forces(); 
   error_code |= collect_cvc_Jacobians();   error_code |= collect_cvc_Jacobians();
   error_code |= calc_colvar_properties();   error_code |= calc_colvar_properties();
  
Line 982
Line 990
     // if (!tasks[task_extended_lagrangian] && (cvm::step_relative() > 0)) {     // if (!tasks[task_extended_lagrangian] && (cvm::step_relative() > 0)) {
    // Disabled check to allow for explicit total force calculation    // Disabled check to allow for explicit total force calculation
     // even with extended Lagrangian     // even with extended Lagrangian
  
     if (cvm::step_relative() > 0) { 
       cvm::increase_depth();       cvm::increase_depth();
       // get from the cvcs the total forces from the PREVIOUS step 
       for (i = first_cvc, cvc_count = 0;       for (i = first_cvc, cvc_count = 0;
           (i < cvcs.size()) && (cvc_count < cvc_max_count);           (i < cvcs.size()) && (cvc_count < cvc_max_count);
           i++) {           i++) {
Line 994
Line 1000
         (cvcs[i])->calc_force_invgrads();         (cvcs[i])->calc_force_invgrads();
       }       }
       cvm::decrease_depth();       cvm::decrease_depth();
     } 
  
     if (cvm::debug())     if (cvm::debug())
       cvm::log("Done calculating total force of colvar \""+this->name+"\".\n");       cvm::log("Done calculating total force of colvar \""+this->name+"\".\n");
Line 1128
Line 1134
   if (is_enabled(f_cv_Jacobian)) {   if (is_enabled(f_cv_Jacobian)) {
     // the instantaneous Jacobian force was not included in the reported total force;     // the instantaneous Jacobian force was not included in the reported total force;
     // instead, it is subtracted from the applied force (silent Jacobian correction)     // instead, it is subtracted from the applied force (silent Jacobian correction)
      // This requires the Jacobian term for the *current* timestep
     if (is_enabled(f_cv_hide_Jacobian))     if (is_enabled(f_cv_hide_Jacobian))
       f -= fj;       f -= fj;
   }   }
  
   if (is_enabled(f_cv_lower_wall) || is_enabled(f_cv_upper_wall)) { 
  
     // Wall force     // Wall force
     colvarvalue fw(x);     colvarvalue fw(x);
     fw.reset();     fw.reset();
  
    if (is_enabled(f_cv_lower_wall) || is_enabled(f_cv_upper_wall)) {
  
     if (cvm::debug())     if (cvm::debug())
       cvm::log("Calculating wall forces for colvar \""+this->name+"\".\n");       cvm::log("Calculating wall forces for colvar \""+this->name+"\".\n");
  
     // For a periodic colvar, both walls may be applicable at the same time     // For a periodic colvar, both walls may be applicable at the same time
     // in which case we pick the closer one     // in which case we pick the closer one
     if ( (!is_enabled(f_cv_upper_wall)) ||     if ( (!is_enabled(f_cv_upper_wall)) ||
          (this->dist2(x_reported, lower_wall) < this->dist2(x_reported, upper_wall)) ) {          (this->dist2(x, lower_wall) < this->dist2(x, upper_wall)) ) {
  
       cvm::real const grad = this->dist2_lgrad(x_reported, lower_wall);       cvm::real const grad = this->dist2_lgrad(x, lower_wall);
       if (grad < 0.0) {       if (grad < 0.0) {
         fw = -0.5 * lower_wall_k * grad;         fw = -0.5 * lower_wall_k * grad;
         f += fw; 
         if (cvm::debug())         if (cvm::debug())
           cvm::log("Applying a lower wall force("+           cvm::log("Applying a lower wall force("+
                     cvm::to_str(fw)+") to \""+this->name+"\".\n");                     cvm::to_str(fw)+") to \""+this->name+"\".\n");
Line 1157
Line 1163
  
     } else {     } else {
  
       cvm::real const grad = this->dist2_lgrad(x_reported, upper_wall);       cvm::real const grad = this->dist2_lgrad(x, upper_wall);
       if (grad > 0.0) {       if (grad > 0.0) {
         fw = -0.5 * upper_wall_k * grad;         fw = -0.5 * upper_wall_k * grad;
         f += fw; 
         if (cvm::debug())         if (cvm::debug())
           cvm::log("Applying an upper wall force("+           cvm::log("Applying an upper wall force("+
                     cvm::to_str(fw)+") to \""+this->name+"\".\n");                     cvm::to_str(fw)+") to \""+this->name+"\".\n");
Line 1168
Line 1173
     }     }
   }   }
  
    // At this point f is the force f from external biases that will be applied to the
    // extended variable if there is one
  
   if (is_enabled(f_cv_extended_Lagrangian)) {   if (is_enabled(f_cv_extended_Lagrangian)) {
  
     cvm::real dt = cvm::dt();     cvm::real dt = cvm::dt();
Line 1175
Line 1183
     f_ext.reset();     f_ext.reset();
  
     // the total force is applied to the fictitious mass, while the     // the total force is applied to the fictitious mass, while the
     // atoms only feel the harmonic force     // atoms only feel the harmonic force + wall force
     // fr: bias force on extended variable (without harmonic spring), for output in trajectory     // fr: bias force on extended variable (without harmonic spring), for output in trajectory
     // f_ext: total force on extended variable (including harmonic spring)     // f_ext: total force on extended variable (including harmonic spring)
     // f: - initially, external biasing force (including wall forces)     // f: - initially, external biasing force
     //    - after this code block, colvar force to be applied to atomic coordinates, ie. spring force     //    - after this code block, colvar force to be applied to atomic coordinates
      //      ie. spring force + wall force
     fr    = f;     fr    = f;
     f_ext = f + (-0.5 * ext_force_k) * this->dist2_lgrad(xr, x);     f_ext = f + (-0.5 * ext_force_k) * this->dist2_lgrad(xr, x);
     f     =     (-0.5 * ext_force_k) * this->dist2_rgrad(xr, x);     f     =     (-0.5 * ext_force_k) * this->dist2_rgrad(xr, x);
Line 1209
Line 1218
     vr  += (0.5 * dt) * f_ext / ext_mass;     vr  += (0.5 * dt) * f_ext / ext_mass;
     xr  += dt * vr;     xr  += dt * vr;
     xr.apply_constraints();     xr.apply_constraints();
     if (this->b_periodic) this->wrap(xr);     if (this->is_enabled(f_cv_periodic)) this->wrap(xr);
   }   }
  
    f += fw;
  
    // Store force to be applied, possibly summed over several timesteps
   f_accumulated += f;   f_accumulated += f;
  
   if (is_enabled(f_cv_fdiff_velocity)) {   if (is_enabled(f_cv_fdiff_velocity)) {


Legend:
Removed in v.1.29 
changed lines
 Added in v.1.30



Made by using version 1.53 of cvs2html