refactor: Single configuration GR node and updated associates

src/classes/configuration.h

@@ -80,7 +80,7 @@ class Configuration : public Serialisable<const CoreData &>
      */
     private:
     // Species populations present in the Configuration
-    std::vector<std::pair<const Species *, int>> speciesPopulations_;
+    std::map<const Species *, int> speciesPopulations_;
     // AtomType populations in the configuration
     AtomTypeMix atomTypePopulations_;
     // Contents version, incremented whenever Configuration content or Atom positions change
@@ -100,7 +100,7 @@ class Configuration : public Serialisable<const CoreData &>
     // Adjust population of specified Species in the Configuration
     void adjustSpeciesPopulation(const Species *sp, int delta);
     // Return Species populations within the Configuration
-    const std::vector<std::pair<const Species *, int>> &speciesPopulations() const;
+    const std::map<const Species *, int> &speciesPopulations() const;
     // Return population of specified species within the Configuration
     int speciesPopulation(const Species *sp) const;
     // Return if the specified Species is present in the Configuration

src/classes/configuration_contents.cpp

@@ -31,35 +31,27 @@ const AtomTypeMix &Configuration::atomTypePopulations() const { return atomTypeP
 // Adjust population of specified Species in the Configuration
 void Configuration::adjustSpeciesPopulation(const Species *sp, int delta)
 {
-    auto it = std::find_if(speciesPopulations_.begin(), speciesPopulations_.end(),
-                           [sp](const auto &data) { return data.first == sp; });
-    if (it == speciesPopulations_.end())
+    if (speciesPopulations_.contains(sp))
+        speciesPopulations_[sp] += delta;
+    else
     {
         if (delta < 0)
             Messenger::exception("Can't decrease population of Species '{}' as it is not in the list.\n", sp->name());
-        speciesPopulations_.emplace_back(sp, delta);
+        speciesPopulations_[sp] = delta;
     }
-    else
-        it->second += delta;
 }
 
 // Return Species populations within the Configuration
-const std::vector<std::pair<const Species *, int>> &Configuration::speciesPopulations() const { return speciesPopulations_; }
+const std::map<const Species *, int> &Configuration::speciesPopulations() const { return speciesPopulations_; }
 
 // Return population of specified species within the Configuration
 int Configuration::speciesPopulation(const Species *sp) const
 {
-    auto it = std::find_if(speciesPopulations_.begin(), speciesPopulations_.end(),
-                           [sp](const auto &spInfo) { return spInfo.first == sp; });
-    return it == speciesPopulations_.end() ? 0 : it->second;
+    return speciesPopulations_.contains(sp) ? speciesPopulations_.at(sp) : 0;
 }
 
 // Return if the specified Species is present in the Configuration
-bool Configuration::containsSpecies(const Species *sp)
-{
-    return std::find_if(speciesPopulations_.begin(), speciesPopulations_.end(),
-                        [sp](const auto &data) { return data.first == sp; }) != speciesPopulations_.end();
-}
+bool Configuration::containsSpecies(const Species *sp) { return speciesPopulations_.contains(sp); }
 
 // Return the total charge of the Configuration
 double Configuration::totalCharge(bool ppIncludeCoulomb) const
@@ -72,11 +64,11 @@ double Configuration::totalCharge(bool ppIncludeCoulomb) const
 // Return the total atomic mass present in the Configuration
 double Configuration::atomicMass() const
 {
-    double mass = 0.0;
+    auto mass = 0.0;
 
     // Get total molar mass in configuration
-    for (const auto &spPop : speciesPopulations_)
-        mass += spPop.first->mass() * spPop.second;
+    for (const auto &[sp, population] : speciesPopulations_)
+        mass += sp->mass() * population;
 
     // Convert to absolute mass
     return mass / DissolveMath::Avogadro;
@@ -324,4 +316,4 @@ bool Configuration::energyIsStable() const { return energyIsStable_; }
 // Energy gradient
 void Configuration::setEnergyGradient(double grad) { energyGradient_ = grad; }
 
-double Configuration::getEnergyGradient() const { return energyGradient_; }
+double Configuration::getEnergyGradient() const { return energyGradient_; }

src/classes/neutronWeights.cpp

@@ -4,6 +4,7 @@
 #include "classes/neutronWeights.h"
 #include "base/lineParser.h"
 #include "classes/atomType.h"
+#include "classes/isotopologueSet.h"
 #include "classes/species.h"
 #include "data/isotopes.h"
 #include "items/deserialisers.h"
@@ -215,6 +216,29 @@ void NeutronWeights::calculateWeightingMatrices()
                             });
 }
 
