diff --git a/lib/colvars/colvaratoms.cpp b/lib/colvars/colvaratoms.cpp index cd47af9512..99f63c5af6 100644 --- a/lib/colvars/colvaratoms.cpp +++ b/lib/colvars/colvaratoms.cpp @@ -315,14 +315,15 @@ int cvm::atom_group::parse(std::string const &conf, #endif if (!b_dummy) { + + // calculate total mass (TODO: this is the step that most needs deferred re-initialization) this->total_mass = 0.0; for (cvm::atom_iter ai = this->begin(); ai != this->end(); ai++) { this->total_mass += ai->mass; } - } - if (!b_dummy) { + // whether these atoms will ever receive forces or not bool enable_forces = true; // disableForces is deprecated if (get_keyval(group_conf, "enableForces", enable_forces, true, mode)) { @@ -330,6 +331,8 @@ int cvm::atom_group::parse(std::string const &conf, } else { get_keyval(group_conf, "disableForces", noforce, false, colvarparse::parse_silent); } + + get_keyval(group_conf, "weights", weights, weights, colvarparse::parse_silent); } // FITTING OPTIONS diff --git a/lib/colvars/colvaratoms.h b/lib/colvars/colvaratoms.h index 5a44d106cb..d950ffb7ac 100644 --- a/lib/colvars/colvaratoms.h +++ b/lib/colvars/colvaratoms.h @@ -139,6 +139,9 @@ public: /// Allocates and populates the sorted list of atom ids int create_sorted_ids(void); + /// List of user-defined weights to be used by certain CVCs + std::vector weights; + /// \brief When updating atomic coordinates, translate them to align with the /// center of mass of the reference coordinates diff --git a/lib/colvars/colvarbias_abf.cpp b/lib/colvars/colvarbias_abf.cpp index d2e35c85fe..2685ea0e9d 100644 --- a/lib/colvars/colvarbias_abf.cpp +++ b/lib/colvars/colvarbias_abf.cpp @@ -118,7 +118,7 @@ colvarbias_abf::colvarbias_abf(std::string const &conf, char const *key) force = new cvm::real [colvars.size()]; // Construct empty grids based on the colvars - samples = new colvar_grid_count (colvars); + samples = new colvar_grid_count(colvars); gradients = new colvar_grid_gradient(colvars); gradients->samples = samples; samples->has_parent_data = true; @@ -126,7 +126,7 @@ colvarbias_abf::colvarbias_abf(std::string const &conf, char const *key) // For shared ABF, we store a second set of grids. // This used to be only if "shared" was defined, // but now we allow calling share externally (e.g. from Tcl). - last_samples = new colvar_grid_count (colvars); + last_samples = new colvar_grid_count(colvars); last_gradients = new colvar_grid_gradient(colvars); last_gradients->samples = last_samples; last_samples->has_parent_data = true; @@ -579,7 +579,7 @@ colvarbias_histogram::colvarbias_histogram(std::string const &conf, char const * cvm::error("User required histogram with zero output frequency"); } - grid = new colvar_grid_count (colvars); + grid = new colvar_grid_count(colvars); bin.assign(colvars.size(), 0); cvm::log("Finished histogram setup.\n"); @@ -588,7 +588,8 @@ colvarbias_histogram::colvarbias_histogram(std::string const &conf, char const * /// Destructor colvarbias_histogram::~colvarbias_histogram() { - if (grid_os.is_open()) grid_os.close(); + if (grid_os.is_open()) + grid_os.close(); if (grid) { delete grid; diff --git a/lib/colvars/colvarbias_meta.cpp b/lib/colvars/colvarbias_meta.cpp index c329c04d5a..812924beed 100644 --- a/lib/colvars/colvarbias_meta.cpp +++ b/lib/colvars/colvarbias_meta.cpp @@ -86,7 +86,7 @@ colvarbias_meta::colvarbias_meta(std::string const &conf, char const *key) colvars[i]->enable(colvar::task_grid); } - hills_energy = new colvar_grid_scalar (colvars); + hills_energy = new colvar_grid_scalar(colvars); hills_energy_gradients = new colvar_grid_gradient(colvars); } else { rebin_grids = false; @@ -959,7 +959,7 @@ void colvarbias_meta::update_replicas_registry() (replicas.back())->comm = multiple_replicas; if (use_grids) { - (replicas.back())->hills_energy = new colvar_grid_scalar (colvars); + (replicas.back())->hills_energy = new colvar_grid_scalar(colvars); (replicas.back())->hills_energy_gradients = new colvar_grid_gradient(colvars); } } diff --git a/lib/colvars/colvarcomp.h b/lib/colvars/colvarcomp.h index fb104aa28f..62bd8a274e 100644 --- a/lib/colvars/colvarcomp.h +++ b/lib/colvars/colvarcomp.h @@ -1139,6 +1139,8 @@ class colvar::cartesian protected: /// Atom group cvm::atom_group atoms; + /// Which Cartesian coordinates to include + std::vector axes; public: cartesian(std::string const &conf); cartesian(); diff --git a/lib/colvars/colvarcomp_distances.cpp b/lib/colvars/colvarcomp_distances.cpp index f67ac87baa..a5d921b34f 100644 --- a/lib/colvars/colvarcomp_distances.cpp +++ b/lib/colvars/colvarcomp_distances.cpp @@ -1296,17 +1296,40 @@ colvar::cartesian::cartesian(std::string const &conf) parse_group(conf, "atoms", atoms); atom_groups.push_back(&atoms); + bool use_x, use_y, use_z; + get_keyval(conf, "useX", use_x, true); + get_keyval(conf, "useY", use_y, true); + get_keyval(conf, "useZ", use_z, true); + + axes.clear(); + if (use_x) axes.push_back(0); + if (use_y) axes.push_back(1); + if (use_z) axes.push_back(2); + + if (axes.size() == 0) { + cvm::error("Error: a \"cartesian\" component was defined with all three axes disabled.\n"); + } + x.type(colvarvalue::type_vector); - x.vector1d_value.resize(atoms.size() * 3); + x.vector1d_value.resize(atoms.size() * axes.size()); } void colvar::cartesian::calc_value() { - int ia, j; + size_t const dim = axes.size(); + size_t ia, j; for (ia = 0; ia < atoms.size(); ia++) { - for (j = 0; j < 3; j++) { - x.vector1d_value[3*ia + j] = atoms[ia].pos[j]; + for (j = 0; j < dim; j++) { + x.vector1d_value[dim*ia + j] = atoms[ia].pos[axes[j]]; + } + } + + if (atoms.weights.size()) { + for (ia = 0; ia < atoms.size(); ia++) { + for (j = 0; j < dim; j++) { + x.vector1d_value[dim*ia + j] *= atoms.weights[ia]; + } } } } @@ -1322,14 +1345,24 @@ void colvar::cartesian::calc_gradients() void colvar::cartesian::apply_force(colvarvalue const &force) { - int ia, j; + size_t const dim = axes.size(); + size_t ia, j; if (!atoms.noforce) { cvm::rvector f; - for (ia = 0; ia < atoms.size(); ia++) { - for (j = 0; j < 3; j++) { - f[j] = force.vector1d_value[3*ia + j]; + if (atoms.weights.size()) { + for (ia = 0; ia < atoms.size(); ia++) { + for (j = 0; j < dim; j++) { + f[axes[j]] = force.vector1d_value[dim*ia + j]; + } + atoms[ia].apply_force(f); + } + } else { + for (ia = 0; ia < atoms.size(); ia++) { + for (j = 0; j < dim; j++) { + f[axes[j]] = force.vector1d_value[dim*ia + j] / atoms.weights[ia]; + } + atoms[ia].apply_force(f); } - atoms[ia].apply_force(f); } } } diff --git a/lib/colvars/colvarparse.cpp b/lib/colvars/colvarparse.cpp index a0f14c3877..84777149c0 100644 --- a/lib/colvars/colvarparse.cpp +++ b/lib/colvars/colvarparse.cpp @@ -482,16 +482,6 @@ bool colvarparse::key_lookup(std::string const &conf, } } - // check it is not between quotes - // if ( (conf.find_last_of ("\"", - // conf.find_last_of (white_space, pos)) != - // std::string::npos) && - // (conf.find_first_of ("\"", - // conf.find_first_of (white_space, pos)) != - // std::string::npos) ) - // return false; - - // save the pointer for a future call (when iterating over multiple // valid instances of the same keyword) save_pos = pos + key.size(); @@ -499,7 +489,7 @@ bool colvarparse::key_lookup(std::string const &conf, // get the remainder of the line size_t pl = conf.rfind("\n", pos); size_t line_begin = (pl == std::string::npos) ? 0 : pos; - size_t nl = conf.find ("\n", pos); + size_t nl = conf.find("\n", pos); size_t line_end = (nl == std::string::npos) ? conf.size() : nl; std::string line(conf, line_begin, (line_end-line_begin)); @@ -569,7 +559,7 @@ bool colvarparse::key_lookup(std::string const &conf, if (data.size() && save_delimiters) { data_begin_pos.push_back(conf.find(data, pos+key.size())); - data_end_pos.push_back (data_begin_pos.back()+data.size()); + data_end_pos.push_back(data_begin_pos.back()+data.size()); // std::cerr << "key = " << key << ", data = \"" // << data << "\", data_begin, data_end = " // << data_begin_pos.back() << ", " << data_end_pos.back()