| version 1.29 | version 1.30 |
|---|
| |
| 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 |
| |
| 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()) |
| |
| | |
| 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(); |
| | |
| |
| // 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++) { |
| |
| (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"); |
| |
| 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"); |
| |
| | |
| } 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"); |
| |
| } | } |
| } | } |
| | |
| | // 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(); |
| |
| 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); |
| |
| 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)) { |