Expose computed density in InterRDF results

[sumitsaraswat362]

Fixes #5362

Changes made in this Pull Request:

• Expose computed average number density as results.density (and alias results.rho) in InterRDF and InterRDF_s when norm is "rdf" or "density".
• Guard against zero-division / zero box volume cases when computing density.
• Add unit tests for exposed density in test_rdf.py and test_rdf_s.py.
• Update package/CHANGELOG and package/AUTHORS.

LLM / AI generated code disclosure

LLMs or other AI-powered tools (beyond simple IDE use cases) were used in this contribution: no

PR Checklist

• Issue raised/referenced?
• Tests updated/added?
• Documentation updated/added?
• package/CHANGELOG file updated?
• Is your name in package/AUTHORS? (If it is not, add it!)
• I have read and understand the current AI Policy
• LLM/AI disclosure was updated.

Developers Certificate of Origin

I certify that I can submit this code contribution as described in the Developer Certificate of Origin, under the MDAnalysis LICENSE.

[read-the-docs-community]

Documentation build overview

📚 MDAnalysis | 🛠️ Build #35213681 | 📁 Comparing 47aac89 against latest (4f6d54a)

  🔍 Preview build  

10 files changed · ± 10 modified

± Modified

• genindex.html
• index.html
• documentation_pages/analysis/rdf.html
• documentation_pages/analysis/wbridge_analysis.html
• documentation_pages/coordinates/TPR.html
• _modules/MDAnalysis/analysis/rdf.html
• _modules/MDAnalysis/coordinates/MMCIF.html
• _modules/MDAnalysis/coordinates/TPR.html
• _modules/MDAnalysis/core/topologyattrs.html
• _modules/MDAnalysis/topology/MMCIFParser.html

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 93.89%. Comparing base (e4884fa) to head (47aac89).
⚠️ Report is 1 commits behind head on develop.

Additional details and impacted files

@@           Coverage Diff            @@
##           develop    #5482   +/-   ##
========================================
  Coverage    93.88%   93.89%           
========================================
  Files          184      184           
  Lines        22640    22657   +17     
  Branches      3217     3221    +4     
========================================
+ Hits         21256    21273   +17     
  Misses         920      920           
  Partials       464      464           

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[sici17]

I think that trying to compute the average density value in the norm="density" case might lead to zero division erros given that a nonzero box_vol is not actually required in these cases in order to compute the distance dependent density. Probably the two cases should be handled separately. I am talking about the line 401 of the diff. I am not an expert or anything though so I would wait someone with more expertise to express his opinion obviously. I was just practicing code reviewing and understing.

[sumitsaraswat362]

Thanks for the great catch, @sici17! You are completely right—for non-periodic trajectories or simulations without box dimensions where norm="density" is used, box_vol can be zero, which would have triggered an unintended ZeroDivisionError.

I've just pushed a commit guarding against this:

1. results.density (and results.rho) is now only computed when box_vol > 0, and set to None if no box volume is available.
2. In the norm="rdf" path, it proceeds with normalization as before.
3. Also exposed results.rho as an alias alongside results.density as originally suggested in Calculating the rdf for norm = "rdf" or "density" at the same time. #5362.

Really appreciate you taking the time to review!

[orbeckst]

Before we can review this PR, please restore the original PR template https://github.com/MDAnalysis/mdanalysis/blob/develop/.github/PULL_REQUEST_TEMPLATE.md and answer all questions.

[sumitsaraswat362]

@orbeckst I have restored the official PR template in the PR description and answered all questions, updated package/CHANGELOG, added my name to package/AUTHORS, and added unit tests in test_rdf.py and test_rdf_s.py.

@sici17 Good catch! I added guards for zero box volume / non-periodic systems (box_vol > 0 and volume_cum > 0) so that results.density and results.rho are cleanly set to None without raising any ZeroDivisionError.

[orbeckst]

Quick initial review

• only include results.density (drop results.rho) — there should be only one way to access the data
• add documentation, for instance, results.density needs to be explained where the other results attributes are listed; add a .. versionadded:: 2.11.0

These are the first things I noticed, there may be more comments.

[orbeckst]

Please also resolve the conflicts with the develop branch.

[sumitsaraswat362]

@orbeckst Thanks for the review!

I have addressed the points:

• Dropped results.rho completely so results.density is the sole way to access the data.
• Added documentation for results.density under the Attributes section of both InterRDF and InterRDF_s with .. versionadded:: 2.11.0.
• Merged the latest develop branch and resolved the merge conflict in package/CHANGELOG.
• Updated test assertions in test_rdf.py and test_rdf_s.py.

[orbeckst]

The computed "density" is $N_a N_b / \langle V \rangle$. This is not a real density and it's not the $\rho$ from the docs at the top. One could write it as $N_a \rho_b = N_b \rho_a$. We cannot store this as an attribute results.density because it's not a density.

(Furthermore, it's confusing because norm="rdf gives results.rdf and then people/LLMs may look for results.density for norm="density".)

We could store it as results.normalization_density or something like that. What would be a good name? Does it still make sense to store this quantity or should we be storing $\rho_a$, $\rho_b$, $\langle V \rangle$ as well?

Whatever we store, we need then document properly, i.e., write down the equations that define the quantity.

We cannot really continue with the PR until we have decided this.

[sumitsaraswat362]

Thank you for pointing this out, Oliver. You are completely right. The quantity $N_a N_b / \langle V \rangle$ is the pair normalization density $N_a \rho_b$ rather than the single-species density $\rho$, and exposing it as results.density is inaccurate and misleading for coordination calculations ($N_{ab} = \rho G_{ab}$). It also collides with norm="density" where users expect the radial profile.

I have updated the implementation to address this:

1. results.normalization_density: Renamed from results.density, storing the pair normalization density $N_{\text{pairs}} / \langle V \rangle = (N_a N_b - N_{\text{excl}}) / \langle V \rangle$ used in the RDF denominator (or $1 / \langle V \rangle$ for InterRDF_s).
2. results.box_volume: Stores the trajectory-averaged box volume $\langle V \rangle$ across analyzed frames.
3. results.density_a and results.density_b: Store the true bulk number densities $\rho_a = N_a / \langle V \rangle$ and $\rho_b = N_b / \langle V \rangle$, where $\rho_b$ is the density $\rho$ appearing in the documentation for $n_{ab}(r) = \rho g_{ab}(r)$ and $N_{ab}(r) = \rho G_{ab}(r)$.
4. If the box volume is zero or unavailable, all four attributes evaluate to None.
5. When norm="none", none of these normalization attributes are set.

The equations and physical definitions have been added to the docstrings of both InterRDF and InterRDF_s with .. versionadded:: 2.11.0, and the test suites in test_rdf.py and test_rdf_s.py have been updated with numerical assertions against the theoretical values.