Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improve detection of hard/mathematical boundaries #731

Merged
merged 3 commits into from
Nov 13, 2024
Merged

Improve detection of hard/mathematical boundaries #731

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
b741a6d
Break up long function
giacomofiorin Oct 16, 2024
a3f2328
Use constexpr to define tolerance parameter in a single place
giacomofiorin Oct 16, 2024
0c60c20
Allow detecting when user-defined boundaries match mathematical bound…
giacomofiorin Oct 16, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
70 changes: 50 additions & 20 deletions src/colvar.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -501,8 +501,6 @@ int colvar::init_grid_parameters(std::string const &conf)
{
int error_code = COLVARS_OK;

colvarmodule *cv = cvm::main();

cvm::real default_width = width;

if (!key_already_set("width")) {
Expand All @@ -528,34 +526,68 @@ int colvar::init_grid_parameters(std::string const &conf)

if (is_enabled(f_cv_scalar)) {

if (is_enabled(f_cv_single_cvc)) {
// Get the default boundaries from the component
// Record the CVC's intrinsic boundaries, and set them as default values for the user's choice
colvarvalue cvc_lower_boundary, cvc_upper_boundary;

if (is_enabled(f_cv_single_cvc)) { // Get the intrinsic boundaries of the CVC

if (cvcs[0]->is_enabled(f_cvc_lower_boundary)) {
enable(f_cv_lower_boundary);
enable(f_cv_hard_lower_boundary);
lower_boundary =
lower_boundary = cvc_lower_boundary =
*(reinterpret_cast<colvarvalue const *>(cvcs[0]->get_param_ptr("lowerBoundary")));
}

if (cvcs[0]->is_enabled(f_cvc_upper_boundary)) {
enable(f_cv_upper_boundary);
enable(f_cv_hard_upper_boundary);
upper_boundary =
*(reinterpret_cast<colvarvalue const *>(cvcs[0]->get_param_ptr("upperBoundary")));
upper_boundary = cvc_upper_boundary =
*(reinterpret_cast<colvarvalue const *>(cvcs[0]->get_param_ptr("upperBoundary")));
}
}

if (get_keyval(conf, "lowerBoundary", lower_boundary, lower_boundary)) {
enable(f_cv_lower_boundary);
// Because this is the user's choice, we cannot assume it is a true
// physical boundary
disable(f_cv_hard_lower_boundary);
if (is_enabled(f_cv_single_cvc) && is_enabled(f_cv_hard_lower_boundary)) {
if (cvm::sqrt(dist2(lower_boundary, cvc_lower_boundary))/width > colvar_boundaries_tol) {
// The user choice is different from the CVC's default
disable(f_cv_hard_lower_boundary);
}
}
}

if (get_keyval(conf, "upperBoundary", upper_boundary, upper_boundary)) {
enable(f_cv_upper_boundary);
disable(f_cv_hard_upper_boundary);
if (is_enabled(f_cv_single_cvc) && is_enabled(f_cv_hard_upper_boundary)) {
if (cvm::sqrt(dist2(upper_boundary, cvc_upper_boundary))/width > colvar_boundaries_tol) {
disable(f_cv_hard_upper_boundary);
}
}
}

get_keyval_feature(this, conf, "hardLowerBoundary", f_cv_hard_lower_boundary,
is_enabled(f_cv_hard_lower_boundary));

get_keyval_feature(this, conf, "hardUpperBoundary", f_cv_hard_upper_boundary,
is_enabled(f_cv_hard_upper_boundary));

get_keyval(conf, "expandBoundaries", expand_boundaries, expand_boundaries);

error_code |= parse_legacy_wall_params(conf);
error_code |= check_grid_parameters();
}

return error_code;
}


int colvar::parse_legacy_wall_params(std::string const &conf)
{
int error_code = COLVARS_OK;
colvarmodule *cv = cvm::main();

if (is_enabled(f_cv_scalar)) {

// Parse legacy wall options and set up a harmonicWalls bias if needed
cvm::real lower_wall_k = 0.0, upper_wall_k = 0.0;
cvm::real lower_wall = 0.0, upper_wall = 0.0;
Expand Down Expand Up @@ -609,13 +641,14 @@ harmonicWalls {\n\
}
}

get_keyval_feature(this, conf, "hardLowerBoundary", f_cv_hard_lower_boundary,
is_enabled(f_cv_hard_lower_boundary));
return error_code;
}

get_keyval_feature(this, conf, "hardUpperBoundary", f_cv_hard_upper_boundary,
is_enabled(f_cv_hard_upper_boundary));

// consistency checks for boundaries and walls
int colvar::check_grid_parameters()
{
int error_code = COLVARS_OK;

if (is_enabled(f_cv_lower_boundary) && is_enabled(f_cv_upper_boundary)) {
if (lower_boundary >= upper_boundary) {
error_code |= cvm::error("Error: the upper boundary, "+
Expand All @@ -626,7 +659,6 @@ harmonicWalls {\n\
}
}

get_keyval(conf, "expandBoundaries", expand_boundaries, expand_boundaries);
if (expand_boundaries && periodic_boundaries()) {
error_code |= cvm::error("Error: trying to expand boundaries that already "
"cover a whole period of a periodic colvar.\n",
Expand Down Expand Up @@ -2190,12 +2222,10 @@ int colvar::set_cvc_param(std::string const &param_name, void const *new_value)
bool colvar::periodic_boundaries(colvarvalue const &lb, colvarvalue const &ub) const
{
if (period > 0.0) {
if ( ((cvm::sqrt(this->dist2(lb, ub))) / this->width)
< 1.0E-10 ) {
if (((cvm::sqrt(this->dist2(lb, ub))) / this->width) < colvar_boundaries_tol) {
return true;
}
}

return false;
}

Expand Down
12 changes: 12 additions & 0 deletions src/colvar.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -263,6 +263,12 @@ class colvar : public colvarparse, public colvardeps {
/// Init defaults for grid options
int init_grid_parameters(std::string const &conf);

/// Consistency check for the grid paramaters
int check_grid_parameters();

/// Read legacy wall keyword (these are biases now)
int parse_legacy_wall_params(std::string const &conf);

/// Init extended Lagrangian parameters
int init_extended_Lagrangian(std::string const &conf);

Expand Down Expand Up @@ -778,4 +784,10 @@ inline void colvar::reset_bias_force() {
fb_actual.reset();
}


namespace {
// Tolerance parameter to decide when two boundaries coincide
constexpr cvm::real colvar_boundaries_tol = 1.0e-10;
}

#endif
Loading