-
Notifications
You must be signed in to change notification settings - Fork 19
refactor: Separate histogram calculation from PartialSet #2216
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.
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
Merged
Merged
Changes from all commits
Commits
Show all changes
4 commits
Select commit
Hold shift + click to select a range
File filter
Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
There are no files selected for viewing
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,128 @@ | ||
| // SPDX-License-Identifier: GPL-3.0-or-later | ||
| // Copyright (c) 2025 Team Dissolve and contributors | ||
|
|
||
| #include "classes/histogramSet.h" | ||
| #include "base/lineParser.h" | ||
| #include "classes/atomType.h" | ||
| #include "items/deserialisers.h" | ||
| #include "math/mathFunc.h" | ||
| #include "templates/algorithms.h" | ||
|
|
||
| HistogramSet::~HistogramSet() | ||
| { | ||
| fullHistograms_.clear(); | ||
| boundHistograms_.clear(); | ||
| unboundHistograms_.clear(); | ||
| } | ||
|
|
||
| /* | ||
| * Data | ||
| */ | ||
|
|
||
| // Initialise histograms | ||
| void HistogramSet::initialise(const AtomTypeMix &atomTypeMix, double rdfRange, double binWidth) | ||
| { | ||
| atomTypeMix_ = atomTypeMix; | ||
|
|
||
| auto nTypes = atomTypeMix_.nItems(); | ||
|
|
||
| fullHistograms_.initialise(nTypes, nTypes, half_); | ||
| boundHistograms_.initialise(nTypes, nTypes, half_); | ||
| unboundHistograms_.initialise(nTypes, nTypes, half_); | ||
|
|
||
| dissolve::for_each_pair( | ||
| ParallelPolicies::par, nTypes, | ||
| [&](int i, int j) | ||
| { | ||
| fullHistograms_[{i, j}].initialise(0.0, rdfRange, binWidth); | ||
| boundHistograms_[{i, j}].initialise(0.0, rdfRange, binWidth); | ||
| unboundHistograms_[{i, j}].initialise(0.0, rdfRange, binWidth); | ||
| }, | ||
| half_); | ||
| } | ||
|
|
||
| // Clear all histogram data | ||
| void HistogramSet::clear() | ||
| { | ||
| fullHistograms_.clear(); | ||
| boundHistograms_.clear(); | ||
| unboundHistograms_.clear(); | ||
| } | ||
|
|
||
| // Zero histogram bins | ||
| void HistogramSet::zeroBins() | ||
| { | ||
| for (auto &histo : fullHistograms_.linearArray()) | ||
| histo.zeroBins(); | ||
| for (auto &histo : boundHistograms_.linearArray()) | ||
| histo.zeroBins(); | ||
| for (auto &histo : unboundHistograms_.linearArray()) | ||
| histo.zeroBins(); | ||
| } | ||
|
|
||
| // Return atom types list | ||
| const AtomTypeMix &HistogramSet::atomTypeMix() const { return atomTypeMix_; } | ||
|
|
||
| // Set new fingerprint | ||
| void HistogramSet::setFingerprint(std::string_view fingerprint) { fingerprint_ = fingerprint; } | ||
|
|
||
| // Return fingerprint of partials | ||
| std::string_view HistogramSet::fingerprint() const { return fingerprint_; } | ||
|
|
||
| // Return full histogram specified | ||
| Histogram1D &HistogramSet::fullHistogram(int i, int j) { return fullHistograms_[{i, j}]; } | ||
|
|
||
| // Return bound histogram specified | ||
| Histogram1D &HistogramSet::boundHistogram(int i, int j) { return boundHistograms_[{i, j}]; } | ||
|
|
||
| // Return unbound histogram specified | ||
| Histogram1D &HistogramSet::unboundHistogram(int i, int j) { return unboundHistograms_[{i, j}]; } | ||
|
|
||
| /* | ||
| * Manipulation | ||
| */ | ||
|
|
||
| // Create partials from stored Histogram data | ||
| void HistogramSet::formPartials(PartialSet &partials, double boxVolume) | ||
| { | ||
| dissolve::for_each_pair( | ||
| ParallelPolicies::seq, atomTypeMix_, | ||
| [&](int n, const AtomTypeData &at1, int m, const AtomTypeData &at2) | ||
| { | ||
| // Calculate RDFs from histogram data | ||
| calculateRDF(partials.partial(n, m), fullHistograms_[{n, m}], boxVolume, at1.population(), at2.population(), | ||
| &at1 == &at2 ? 2.0 : 1.0); | ||
| calculateRDF(partials.boundPartial(n, m), boundHistograms_[{n, m}], boxVolume, at1.population(), at2.population(), | ||
| &at1 == &at2 ? 2.0 : 1.0); | ||
| calculateRDF(partials.unboundPartial(n, m), unboundHistograms_[{n, m}], boxVolume, at1.population(), | ||
| at2.population(), &at1 == &at2 ? 2.0 : 1.0); | ||
|
|
||
| // Set flags for bound partials specifying if they are empty (i.e. there are no | ||
| // contributions of that type) | ||
| partials.emptyBoundPartial(n, m) = boundHistograms_[{n, m}].nBinned() == 0; | ||
| }, | ||
| half_); | ||
| } | ||
|
|
||
| // Calculate RDF from supplied Histogram and normalisation data | ||
| void HistogramSet::calculateRDF(Data1D &destination, const Histogram1D &histogram, double boxVolume, int nCentres, | ||
| int nSurrounding, double multiplier) | ||
| { | ||
| auto nBins = histogram.nBins(); | ||
| double delta = histogram.binWidth(); | ||
| const auto &bins = histogram.bins(); | ||
|
|
||
| destination.clear(); | ||
|
|
||
| double shellVolume, factor, r = 0.5 * delta, lowerShellLimit = 0.0, numberDensity = nSurrounding / boxVolume; | ||
| for (auto n = 0; n < nBins; ++n) | ||
| { | ||
| shellVolume = (4.0 / 3.0) * M_PI * (pow(lowerShellLimit + delta, 3.0) - pow(lowerShellLimit, 3.0)); | ||
| factor = nCentres * (shellVolume * numberDensity); | ||
|
|
||
| destination.addPoint(r, bins[n] * (multiplier / factor)); | ||
|
|
||
| r += delta; | ||
| lowerShellLimit += delta; | ||
| } | ||
| } | ||
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,64 @@ | ||
| // SPDX-License-Identifier: GPL-3.0-or-later | ||
| // Copyright (c) 2025 Team Dissolve and contributors | ||
|
|
||
| #pragma once | ||
|
|
||
| #include "classes/atomTypeMix.h" | ||
| #include "classes/partialSet.h" | ||
| #include "math/histogram1D.h" | ||
| #include "templates/array2D.h" | ||
|
|
||
| // Set of Histograms | ||
| class HistogramSet | ||
| { | ||
| public: | ||
| HistogramSet() = default; | ||
| ~HistogramSet(); | ||
|
|
||
| /* | ||
| * Data | ||
| */ | ||
| private: | ||
| // AtomTypeMix used to generate matrices | ||
| AtomTypeMix atomTypeMix_; | ||
| // Fingerprint for these partials (e.g. reflecting Configuration indices at which they were calculated) | ||
| std::string fingerprint_; | ||
| // Histograms used for calculating full atom-atom partials in r | ||
| Array2D<Histogram1D> fullHistograms_; | ||
| // Histograms used for calculating bound atom-atom partials in r | ||
| Array2D<Histogram1D> boundHistograms_; | ||
| // Histograms used for deriving unbound atom-atom partials in r | ||
| Array2D<Histogram1D> unboundHistograms_; | ||
| // Check for full or half matrix | ||
| bool half_{true}; | ||
|
|
||
| public: | ||
| // Set up histograms | ||
| void initialise(const AtomTypeMix &atomTypeMix, double rdfRange, double binWidth); | ||
| // Clear all histogram data | ||
| void clear(); | ||
| // Zero histogram bins | ||
| void zeroBins(); | ||
| // Return atom types mis | ||
| const AtomTypeMix &atomTypeMix() const; | ||
| // Set new fingerprint | ||
| void setFingerprint(std::string_view fingerprint); | ||
|
trisyoungs marked this conversation as resolved.
|
||
| // Return fingerprint of partials | ||
| std::string_view fingerprint() const; | ||
| // Return full histogram specified | ||
| Histogram1D &fullHistogram(int i, int j); | ||
| // Return bound histogram specified | ||
| Histogram1D &boundHistogram(int i, int j); | ||
| // Return unbound histogram specified | ||
| Histogram1D &unboundHistogram(int i, int j); | ||
|
|
||
| /* | ||
| * Manipulation | ||
| */ | ||
| public: | ||
| // Form partials from stored Histogram data | ||
| void formPartials(PartialSet &partials, double boxVolume); | ||
| // Calculate RDF from supplied Histogram and normalisation data | ||
| static void calculateRDF(Data1D &destination, const Histogram1D &histogram, double boxVolume, int nCentres, | ||
| int nSurrounding, double multiplier); | ||
| }; | ||
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Oops, something went wrong.
Oops, something went wrong.
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
What do you think of a static method in
Array2D,flatten, that takes a couple ofArray2D<T>'s (like the abovefull/bound/unboundHistograms_) and flattens theirlinearArray()s into a single vector?I believe the helper method could be implemented like:
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I see what you mean - in fact Adam's
zip()operator would allow us to do this in one loop, albeit with three variables. Your suggestion makes a copy of the data, so unfortunately wouldn't have the intended result. I think I will leave this as-is for now, but your suggestion makes me think that some kind of function which takes a lambda and operates on every element in theArray2Dwould be useful to implement in the future (e.g.fullHistograms_.operate([](auto &histo) { histo.zeroBins(); }))