Repository navigation
Allow reusing computation and data between CVC objects - #644
giacomofiorin wants to merge 18 commits into
Conversation
14dbb58 to
d1333a7
Compare
1b984a4 to
5ec6f56
Compare
|
Hi @giacomofiorin! It seems there is no change to the SMP related code, but maybe I miss somethign. If a cvc reuses other ones, how does the SMP work? |
See the last paragraph description. What I would like to do is keep a list of reusable CVCs, and run those first over SMP. |
Thanks! I am working on |
That's one option, but then all CVC classes would then need to have the necessary API changes. Okay with you if we discuss more next week? |
OK. |
|
Also, So it would make sense to make them special cases. |
|
I just doubt if it is really necessary to distribute the CVCs to different threads. On one hand, if the thread for CVC_A takes longer time to finish, then the thread for CVC_B has to wait. On the other hand, it makes the reusable CVCs complicated. I still think it would be better if Colvars can give up distributing CVCs among threads completely, and use the threads in fine-grained cases, such as calculating correlation and overlapping matrices in optimal rotation, rotating atoms, and projecting a hill in metadynamics... |
|
The main use case for distributing CVCs over threads is that of multiple, expensive CVCs of similar computational cost. It is hard to know how widespread it is, but in that case, the current implementation has the benefit of improving performance transparently without user intervention. |
If the expensive CVCs refer to RMSD, anything related to optimal alignment and coordNum, then it is still better to distribute the expensive loops in a single CVC over threads instead of distribute CVCs themselves over threads. The former approach is more friendly to CPU caches, especially the shared L3 cache, because the CPU could fetch a large chunk of continuous data, instead of many small blocks, from the main memory to L3 directly for all threads. Moreover, distributing the inner loops in a CVC over threads can guarantee to accelerate the calculation of the most expensive CVC. Distributing CVCs over threads means that the slowest CVC uses always a single thread, and if there are other fast CVs (for example, using RMSD, gyration and a few distances), then the other threads have to wait for the slowest done. |
|
I am looking into ways to parallelize the computation of CVCs even if there are interdependencies. With the dependencies the CVCs are in a compute graph. A way to parallelize the computation could be:
For example, in the following compute graph, Now the major issue is that Colvars actually does not distribute CVCs over threads, and instead, the parallelization scheme is worse than I expect. It distributes the outermost |
|
There should be an additional reusability improvement about applying biasing force only to a single value of a vector CV, so that we can safely compute some CVs in a single CVC and output a vector, and then bias only a scalar part of the vector. |
Hi @HanatoK the loop is done over colvar objects to avoid having the top-level class reach too deep into the CVCs. However, in the case of a colvar with multiple CVCs, that colvar will appear in the loop multiple times to achieve parallelism over the entire collection of CVC objects. Lines 273 to 277 in 0f2d682 Your other comments all make sense, especially regarding building a dependency tree and do multiple parallel loops, one for each level. This PR is meant to introduce only one level to start with. |
d05529f to
32a3aca
Compare
There was a problem hiding this comment.
🟡 Changes recommended
Reuse validation is bypassed, force signs can be incorrect, and SMP execution has an unresolved data race.
Once you've addressed the issues Copilot identified, you can request another Copilot review.
Pull request overview
Adds reusable orientation computation across CVC objects, implementing the first phase of issue #232.
Changes:
- Adds named CVC registration and reuse dependencies.
- Refactors quaternion rotation calculations and force propagation.
- Adds spin-angle reuse regression fixtures.
File summaries
| File | Description |
|---|---|
tests/spinangle-reuse.in |
Adds reuse configuration test. |
tests/input_files/spinangle-reuse_harmonic-fixed/test.in |
Adds harmonic reuse test. |
tests/build_tests.sh |
Registers the new test. |
src/colvartypes.h |
Moves rotation operations onto quaternions. |
src/colvartypes.cpp |
Implements quaternion rotation operations. |
src/colvarmodule.h |
Adds CVC registry and smart-pointer ownership. |
src/colvarmodule.cpp |
Implements registry, lifecycle, and output updates. |
src/colvardeps.h |
Defines CVC reuse features. |
src/colvarcomp.h |
Adds reusable-CVC state and orientation helpers. |
src/colvarcomp.cpp |
Implements registration and reuse validation. |
src/colvarcomp_rotations.cpp |
Reuses orientation values and forwards forces. |
src/colvarbias.cpp |
Adjusts trajectory spacing. |
src/colvar.h |
Exposes reference cleanup. |
src/colvar.cpp |
Registers and dereferences named components. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/test.in |
Adds NAMD reuse configuration. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/notes.txt |
Documents fixture provenance. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/namd-version.txt |
Records fixture environment. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/AutoDiff/test.restart.colvars.traj |
Adds restart trajectory reference. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/AutoDiff/test.restart.colvars.state.stripped |
Adds restart state reference. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/AutoDiff/test.restart.colvars.out |
Adds restart output reference. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/AutoDiff/test.colvars.traj |
Adds trajectory reference. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/AutoDiff/test.colvars.state.stripped |
Adds state reference. |
namd/tests/library/000_spinangle-reuse_harmonic-fixed/AutoDiff/test.colvars.out |
Adds output reference. |
Review details
Suppressed comments (1)
src/colvarcomp_rotations.cpp:285
dx_dqis the derivative with respect to the raw rotationq, but the source orientation'sapply_force()interprets its argument as a force on its sign-adjusted output and multiplies by theclosestToQuaternionsign. When the source chooses-q, this forwarding reverses every scalar consumer's atomic force, making an output-only option affect dynamics. Forward through an API that applies a derivative with respect to rawq, or compensate for the source sign before callingapply_force().
if (has_precomputed_cvc("orientation")) {
get_orientation_cvc()->apply_force(force.real_value * dx_dq);
- Files reviewed: 23/23 changed files
- Comments generated: 4
- Review effort level: Balanced
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
3b8e7d0 to
012f044
Compare
|
In one of the recent commits I work around the race condition in SMP parallelism (also correctly identified by the AI). Putting it up for review. Note: I have not added documentation changes, since I'd first like feedback on the |
6323cf4 to
147af9e
Compare
There was a problem hiding this comment.
Copilot review overview
🟡 Changes recommended
Reused rotations currently have initialization, scheduling, and force-sign correctness defects.
Get a fresh assessment by requesting another Copilot review.
Review effort: Balanced
Findings: 4
Open (6)
Regenerate golden output with Colvars SMP disabled · New Regenerate restart output with Colvars SMP disabled · New Avoid constructing derivative with uninitialized position counts · New This reference does not establish any scheduling dependency on the source colvar. If the source has… Assign component names before setup diagnostics · New Reject only unsupported SMP scheduling modes · New
| colvars: Protein_C-alpha_10 (1 atoms) | ||
| colvars: group10 (4 atoms) | ||
| colvars: heavy_atoms (51 atoms) | ||
| colvars: # smp = on [default] |
| colvars: Protein_C-alpha_10 (1 atoms) | ||
| colvars: group10 (4 atoms) | ||
| colvars: heavy_atoms (51 atoms) | ||
| colvars: # smp = on [default] |
| colvar::orientation::orientation() | ||
| { | ||
| set_function_type("orientation"); | ||
| rot_deriv_impl = std::unique_ptr<rotation_derivative_impl_>(new rotation_derivative_impl_(this)); |
| // Update properties and check that dynamic dependencies are fulfilled | ||
| error_code_this |= cvcs.back()->setup(); |
| if (cvmodule->proxy->get_smp_mode() != colvarproxy_smp::smp_mode_t::none) { | ||
| error_code |= cvmodule->error( | ||
| "reusable components are currently not compatible with SMP object-based parallelism" | ||
| "(\"smp on\" or \"smp cvcs\"); please disable it by using \"smp off\" as a global " | ||
| "keyword", | ||
| COLVARS_INPUT_ERROR); |
…ber (the latter may be unused)
Co-authored-by: Copilot Autofix powered by AI <175728472+Copilot@users.noreply.github.com>
147af9e to
8e94740
Compare
| inline cvm::real spin_angle(cvm::rvector const &axis) const | ||
| { | ||
| cvm::rvector const q_vec = q.get_vector(); | ||
| cvm::real alpha = (180.0/PI) * 2.0 * cvm::atan2(axis * q_vec, q.q0); | ||
| while (alpha > 180.0) alpha -= 360; | ||
| while (alpha < -180.0) alpha += 360; | ||
| return alpha; | ||
| } | ||
|
|
||
| /// \brief Return the derivative of the spin angle with respect to | ||
| /// the quaternion | ||
| inline cvm::quaternion dspin_angle_dq(cvm::rvector const &axis) const | ||
| { | ||
| cvm::rvector const q_vec = q.get_vector(); | ||
| cvm::real const iprod = axis * q_vec; | ||
|
|
||
| if (q.q0 != 0.0) { | ||
|
|
||
| cvm::real const dspindx = | ||
| (180.0/PI) * 2.0 * (1.0 / (1.0 + (iprod*iprod)/(q.q0*q.q0))); | ||
|
|
||
| return cvm::quaternion( dspindx * (iprod * (-1.0) / (q.q0*q.q0)), | ||
| dspindx * ((1.0/q.q0) * axis.x), | ||
| dspindx * ((1.0/q.q0) * axis.y), | ||
| dspindx * ((1.0/q.q0) * axis.z)); | ||
| } else { | ||
| // (1/(1+x^2)) ~ (1/x)^2 | ||
| // The documentation of spinAngle discourages its use when q_vec and | ||
| // axis are not close | ||
| return cvm::quaternion((180.0/PI) * 2.0 * ((-1.0)/iprod), 0.0, 0.0, 0.0); | ||
| } | ||
| } | ||
|
|
||
| /// \brief Return the projection of the orientation vector onto a | ||
| /// predefined axis | ||
| inline cvm::real cos_theta(cvm::rvector const &axis) const | ||
| { | ||
| cvm::rvector const q_vec = q.get_vector(); | ||
| cvm::real const alpha = | ||
| (180.0/PI) * 2.0 * cvm::atan2(axis * q_vec, q.q0); | ||
|
|
||
| cvm::real const cos_spin_2 = cvm::cos(alpha * (PI/180.0) * 0.5); | ||
| cvm::real const cos_theta_2 = ( (cos_spin_2 != 0.0) ? | ||
| (q.q0 / cos_spin_2) : | ||
| (0.0) ); | ||
| // cos(2t) = 2*cos(t)^2 - 1 | ||
| return 2.0 * (cos_theta_2*cos_theta_2) - 1.0; | ||
| } | ||
|
|
||
| /// Return the derivative of the tilt wrt the quaternion | ||
| inline cvm::quaternion dcos_theta_dq(cvm::rvector const &axis) const | ||
| { | ||
| cvm::rvector const q_vec = q.get_vector(); | ||
| cvm::real const iprod = axis * q_vec; | ||
|
|
||
| cvm::real const cos_spin_2 = cvm::cos(cvm::atan2(iprod, q.q0)); | ||
|
|
||
| if (q.q0 != 0.0) { | ||
|
|
||
| cvm::real const d_cos_theta_dq0 = | ||
| (4.0 * q.q0 / (cos_spin_2*cos_spin_2)) * | ||
| (1.0 - (iprod*iprod)/(q.q0*q.q0) / (1.0 + (iprod*iprod)/(q.q0*q.q0))); | ||
|
|
||
| cvm::real const d_cos_theta_dqn = | ||
| (4.0 * q.q0 / (cos_spin_2*cos_spin_2) * | ||
| (iprod/q.q0) / (1.0 + (iprod*iprod)/(q.q0*q.q0))); | ||
|
|
||
| return cvm::quaternion(d_cos_theta_dq0, | ||
| d_cos_theta_dqn * axis.x, | ||
| d_cos_theta_dqn * axis.y, | ||
| d_cos_theta_dqn * axis.z); | ||
| } else { | ||
|
|
||
| cvm::real const d_cos_theta_dqn = | ||
| (4.0 / (cos_spin_2*cos_spin_2 * iprod)); | ||
|
|
||
| return cvm::quaternion(0.0, | ||
| d_cos_theta_dqn * axis.x, | ||
| d_cos_theta_dqn * axis.y, | ||
| d_cos_theta_dqn * axis.z); | ||
| } | ||
| } |
There was a problem hiding this comment.
Why are these functions moved to the .cpp file? I expect them to be inlined and can be used in future GPU porting.
There was a problem hiding this comment.
You're right: they probably belong in a *_impl.h header though.
| cvcp->name = s.str(); | ||
| } | ||
| // Record the name of the parent colvar to make it globally searchable | ||
| cvcp->set_parent_name(name); |
There was a problem hiding this comment.
A single parent name seems strange. What would happen if the CVC has multiple parents? Are we going to reuse the whole CVCs, or only reuse some specific computational details (like q in some of them) of a few specific CVCs?
| if (cvc_ptr) { | ||
| precomputed_cvcs[id] = cvc_ptr; | ||
| disable(f_cvc_explicit_atom_groups); | ||
| disable(f_cvc_explicit_gradient); |
There was a problem hiding this comment.
If so, then is debugGradients on still possible?
| std::vector<std::string> function_types; | ||
|
|
||
| /// Role identifiers and pointers to the external CVCs whose data is being reused | ||
| std::map<std::string, std::shared_ptr<cvc>> precomputed_cvcs; |
There was a problem hiding this comment.
std::unordered_map might have better performance.
|
@HanatoK at the moment I'm simply trying to keep this branch up to date. I'm going to test out the "Watch and Fix" Copilot feature, so that it keeps an eye out for merge conflicts as more PRs are merged into master. |



This PR contains an implementation of the code to reuse the output of certain computations between colvar component (CVC) objects. The scheme aims to provide support a broader set of use cases than the ones described in #232, i.e. it aims to share specific pieces of computation between CVCs of different types, as well as allowing concurrent application of forces.
The proposed scheme works as follows. In the configuration string of CVC_B, a reference is provided to CVC_A, which is an object been defined previously in the Colvars config. The type of CVC_A may be the same as CVC_B (in which case its value is simply copied), or it may be a base class of CVC_B, or the same as one of its members. In the latter two cases, It's assumed that the computation of CVC_B can be written as follows:
When CVC_A is an object that was already computed, CVC_B can skip the first two steps altogether, and just replace the third step with an assignment.
This PR adds the ability to reuse an
orientationcomponent in atilt,spinAngle,euler{Phi,Theta,Psi}, component. Other potential use cases, not covered in this PR, are:distanceVeccomponent as input todistance,distanceZ,distanceXY, etc.gyration,inertia,inertiaZ, etcIn all of the above, the speedup would be either 2 or 3, assuming that the explicit loop over atoms is speed-limiting. (The inertia tensor has six independent components, but I assume that very few would be interested in biasing them all!)
Because the SMP loop over objects is currently unaware of the mutual dependencies, the
reuseComponentkeyword is not be allowed in the default case (smp onorsmp cvcs). Using it would require thatsmp off,smp inner_looporsmp gpuare used instead. Future improvements may lift this constraint, e.g. by leveraging the fact that currently only one level of reusal is allowed.Implements #232