+// Create from species populations and isotopologues
+void NeutronWeights::create(const std::map<const Species *, double> &populations, const IsotopologueSet &isotopologues,
+                            const std::vector<std::shared_ptr<AtomType>> &exchangeableTypes)
+{
+    clear();
+
+    for (auto &[sp, pop] : populations)
+    {
+        // Find the defined Isotopologue for this Species - if it doesn't exist, use the Natural one
+        auto isoRef = isotopologues.getIsotopologues(sp);
+        if (isoRef)
+        {
+            const Isotopologues &topes = *isoRef;
+            for (const auto &isoWeight : topes.mix())
+                addIsotopologue(sp, pop, isoWeight.isotopologue(), isoWeight.weight());
+        }
+        else
+            addIsotopologue(sp, pop, sp->naturalIsotopologue(), 1.0);
+    }
+
+    createFromIsotopologues(exchangeableTypes);
+}
+
 // Create AtomType list and matrices based on stored Isotopologues information
 void NeutronWeights::createFromIsotopologues(const std::vector<std::shared_ptr<AtomType>> &exchangeableTypes)
 {

src/classes/neutronWeights.h

@@ -8,6 +8,9 @@
 #include "templates/array2D.h"
 #include <vector>
 
+// Forward Declarations
+class IsotopologueSet;
+
 // Neutron Weights Container
 class NeutronWeights
 {
@@ -60,6 +63,9 @@ class NeutronWeights
     void calculateWeightingMatrices();
 
     public:
+    // Create from species populations and isotopologues
+    void create(const std::map<const Species *, double> &populations, const IsotopologueSet &isotopologues,
+                const std::vector<std::shared_ptr<AtomType>> &exchangeableTypes);
     // Create AtomType list and matrices based on stored Isotopologues information
     void createFromIsotopologues(const std::vector<std::shared_ptr<AtomType>> &exchangeableTypes);
     // Reduce data to be naturally-weighted

src/classes/partialSet.cpp

@@ -12,7 +12,10 @@
 #include "math/mathFunc.h"
 #include "templates/algorithms.h"
 
-PartialSet::PartialSet(const SpeciesPopulations &speciesPopulations) : speciesPopulations_(speciesPopulations) {}
+PartialSet::PartialSet(const std::map<const Species *, double> &realSpeciesPopulations)
+    : realSpeciesPopulations_(realSpeciesPopulations)
+{
+}
 
 PartialSet::~PartialSet()
 {
@@ -268,7 +271,7 @@ Data1D &PartialSet::unboundTotal() { return unboundTotal_; }
 const Data1D &PartialSet::unboundTotal() const { return unboundTotal_; }
 
 // Species populations
-SpeciesPopulations &PartialSet::speciesPopulations() { return speciesPopulations_; }
+const std::map<const Species *, double> &PartialSet::realSpeciesPopulations() const { return realSpeciesPopulations_; }
 
 // Save all partials and total
 bool PartialSet::save(std::string_view prefix, std::string_view tag, std::string_view suffix,

src/classes/partialSet.h

@@ -13,14 +13,12 @@
 class Configuration;
 class Interpolator;
 
-using SpeciesPopulations = std::vector<std::pair<const Species *, double>>;
-
 // Set of Partials
 class PartialSet
 {
     public:
     PartialSet() = default;
-    PartialSet(const SpeciesPopulations &speciesPopulations);
+    PartialSet(const std::map<const Species *, double> &realSpeciesPopulations);
     ~PartialSet();
 
     /*
@@ -52,7 +50,7 @@ class PartialSet
     // Effective density
     double rho_;
     // Species populations
-    std::vector<std::pair<const Species *, double>> speciesPopulations_;
+    std::map<const Species *, double> realSpeciesPopulations_;
 
     public:
     // Set up PartialSet, including initialising histograms for g(r) use
@@ -104,8 +102,8 @@ class PartialSet
     // Return total unbound function
     Data1D &unboundTotal();
     const Data1D &unboundTotal() const;
-    // Species populations
-    SpeciesPopulations &speciesPopulations();
+    // Return real species populations
+    const std::map<const Species *, double> &realSpeciesPopulations() const;
     // Save all partials and total
     bool save(std::string_view prefix, std::string_view tag, std::string_view suffix, std::string_view abscissaUnits) const;
     // Name all object based on the supplied prefix

src/modules/gr/process.cpp

@@ -35,7 +35,7 @@ Module::ExecutionResult GRModule::process(Dissolve &dissolve)
 
     /*
      * Regardless of whether we are a main processing task (summing some combination of Configuration's partials) or
-     * multiple independent Configurations, we must loop over the specified targetConfigurations_ and calculate the partials
+     * multiple independent Configurations, we must loop over the specified targetConfiguration_ and calculate the partials
      * for each.
      */
 

src/nodes/gr/gr.cpp

@@ -5,8 +5,7 @@
 
 GRNode::GRNode(Graph *parentGraph) : Node(parentGraph)
 {
-    addInput<std::vector<Configuration *>>("Configurations", "Set target configuration(s) for the module",
-                                           targetConfigurations_)
+    addInput<Configuration *>("Configuration", "Set target configuration for the module", targetConfiguration_)
         ->setFlags({ParameterBase::Required, ParameterBase::ClearData});
     addOption<Number>("BinWidth", "Bin width (spacing in r) to use", binWidth_);
     addOption<std::optional<Number>>("Range", "Maximum r to calculate g(r) out to", requestedRange_);
@@ -17,8 +16,7 @@ GRNode::GRNode(Graph *parentGraph) : Node(parentGraph)
     addOption<Function1DWrapper>("IntraBroadening", "Type of broadening to apply to intramolecular g(r)", intraBroadening_);
     addOption<std::optional<Number>>("Smoothing", "Specifies the degree of smoothing to apply to calculated g(r)", nSmooths_);
     addOption<bool>("Save", "Whether to save partials and total functions to disk", save_);
-    addOption<bool>("SaveOriginal", "Whether to save original (unbroadened) partials and total functions to disk",
-                    saveOriginal_);
+    addOption<bool>("SaveRaw", "Whether to save raw simulation partial and total functions to disk", saveRaw_);
     addOption<bool>(
         "InternalTest",
         "Perform internal check of calculated partials against a set calculated by a simple unoptimised double-loop",

src/nodes/gr/gr.h

@@ -44,13 +44,11 @@ class GRNode : public Node
 
     private:
     // Target configurations
-    std::vector<Configuration *> targetConfigurations_;
-    // Original g(r)
-    std::optional<PartialSet> originalgr_;
+    Configuration *targetConfiguration_{nullptr};
+    // Raw simulation g(r)
+    std::optional<PartialSet> rawGR_;
     // Unweighted g(r)
     std::optional<PartialSet> unweightedGR_;
-    // Summed unweighted g(r)
-    std::optional<PartialSet> summedUnweightedGR_;
     // Number of historical partial sets to combine into final partials
     std::optional<Number> averagingLength_{5};
     // Weighting scheme to use when averaging partials
@@ -71,40 +69,25 @@ class GRNode : public Node
     std::optional<Number> requestedRange_;
     // Whether to save partials and total functions to disk
     bool save_{false};
-    // Whether to save original (unbroadened) partials and total functions to disk
-    bool saveOriginal_{false};
+    // Whether to save raw partials and total functions to disk
+    bool saveRaw_{false};
 
     /*
      * Functions
      */
     private:
     // Calculate partial g(r) in serial with simple double-loop
-    bool calculateGRTestSerial(Configuration *cfg, PartialSet &partialSet);
+    bool calculateGRTestSerial();
     // Calculate partial g(r) with optimised double-loop
-    bool calculateGRSimple(Configuration *cfg, PartialSet &partialSet, const double rdfRange);
+    bool calculateGRSimple();
     // Calculate partial g(r) utilising Cell neighbour lists
-    bool calculateGRCells(Configuration *cfg, PartialSet &partialSet, const double binWidth);
+    bool calculateGRCells(double grRange);
 
     public:
-    // Get original g(r), constructing if empty
-    PartialSet &originalGR(Configuration *cfg, const double rdfRange, const double rdfBinWidth);
-    // Get unweighted g(r), constructing if empty
-    PartialSet &unweightedGR();
-    // Get summed unweighted g(r), constructing if empty
-    PartialSet &summedUnweightedGR();
-    // Calculate and return effective density based on target Configurations
-    double effectiveDensity() const;
-    // Calculate and return used species populations based on target Configurations
-    SpeciesPopulations speciesPopulations() const;
-    // (Re)calculate partial g(r) for the specified Configuration
-    bool calculateGR(Configuration *cfg, PartialSet &originalgr, PartialsMethod method, const double rdfRange,
-                     const double rdfBinWidth, bool &alreadyUpToDate);
+    // Calculate raw partials
+    bool calculateRawGR(const double grRange, bool &alreadyUpToDate);
     // Calculate smoothed/broadened partial g(r) from supplied partials
-    bool calculateUnweightedGR(Configuration *cfg, const PartialSet &originalgr, PartialSet &weightedgr,
-                               const Function1DWrapper intraBroadening, int smoothing);
-    // Sum unweighted g(r) over the supplied Module's target Configurations
-    bool sumUnweightedGR(std::string_view targetPrefix, std::string_view parentPrefix,
-                         const std::vector<Configuration *> &parentCfgs, PartialSet &summedUnweightedGR);
+    bool calculateUnweightedGR();
     // Test supplied PartialSets against each other
     bool testReferencePartials(PartialSet &setA, PartialSet &setB, double testThreshold);
     // Test calculated partial against supplied reference data