diff --git a/src/classes/partialSet.cpp b/src/classes/partialSet.cpp index 69c77721b6..f2f0d47ded 100644 --- a/src/classes/partialSet.cpp +++ b/src/classes/partialSet.cpp @@ -91,6 +91,10 @@ void PartialSet::initialise(const PartialSet &partialSet) triangular_ = partialSet.triangular_; rho_ = partialSet.rho_; + partials_.clear(triangular_); + boundPartials_.clear(triangular_); + unboundPartials_.clear(triangular_); + // Template data from source PartialSet and set tags dissolve::for_each_pair( ParallelPolicies::seq, atomTypeFractions(), @@ -289,11 +293,7 @@ bool PartialSet::save(std::string_view prefix, std::string_view tag, std::string std::string filename{std::format("{}-{}-{}-{}.{}", prefix, tag, popI.first->name(), popJ.first->name(), suffix)}; Messenger::printVerbose("Writing partial file '{}'...\n", filename); - auto cwd = std::filesystem::current_path(); - auto path = cwd.parent_path().parent_path() / "tests" / "nodes" / "output" / filename; - auto fullPath = path.string(); - - parser.openOutput(fullPath, true); + parser.openOutput(filename, true); if (!parser.isFileGoodForWriting()) return Messenger::error("Couldn't open file '{}' for writing.\n", filename); diff --git a/src/classes/species.h b/src/classes/species.h index d8b8bf3b4b..3aba19f085 100644 --- a/src/classes/species.h +++ b/src/classes/species.h @@ -71,7 +71,7 @@ class Species : public Serialisable // Add a new atom to the Species, returning its index int addAtom(Elements::Element Z, Vector3 r, double q = 0.0, std::shared_ptr atomType = nullptr); // Add new atom type to atom types - const std::shared_ptr addAtomType(Elements::Element Z); + const std::shared_ptr addAtomType(Elements::Element Z, std::string_view name = ""); // Remove the specified atom from the species void removeAtom(int index); // Remove set of atom indices diff --git a/src/classes/species_atomic.cpp b/src/classes/species_atomic.cpp index 5634852dc3..8cb7a75a34 100644 --- a/src/classes/species_atomic.cpp +++ b/src/classes/species_atomic.cpp @@ -40,17 +40,15 @@ int Species::addAtom(Elements::Element Z, Vector3 r, double q, std::shared_ptr Species::addAtomType(Elements::Element Z) +const std::shared_ptr Species::addAtomType(Elements::Element Z, std::string_view name) { - auto newAtomType = std::make_shared(); - atomTypes_.push_back(newAtomType); - // Create a suitable unique name - newAtomType->setName(DissolveSys::uniqueName(Elements::symbol(Z), atomTypes_, - [&](const auto &at) { return newAtomType == at ? "" : at->name(); })); + auto uniqueName = DissolveSys::uniqueName(name == "" ? Elements::symbol(Z) : name, atomTypes_, + [&](const auto &at) { return at->name(); }); - // Set data - newAtomType->setZ(Z); + // Create atom type and set data + auto newAtomType = std::make_shared(Z, uniqueName); + atomTypes_.push_back(newAtomType); newAtomType->setIndex(atomTypes_.size() - 1); return newAtomType; diff --git a/src/io/import/coordinates.h b/src/io/import/coordinates.h index 0b4014f990..e7b309cfe0 100644 --- a/src/io/import/coordinates.h +++ b/src/io/import/coordinates.h @@ -6,8 +6,6 @@ #include "base/enumOptions.h" #include "io/fileAndFormat.h" -// Forward Declarations - // Coordinate Import Formats class CoordinateImportFileFormat : public FileAndFormat { diff --git a/src/math/mathFunc.cpp b/src/math/mathFunc.cpp index b8be283994..8c4e239a53 100644 --- a/src/math/mathFunc.cpp +++ b/src/math/mathFunc.cpp @@ -78,6 +78,16 @@ double sgn(double a, double signOf) { return signOf >= 0.0 ? fabs(a) : -fabs(a); // Return the cyclic permutation of the integer 'i', span 3 int cp3(int i) { return (i % 3); } +// Return the integer wrapped into the specified range +int wrap(int i, int lower, int upper) +{ + const auto range = upper - lower + 1; + + if (i < 0) + i += range * ((lower - i) / range + 1); + return lower + (i - lower) % range; +} + /* * Conversion */ diff --git a/src/math/mathFunc.h b/src/math/mathFunc.h index be3ed6ea6e..49a41b09bb 100644 --- a/src/math/mathFunc.h +++ b/src/math/mathFunc.h @@ -39,6 +39,8 @@ int sgn(double x); double sgn(double a, double signOf); // Return the cyclic permutation of the integer 'i', span 3 int cp3(int i); +// Return the integer wrapped into the specified range +int wrap(int i, int lower, int upper); /* * Conversion diff --git a/src/nodes/atomicMC/process.cpp b/src/nodes/atomicMC/process.cpp index 1a90faf8f5..32c356c422 100644 --- a/src/nodes/atomicMC/process.cpp +++ b/src/nodes/atomicMC/process.cpp @@ -32,7 +32,7 @@ NodeConstants::ProcessResult AtomicMCNode::process() message("\n"); // Prepare for energy calculation, generate kernel - auto kernel = dissolveGraph()->prepareEnergyCalculation(targetConfiguration_); + auto kernel = dissolveGraph()->createEnergyCalculation(targetConfiguration_); auto nAttempts = 0, nAccepted = 0; auto totalDelta = 0.0; diff --git a/src/nodes/dissolve.cpp b/src/nodes/dissolve.cpp index 227d375fc0..cb76c06e88 100644 --- a/src/nodes/dissolve.cpp +++ b/src/nodes/dissolve.cpp @@ -32,30 +32,32 @@ const DoubleKeyedMap &DissolveGraph::pairPotentialStore() { retur * Functions */ -// Return maximum distance for tabulated PairPotentials -const double DissolveGraph::pairPotentialRange() const { return pairPotentialRange_; } +// Return range of tabulated pari potentials +double DissolveGraph::pairPotentialRange() const { return pairPotentialRange_; } -// Return energy kernel containing potential map -std::unique_ptr DissolveGraph::prepareEnergyCalculation(Configuration *cfg) +// Ensure that the specified Configuration has updated type indexing, cells etc. +void DissolveGraph::updateIndexingAndCells(Configuration *cfg) const { - auto atomTypes = cfg->atomTypeVector(); - // Update atom type indexing cfg->updateTypeIndexing(); + // Regenerate cells in the configuration if necessary + cfg->updateCells(pairPotentialRange()); +} + +// Create an energy kernel suitable for the supplied Configuration +std::unique_ptr DissolveGraph::createEnergyCalculation(Configuration *cfg) +{ + // Update types and cells in Configuration + updateIndexingAndCells(cfg); + auto atomTypes = cfg->atomTypeVector(); + // Update pair potentials dissolve::for_each_pair(ParallelPolicies::seq, atomTypes, [&](int i, const auto &atI, int j, const auto &atJ) { updatePairPotentials(*atI, *atJ); }); - // Generate configuration potential map - PotentialMap potentialMap(atomTypes, pairPotentialStore(), pairPotentialRange()); - - // Regenerate cells in the configuration if necessary - cfg->updateCells(potentialMap.range()); - - auto kernel = KernelProducer::energyKernel(cfg, potentialMap); - - return kernel; + // Generate and return kernel + return KernelProducer::energyKernel(cfg, PotentialMap(atomTypes, pairPotentialStore(), pairPotentialRange())); } // Update pair potential store diff --git a/src/nodes/dissolve.h b/src/nodes/dissolve.h index 0bd80b2da7..c991a486c9 100644 --- a/src/nodes/dissolve.h +++ b/src/nodes/dissolve.h @@ -56,10 +56,12 @@ class DissolveGraph : public Graph * Functions */ public: - // Return maximum distance for tabulated PairPotentials - const double pairPotentialRange() const; - // Return suitable energy kernel for the supplied Configuration - std::unique_ptr prepareEnergyCalculation(Configuration *cfg); + // Return range of tabulated pari potentials + double pairPotentialRange() const; + // Ensure that the specified Configuration has updated type indexing, cells etc. + void updateIndexingAndCells(Configuration *cfg) const; + // Create an energy kernel suitable for the supplied Configuration + std::unique_ptr createEnergyCalculation(Configuration *cfg); private: // Update pair potential store diff --git a/src/nodes/edge.cpp b/src/nodes/edge.cpp index b77172bbef..1ca8a7259c 100644 --- a/src/nodes/edge.cpp +++ b/src/nodes/edge.cpp @@ -32,7 +32,7 @@ class EdgeConstructor : public Edge std::unique_ptr Edge::create(Graph *parent, const EdgeDefinition &definition) { // Get source node and output - auto sourceNode = parent->node(definition.sourceNode); + auto sourceNode = parent->findNode(definition.sourceNode); if (!sourceNode) { Messenger::error("Source node '{}' does not exist in the graph.\n", definition.sourceNode); @@ -54,7 +54,7 @@ std::unique_ptr Edge::create(Graph *parent, const EdgeDefinition &definiti } // Get target node and input - auto targetNode = parent->node(definition.targetNode); + auto targetNode = parent->findNode(definition.targetNode); if (!targetNode) { Messenger::error("Target node '{}' does not exist in the graph.\n", definition.targetNode); diff --git a/src/nodes/energy/process.cpp b/src/nodes/energy/process.cpp index 2fde9a9124..85951a0b99 100644 --- a/src/nodes/energy/process.cpp +++ b/src/nodes/energy/process.cpp @@ -24,8 +24,8 @@ NodeConstants::ProcessResult EnergyNode::process() * This is a serial routine (subroutines called from within are parallel). */ - auto kernel = dissolveGraph()->prepareEnergyCalculation(targetConfiguration_); - auto potentialMap = kernel->potentialMap(); + auto kernel = dissolveGraph()->createEnergyCalculation(targetConfiguration_); + auto &potentialMap = kernel->potentialMap(); // Calculate pair potential energy Timer interTimer; diff --git a/src/nodes/gr/gr.cpp b/src/nodes/gr/gr.cpp index c7fa9656e0..200a2eb5f4 100644 --- a/src/nodes/gr/gr.cpp +++ b/src/nodes/gr/gr.cpp @@ -29,10 +29,6 @@ GRNode::GRNode(Graph *parentGraph) addOption>("Smoothing", "Specifies the degree of smoothing to apply to calculated g(r)", nSmooths_); addOption("Save", "Whether to save partials and total functions to disk", save_); addOption("SaveRaw", "Whether to save raw simulation partial and total functions to disk", saveRaw_); - addOption( - "InternalTest", - "Perform internal check of calculated partials against a set calculated by a simple unoptimised double-loop", - internalTest_); addOption("Method", "Calculation method for partial radial distribution functions", partialsMethod_); diff --git a/src/nodes/gr/gr.h b/src/nodes/gr/gr.h index 4c3194fc67..44a4209157 100644 --- a/src/nodes/gr/gr.h +++ b/src/nodes/gr/gr.h @@ -53,8 +53,6 @@ class GRNode : public Node std::optional averagingLength_{5}; // Bin width (spacing in r) to use Number binWidth_{0.001}; - // Perform internal check of calculated partials against a set calculated by a simple unoptimised double-loop - bool internalTest_{false}; // Type of broadening to apply to intramolecular g(r) Function1DWrapper intraBroadening_{Functions1D::Form::Gaussian, {0.18}}; // Degree of smoothing to apply @@ -89,12 +87,6 @@ class GRNode : public Node bool calculateRawGR(const double grRange, bool &alreadyUpToDate); // Calculate smoothed/broadened partial g(r) from supplied partials bool calculateUnweightedGR(); - // Test supplied PartialSets against each other - bool testReferencePartials(const std::vector &types, PartialSet &setA, PartialSet &setB, - double testThreshold); - // Test calculated partial against supplied reference data - bool testReferencePartial(const PartialSet &partials, double testThreshold, const Data1D &testData, - std::string_view typeIorTotal, std::string_view typeJ = "", std::string_view target = ""); /* * Processing diff --git a/src/nodes/gr/helpers.cpp b/src/nodes/gr/helpers.cpp index ec5d68bd04..561299a723 100644 --- a/src/nodes/gr/helpers.cpp +++ b/src/nodes/gr/helpers.cpp @@ -8,9 +8,9 @@ #include "classes/cell.h" #include "main/dissolve.h" #include "math/combinations.h" -#include "math/error.h" #include "math/filters.h" #include "module/group.h" +#include "nodes/dissolve.h" #include "nodes/gr/gr.h" #include "templates/algorithms.h" #include "templates/combinable.h" @@ -241,6 +241,7 @@ bool GRNode::calculateGRCells(double grRange, const Array2DupdateIndexingAndCells(targetConfiguration_); calculateGRCells(grRange, fullLUT); + } else if (partialsMethod_ == PartialsMethod::AutoMethod) { - targetConfiguration_->nAtoms() > 10000 ? calculateGRCells(grRange, fullLUT) : calculateGRSimple(fullLUT); + if (targetConfiguration_->nAtoms() > 10000) + { + dissolveGraph()->updateIndexingAndCells(targetConfiguration_); + calculateGRCells(grRange, fullLUT); + } + else + calculateGRSimple(fullLUT); } timer.stop(); message("Finished calculation of partials ({} elapsed).\n", timer.totalTimeString()); @@ -444,116 +454,3 @@ bool GRNode::calculateUnweightedGR() return true; } - -// Test supplied PartialSets against each other -bool GRNode::testReferencePartials(const std::vector &types, PartialSet &setA, PartialSet &setB, - double testThreshold) -{ - for_each_pair_early(types, - [&](int n, const auto *typeI, int m, const auto *typeJ) -> EarlyReturn - { - DoubleKeyedMapKey key{typeI->name(), typeJ->name()}; - - // Full partial - auto errorReport = Error::percent(setA.partials().get(key), setB.partials().get(key)); - message("{}", Error::errorReportString(errorReport)); - message("Test reference full partial '{}-{}' has {} error of {:7.3f}{} with calculated data and is " - "{} (threshold is {:6.3f}%)\n\n", - typeI->name(), typeJ->name(), Error::errorTypes().keyword(errorReport.errorType), - errorReport.error, errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", - errorReport.error <= testThreshold ? "OK" : "NOT OK", testThreshold); - if (errorReport.error > testThreshold) - return false; - - // Bound partial - errorReport = Error::percent(setA.boundPartials().get(key), setB.boundPartials().get(key)); - message("{}", Error::errorReportString(errorReport)); - message("Test reference bound partial '{}-{}' has {} error of {:7.3f}{} with calculated data and " - "is {} (threshold is {:6.3f}%)\n\n", - typeI->name(), typeJ->name(), Error::errorTypes().keyword(errorReport.errorType), - errorReport.error, errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", - errorReport.error <= testThreshold ? "OK" : "NOT OK", testThreshold); - if (errorReport.error > testThreshold) - return false; - - // Unbound reference - errorReport = Error::percent(setA.unboundPartials().get(key), setB.unboundPartials().get(key)); - message("{}", Error::errorReportString(errorReport)); - message("Test reference unbound partial '{}-{}' has {} error of {:7.3f}{} with calculated data and " - "is {} (threshold is {:6.3f}%)\n\n", - typeI->name(), typeJ->name(), Error::errorTypes().keyword(errorReport.errorType), - errorReport.error, errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", - errorReport.error <= testThreshold ? "OK" : "NOT OK", testThreshold); - if (errorReport.error > testThreshold) - return false; - - return EarlyReturn::Continue; - }); - - // Total reference data supplied? - auto errorReport = Error::percent(setA.total(), setB.total()); - message("{}", Error::errorReportString(errorReport)); - message("Test reference total has {} error of {:7.3f}{} with calculated data and is {} (threshold is {:6.3f}%)\n\n", - Error::errorTypes().keyword(errorReport.errorType), errorReport.error, - errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", - errorReport.error <= testThreshold ? "OK" : "NOT OK", testThreshold); - if (errorReport.error > testThreshold) - return false; - - return true; -} - -// Test calculated partial against supplied reference data -bool GRNode::testReferencePartial(const PartialSet &partials, double testThreshold, const Data1D &testData, - std::string_view typeIorTotal, std::string_view typeJ, std::string_view target) -{ - // We either expect two AtomType names and a target next, or the target 'total' - auto testResult = false; - if (DissolveSys::sameString(typeIorTotal, "total") && typeJ.empty() && target.empty()) - { - auto errorReport = Error::percent(partials.total(), testData); - message("{}", Error::errorReportString(errorReport)); - testResult = (errorReport.error <= testThreshold); - message("Test reference data '{}' has {} error of {:7.3f}{} with calculated data and is {} (threshold is " - "{:6.3f}%)\n\n", - testData.tag(), Error::errorTypes().keyword(errorReport.errorType), errorReport.error, - errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", testResult ? "OK" : "NOT OK", - testThreshold); - } - else - { - DoubleKeyedMapKey key{typeIorTotal, typeJ}; - Error::ErrorReport errorReport; - if (DissolveSys::sameString(target, "bound")) - { - errorReport = Error::percent(partials.boundPartials().get(key), testData); - message("{}", Error::errorReportString(errorReport)); - } - - else if (DissolveSys::sameString(target, "unbound")) - { - errorReport = Error::percent(partials.unboundPartials().get(key), testData); - message("{}", Error::errorReportString(errorReport)); - } - else if (DissolveSys::sameString(target, "full")) - { - errorReport = Error::percent(partials.partials().get(key), testData); - message("{}", Error::errorReportString(errorReport)); - } - - else - { - error("Unrecognised test data name '{}'.\n", testData.tag()); - return false; - } - - testResult = (errorReport.error <= testThreshold); - message("Test reference data '{}' has {} error of {:7.3f}{} with calculated data and is {} (threshold is " - "{:6.3f}%)\n\n", - testData.tag(), Error::errorTypes().keyword(errorReport.errorType), errorReport.error, - errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", testResult ? "OK" : "NOT OK", - testThreshold); - } - - return testResult; -} diff --git a/src/nodes/gr/process.cpp b/src/nodes/gr/process.cpp index 844c2fc036..048907ab7d 100644 --- a/src/nodes/gr/process.cpp +++ b/src/nodes/gr/process.cpp @@ -85,19 +85,6 @@ NodeConstants::ProcessResult GRNode::process() if ((averagingLength_.value_or(1) > 1) && (!alreadyUpToDate)) (*rawGR_) = rawGRHistory_.push((*rawGR_), averagingLength_.value().asInteger()); - /* - // Perform internal test of original g(r)? - if (internalTest_) - { - // Copy the already-calculated g(r), then calculate a new set using the Test method - PartialSet referencePartials = originalgr; - calculateGR(dissolve.processingModuleData(), moduleContext.processPool(), cfg, GRModule::TestMethod, - grRange, binWidth_, alreadyUpToDate); - if (!testReferencePartials(referencePartials, originalgr, 1.0e-6)) - return ExecutionResult::Failed; - } - */ - // Form unweighted g(r) from original g(r), applying any requested smoothing and/or intramolecular broadening calculateUnweightedGR(); diff --git a/src/nodes/graph.cpp b/src/nodes/graph.cpp index 0e2bd0612b..0e2f0958e9 100644 --- a/src/nodes/graph.cpp +++ b/src/nodes/graph.cpp @@ -224,7 +224,7 @@ Edge *Graph::findEdge(const EdgeDefinition &definition) const } // Return named node, if it exists -Node *Graph::node(std::string_view nodeName) +Node *Graph::findNode(std::string_view nodeName) { // Return ourself if this is our name if (name() == nodeName) diff --git a/src/nodes/graph.h b/src/nodes/graph.h index 7565e8615e..018d6d395c 100644 --- a/src/nodes/graph.h +++ b/src/nodes/graph.h @@ -99,7 +99,7 @@ class Graph : public Node // Find edge between nodes Edge *findEdge(const EdgeDefinition &definition) const; // Return named node, if it exists - Node *node(std::string_view nodeName); + Node *findNode(std::string_view nodeName); // Return container of nodes Nodes &nodes(); // Return container of edges between nodes diff --git a/src/nodes/importConfigurationCoordinates.cpp b/src/nodes/importConfigurationCoordinates.cpp index e946a88fdd..503a095dbd 100644 --- a/src/nodes/importConfigurationCoordinates.cpp +++ b/src/nodes/importConfigurationCoordinates.cpp @@ -25,8 +25,5 @@ NodeConstants::ProcessResult ImportConfigurationCoordinatesNode::process() { CoordinateImportFileFormat fileSource(filePath_, format_); - if (fileSource.importData(configuration_)) - return NodeConstants::ProcessResult::Success; - - return NodeConstants::ProcessResult::Unchanged; + return fileSource.importData(configuration_) ? NodeConstants::ProcessResult::Success : NodeConstants::ProcessResult::Failed; } diff --git a/src/nodes/iterableGraph.cpp b/src/nodes/iterableGraph.cpp index 8bc74b32a9..6596186904 100644 --- a/src/nodes/iterableGraph.cpp +++ b/src/nodes/iterableGraph.cpp @@ -55,7 +55,7 @@ void IterableGraph::releaseLoopBack(const std::string &name) // Add edge between nodes bool IterableGraph::addEdge(const EdgeDefinition &definition) { - if (dynamic_cast(parentGraph()->node(definition.sourceNode))) + if (dynamic_cast(parentGraph()->findNode(definition.sourceNode))) setLoopBacks(); else if (loopBacks_->findInput(definition.targetInput)) { diff --git a/src/nodes/md/process.cpp b/src/nodes/md/process.cpp index 95f6c13dca..d322d905e3 100644 --- a/src/nodes/md/process.cpp +++ b/src/nodes/md/process.cpp @@ -48,7 +48,7 @@ NodeConstants::ProcessResult MDNode::process() joinStrings(restrictToSpecies_, " ", [](const auto &sp) { return sp->name(); })); message("\n"); - auto kernel = dissolveGraph()->prepareEnergyCalculation(targetConfiguration_); + auto kernel = dissolveGraph()->createEnergyCalculation(targetConfiguration_); /* if (onlyWhenEnergyStable_) diff --git a/tests/graphData.h b/tests/graphData.h new file mode 100644 index 0000000000..473eb70258 --- /dev/null +++ b/tests/graphData.h @@ -0,0 +1,212 @@ +// SPDX-License-Identifier: GPL-3.0-or-later +// Copyright (c) 2026 Team Dissolve and contributors + +#pragma once + +#include "main/dissolve.h" +#include "nodes/atomicSpecies.h" +#include "nodes/configuration.h" +#include "nodes/dissolve.h" +#include "nodes/insert.h" +#include "tests/speciesData.h" +#include + +namespace UnitTest +{ +// Basic object setup for any Graph-based test +class GraphTestData +{ + public: + GraphTestData() : dissolve(coreData), graphRoot(dissolve) { Node::echo_ = true; } + CoreData coreData; + Dissolve dissolve; + DissolveGraph graphRoot; +}; + +// Create an Argon graph in the supplied root node +inline void createArgonGraph(Graph *root, int population = 1000) +{ + /* + * Configuration (Bulk) + * ------------------ + * - Configuration-o ---+ + * - | \ Insert (Insert) + * -----------------/ \ ------------------ + * +---- o-Configuration | + * AtomicSpecies (Ar) +------ o-Species | + * ------------------ / -----------------/ + * - Species-o --+ + * - | + * -----------------/ + */ + + // Create nodes + auto arNode = root->addNode(std::make_unique(root, Elements::Ar), "Ar"); + ASSERT_TRUE(arNode); + auto configurationNode = root->createNode("Configuration", "Bulk"); + ASSERT_TRUE(configurationNode); + auto insertNode = root->createNode("Insert", "Insert"); + ASSERT_TRUE(insertNode); + + // Create edges + ASSERT_TRUE(root->addEdge({"Ar", "Species", "Insert", "Species"})); + ASSERT_TRUE(root->addEdge({"Bulk", "Configuration", "Insert", "Configuration"})); + + // Set configuration contents + ASSERT_TRUE(insertNode->setInput("Population", population)); + ASSERT_TRUE(insertNode->setInput("Density", 0.0213)); + ASSERT_TRUE(insertNode->setOption("DensityUnits", Units::DensityUnits::AtomsPerAngstromUnits)); +} + +// Create a water graph in the supplied root node +inline void createWater1000Graph(Graph *root, CoordinateImportFileFormat initialCoordinates) +{ + // Create species and configuration + auto waterNode = createWater(root); + ASSERT_TRUE(waterNode); + auto configurationNode = root->createNode("Configuration", "Bulk"); + ASSERT_TRUE(configurationNode); + auto insertNode = root->createNode("Insert"); + ASSERT_TRUE(insertNode); + ASSERT_TRUE(root->addEdge({"Water", "Species", "Insert", "Species"})); + ASSERT_TRUE(root->addEdge({"Bulk", "Configuration", "Insert", "Configuration"})); + ASSERT_TRUE(insertNode->setInput("Population", 1000)); + ASSERT_TRUE(insertNode->setInput("Density", 0.1)); + ASSERT_TRUE(insertNode->setOption("DensityUnits", Units::DensityUnits::AtomsPerAngstromUnits)); + + // Import reference coordinates + auto importCoordinates = root->createNode("ImportConfigurationCoordinates", "Import"); + ASSERT_TRUE(importCoordinates->setOption("FilePath", std::string(initialCoordinates.filename()))); + ASSERT_TRUE(importCoordinates->setOption( + "FileFormat", + CoordinateImportFileFormat::coordinateImportFileFormat().enumerationByIndex(initialCoordinates.formatIndex()))); + ASSERT_TRUE(root->addEdge({"Insert", "Configuration", "Import", "Configuration"})); + + // Add GR node and link to the import node + auto grNode = root->createNode("GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(root->addEdge({"Import", "Configuration", "GR", "Configuration"})); + + // Create the SQ node + auto sqNode = root->createNode("SQ"); + ASSERT_TRUE(sqNode); + ASSERT_TRUE(root->addEdge({"GR", "UnweightedGR", "SQ", "UnweightedGR"})); + + // Add in NeutronSQ + auto h2o = root->createNode("NeutronSQ", "H2O"); + ASSERT_TRUE(h2o); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "H2O", "UnweightedGR"})); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "H2O", "UnweightedSQ"})); + + auto d2o = root->createNode("NeutronSQ", "D2O"); + ASSERT_TRUE(d2o); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "D2O", "UnweightedGR"})); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "D2O", "UnweightedSQ"})); + + auto hdo = root->createNode("NeutronSQ", "HDO"); + ASSERT_TRUE(hdo); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "HDO", "UnweightedGR"})); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "HDO", "UnweightedSQ"})); + + // Add in XRaySQ + // auto h2ox = root->createNode("XRaySQ", "H2OX"); + // ASSERT_TRUE(h2ox); + // ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "H2OX", "UnweightedGR"})); + // ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "H2OX", "UnweightedSQ"})); +} + +// Create a water graph in the supplied root node +inline void createWaterMethanolGraph(Graph *root) +{ + // Create species and configuration + auto waterNode = createWater(root); + ASSERT_TRUE(waterNode); + auto methanolNode = createMethanol(root); + ASSERT_TRUE(methanolNode); + + auto configurationNode = root->createNode("Configuration", "Bulk"); + ASSERT_TRUE(configurationNode); + auto insertWaterNode = root->createNode("Insert", "InsertWater"); + ASSERT_TRUE(insertWaterNode); + ASSERT_TRUE(root->addEdge({"Water", "Species", "InsertWater", "Species"})); + ASSERT_TRUE(root->addEdge({"Bulk", "Configuration", "InsertWater", "Configuration"})); + ASSERT_TRUE(insertWaterNode->setInput("Population", 300)); + ASSERT_TRUE(insertWaterNode->setInput("Density", 0.1)); + ASSERT_TRUE(insertWaterNode->setOption("DensityUnits", Units::DensityUnits::AtomsPerAngstromUnits)); + auto insertMethanolNode = root->createNode("Insert", "InsertMethanol"); + ASSERT_TRUE(insertMethanolNode); + ASSERT_TRUE(root->addEdge({"Methanol", "Species", "InsertMethanol", "Species"})); + ASSERT_TRUE(root->addEdge({"InsertWater", "Configuration", "InsertMethanol", "Configuration"})); + ASSERT_TRUE(insertMethanolNode->setInput("Population", 600)); + ASSERT_TRUE(insertMethanolNode->setInput("Density", 0.1)); + ASSERT_TRUE(insertMethanolNode->setOption("DensityUnits", Units::DensityUnits::AtomsPerAngstromUnits)); + + // Import reference coordinates + auto importCoordinates = root->createNode("ImportConfigurationCoordinates", "Import"); + ASSERT_TRUE(importCoordinates->setOption("FilePath", "epsr25/water300methanol600/watermeth.ato")); + ASSERT_TRUE(importCoordinates->setOption( + "FileFormat", CoordinateImportFileFormat::CoordinateImportFormat::EPSR)); + ASSERT_TRUE(root->addEdge({"InsertMethanol", "Configuration", "Import", "Configuration"})); + + // Add GR node and link to the import node + auto grNode = root->createNode("GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(root->addEdge({"Import", "Configuration", "GR", "Configuration"})); +} + +// Create a benzene graph in the supplied root node +inline void createBenzeneGraph(Graph *root) +{ + // Create species and configuration + auto benzeneNode = createBenzene(root); + ASSERT_TRUE(benzeneNode); + auto configurationNode = root->createNode("Configuration", "Bulk"); + ASSERT_TRUE(configurationNode); + auto insertNode = root->createNode("Insert"); + ASSERT_TRUE(insertNode); + ASSERT_TRUE(root->addEdge({"Benzene", "Species", "Insert", "Species"})); + ASSERT_TRUE(root->addEdge({"Bulk", "Configuration", "Insert", "Configuration"})); + ASSERT_TRUE(insertNode->setInput("Population", 200)); + ASSERT_TRUE(insertNode->setInput("Density", 0.876)); + + // Import reference coordinates + auto importCoordinates = root->createNode("ImportConfigurationCoordinates", "Import"); + ASSERT_TRUE(importCoordinates->setOption("FilePath", "epsr25/benzene200-neutron/boxbenz.ato")); + ASSERT_TRUE(importCoordinates->setOption( + "FileFormat", CoordinateImportFileFormat::CoordinateImportFormat::EPSR)); + ASSERT_TRUE(root->addEdge({"Insert", "Configuration", "Import", "Configuration"})); + + // Add GR node and link to the import node + auto grNode = root->createNode("GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(root->addEdge({"Import", "Configuration", "GR", "Configuration"})); + + // Create the SQ node + auto sqNode = root->createNode("SQ"); + ASSERT_TRUE(sqNode); + ASSERT_TRUE(root->addEdge({"GR", "UnweightedGR", "SQ", "UnweightedGR"})); + + // Add in NeutronSQ + auto H = root->createNode("NeutronSQ", "H"); + ASSERT_TRUE(H); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "H", "UnweightedGR"})); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "H", "UnweightedSQ"})); + + auto D = root->createNode("NeutronSQ", "D"); + ASSERT_TRUE(D); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "D", "UnweightedGR"})); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "D", "UnweightedSQ"})); + + auto HD = root->createNode("NeutronSQ", "HD"); + ASSERT_TRUE(HD); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "HD", "UnweightedGR"})); + ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "HD", "UnweightedSQ"})); + + // Add in XRaySQ? + // auto X = root->createNode("XRaySQ", "X"); + // ASSERT_TRUE(X); + // ASSERT_TRUE(root->addEdge({"SQ", "UnweightedGR", "X", "UnweightedGR"})); + // ASSERT_TRUE(root->addEdge({"SQ", "UnweightedSQ", "X", "UnweightedSQ"})); +} + +} // namespace UnitTest diff --git a/tests/math/CMakeLists.txt b/tests/math/CMakeLists.txt index e3da128a4d..04dff89f8b 100644 --- a/tests/math/CMakeLists.txt +++ b/tests/math/CMakeLists.txt @@ -13,3 +13,4 @@ dissolve_add_test(SRC polynomial.cpp) dissolve_add_test(SRC sampledValues.cpp) dissolve_add_test(SRC svd.cpp) dissolve_add_test(SRC vector3.cpp) +dissolve_add_test(SRC wrap.cpp) diff --git a/tests/math/wrap.cpp b/tests/math/wrap.cpp new file mode 100644 index 0000000000..4eadf5b05a --- /dev/null +++ b/tests/math/wrap.cpp @@ -0,0 +1,30 @@ +// SPDX-License-Identifier: GPL-3.0-or-later +// Copyright (c) 2026 Team Dissolve and contributors + +#include "math/vector3.h" +#include "tests/testData.h" +#include + +namespace UnitTest +{ +TEST(WrapIntegers, WrapIntegers) +{ + EXPECT_EQ(DissolveMath::wrap(2, 0, 4), 2); + EXPECT_EQ(DissolveMath::wrap(0, 0, 4), 0); + EXPECT_EQ(DissolveMath::wrap(-5, 0, 4), 0); + EXPECT_EQ(DissolveMath::wrap(-15, 0, 4), 0); + EXPECT_EQ(DissolveMath::wrap(20, 0, 4), 0); + EXPECT_EQ(DissolveMath::wrap(4, 0, 4), 4); + EXPECT_EQ(DissolveMath::wrap(8, 0, 4), 3); + EXPECT_EQ(DissolveMath::wrap(9, 0, 4), 4); + + EXPECT_EQ(DissolveMath::wrap(2, -5, 5), 2); + EXPECT_EQ(DissolveMath::wrap(0, -5, 5), 0); + EXPECT_EQ(DissolveMath::wrap(-5, -5, 5), -5); + EXPECT_EQ(DissolveMath::wrap(-16, -5, 5), -5); + EXPECT_EQ(DissolveMath::wrap(5, -5, 5), 5); + EXPECT_EQ(DissolveMath::wrap(16, -5, 5), 5); + EXPECT_EQ(DissolveMath::wrap(6, -5, 5), -5); +} + +} // namespace UnitTest diff --git a/tests/modules/CMakeLists.txt b/tests/modules/CMakeLists.txt index 70169dd9bb..fbfa2dd07a 100644 --- a/tests/modules/CMakeLists.txt +++ b/tests/modules/CMakeLists.txt @@ -10,7 +10,6 @@ dissolve_add_test(SRC dAngle.cpp) dissolve_add_test(SRC energy.cpp) dissolve_add_test(SRC epsr.cpp) dissolve_add_test(SRC forces.cpp) -dissolve_add_test(SRC gr.cpp) dissolve_add_test(SRC histogramCN.cpp) dissolve_add_test(SRC intraAngle.cpp) dissolve_add_test(SRC intraDistance.cpp) diff --git a/tests/modules/gr.cpp b/tests/modules/gr.cpp deleted file mode 100644 index a6894dd7ad..0000000000 --- a/tests/modules/gr.cpp +++ /dev/null @@ -1,250 +0,0 @@ -// SPDX-License-Identifier: GPL-3.0-or-later -// Copyright (c) 2026 Team Dissolve and contributors - -#include "modules/gr/gr.h" -#include "tests/testData.h" -#include -#include - -namespace UnitTest -{ -class GRModuleTest : public ::testing::Test -{ - protected: - DissolveSystemTest systemTest; -}; - -TEST_F(GRModuleTest, Methods) -{ - ASSERT_NO_THROW_VERBOSE(systemTest.setUp("dissolve/input/rdfMethod.txt")); - auto *grModule = systemTest.getModule("GR01"); - - // Simple method - grModule->keywords().setEnumeration("Method", GRModule::PartialsMethod::SimpleMethod); - ASSERT_TRUE(systemTest.dissolve().iterate(1)); - - // Cells method - systemTest.coreData().configurations().front()->notifyAtomicPositionsChanged(); - grModule->keywords().setEnumeration("Method", GRModule::PartialsMethod::CellsMethod); - ASSERT_TRUE(systemTest.dissolve().iterate(1)); -} - -TEST_F(GRModuleTest, Water) -{ - ASSERT_NO_THROW_VERBOSE(systemTest.setUp("dissolve/input/correlations-water.txt")); - ASSERT_TRUE(systemTest.dissolve().iterate(1)); - - // Partial g(r) (unbound terms) - EXPECT_TRUE(systemTest.checkData1D( - "GR01//Bulk//OriginalGR//OW-OW//Unbound", - {"epsr25/water1000-neutron/water.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 6.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GR01//Bulk//OriginalGR//OW-HW//Unbound", - {"epsr25/water1000-neutron/water.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 2.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GR01//Bulk//OriginalGR//HW-HW//Unbound", - {"epsr25/water1000-neutron/water.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 2.0e-2)); - - // Partial g(r) (intramolecular terms) - EXPECT_TRUE(systemTest.checkData1D( - "GR01//Bulk//OriginalGR//OW-HW//Bound", - {"epsr25/water1000-neutron/water.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.1)); - EXPECT_TRUE(systemTest.checkData1D( - "GR01//Bulk//OriginalGR//HW-HW//Bound", - {"epsr25/water1000-neutron/water.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 1.5e-2)); - - // Partial g(r) (intramolecular terms) - EXPECT_TRUE(systemTest.checkData1D( - "GR01//Bulk//OriginalGR//OW-OW//Bound", - {"epsr25/water1000-neutron/water.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 1.0e-5, - Error::ErrorType::RMSEError)); -} - -TEST_F(GRModuleTest, WaterMethanol) -{ - ASSERT_NO_THROW_VERBOSE(systemTest.setUp("dissolve/input/correlations-waterMethanol.txt")); - ASSERT_TRUE(systemTest.dissolve().iterate(1)); - - /* - * Partial Radial Distribution Functions - * Order of partials in EPSR files is: - * 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 - * OW-OW OW-HW OW-CT OW-HC OW-OH OW-HO HW-HW HW-CT HW-HC HW-OH HW-HO CT-CT CT-HC CT-OH CT-HO HC-HC HC-OH - * 36 38 40 42 - * HC-HO OH-OH OH-HO HO-HO - */ - - // Partial g(r) (unbound terms - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-OW//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 1.0)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-HW//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.5)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-CT//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 0.2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-HC//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 8}, 7.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-OH//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 10}, 0.2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-HO//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 12}, 0.3)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-HW//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 14}, 0.4)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-CT//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 16}, 0.1)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-HC//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 18}, 4.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-OH//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 20}, 0.2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-HO//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 22}, 0.2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-CT//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 24}, 0.2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-HC//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 26}, 4.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-OH//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 28}, 0.1)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-HO//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 30}, 0.1)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HC-HC//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 32}, 4.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HC-OH//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 34}, 4.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HC-HO//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 36}, 5.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OH-OH//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 38}, 0.3)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OH-HO//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 40}, 0.1)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HO-HO//Unbound", - {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 42}, 0.3)); - - // Partial g(r) (intramolecular terms) - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-HW//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.8)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-HW//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 14}, 0.5)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-HC//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 26}, 0.3)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-OH//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 28}, 0.5)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-HO//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 30}, 0.2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HC-HC//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 32}, 0.06)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HC-OH//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 34}, 0.08)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HC-HO//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 36}, 0.5)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OH-HO//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 40}, 0.5)); - - // Partial g(r) (intramolecular terms, zero) - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-OW//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-CT//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-HC//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 8}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-OH//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 10}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OW-HO//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 12}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-CT//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 16}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-HC//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 18}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-OH//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 20}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HW-HO//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 22}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//CT-CT//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 24}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//OH-OH//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 38}, 1.0e-5, - Error::ErrorType::RMSEError)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Mix//OriginalGR//HO-HO//Bound", - {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 42}, 1.0e-5, - Error::ErrorType::RMSEError)); -} - -TEST_F(GRModuleTest, Benzene) -{ - ASSERT_NO_THROW_VERBOSE(systemTest.setUp("dissolve/input/epsr-benzene-3n.txt")); - ASSERT_TRUE(systemTest.dissolve().iterate(1)); - - // Partial g(r) (unbound terms) - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Liquid//OriginalGR//CA-CA//Unbound", - {"epsr25/benzene200-neutron/benzene.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 3.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Liquid//OriginalGR//CA-HA//Unbound", - {"epsr25/benzene200-neutron/benzene.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 2.0e-2)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Liquid//OriginalGR//HA-HA//Unbound", - {"epsr25/benzene200-neutron/benzene.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 4.0e-2)); - - // Partial g(r) (intramolecular terms) - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Liquid//OriginalGR//CA-CA//Bound", - {"epsr25/benzene200-neutron/benzene.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 0.12)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Liquid//OriginalGR//CA-HA//Bound", - {"epsr25/benzene200-neutron/benzene.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.18)); - EXPECT_TRUE(systemTest.checkData1D( - "GRs//Liquid//OriginalGR//HA-HA//Bound", - {"epsr25/benzene200-neutron/benzene.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 9.0e-2)); -} - -} // namespace UnitTest diff --git a/tests/nodes/CMakeLists.txt b/tests/nodes/CMakeLists.txt index 263f037f2e..4c9889d23c 100644 --- a/tests/nodes/CMakeLists.txt +++ b/tests/nodes/CMakeLists.txt @@ -1,5 +1,6 @@ dissolve_add_test(SRC atomicMC.cpp) dissolve_add_test(SRC flow.cpp) +dissolve_add_test(SRC gr.cpp) dissolve_add_test(SRC graph.cpp) dissolve_add_test(SRC graph_argon.cpp) dissolve_add_test(SRC loop.cpp) diff --git a/tests/nodes/atomicMC.cpp b/tests/nodes/atomicMC.cpp index 16b2231d92..4ce43a527f 100644 --- a/tests/nodes/atomicMC.cpp +++ b/tests/nodes/atomicMC.cpp @@ -12,6 +12,7 @@ #include "nodes/iterableGraph.h" #include "nodes/numberNode.h" #include "nodes/species.h" +#include "tests/speciesData.h" #include "tests/testData.h" #include #include diff --git a/tests/nodes/gr.cpp b/tests/nodes/gr.cpp new file mode 100644 index 0000000000..d91134d85f --- /dev/null +++ b/tests/nodes/gr.cpp @@ -0,0 +1,301 @@ +// SPDX-License-Identifier: GPL-3.0-or-later +// Copyright (c) 2026 Team Dissolve and contributors + +#include "nodes/gr/gr.h" +#include "math/windowFunction.h" +#include "tests/graphData.h" +#include "tests/testData.h" +#include + +namespace UnitTest +{ +TEST(GRNodeTest, Methods) +{ + GraphTestData data; + createArgonGraph(&data.graphRoot, 5000); + + // Add GR node and link in configuration + auto grNode = data.graphRoot.createNode("GR", "GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(grNode->setOption("Averaging", std::optional())); + ASSERT_TRUE(grNode->setOption("IntraBroadening", Function1DWrapper())); + ASSERT_TRUE(data.graphRoot.addEdge({"Insert", "Configuration", "GR", "Configuration"})); + + // Calculate baseline GR with the "Test" method, a simple, serial double-loop + ASSERT_TRUE(grNode->setOption("Method", GRNode::PartialsMethod::TestMethod)); + ASSERT_EQ(grNode->run(), NodeConstants::ProcessResult::Success); + ASSERT_EQ(grNode->versionIndex(), 0); + auto rawGRBaseline = *grNode->getOutputValue("RawGR"); + + // Test against simple method + ASSERT_TRUE(grNode->setOption("Method", GRNode::PartialsMethod::SimpleMethod)); + ASSERT_EQ(grNode->run(), NodeConstants::ProcessResult::Success); + auto rawGRSimple = *grNode->getOutputValue("RawGR"); + ASSERT_EQ(grNode->versionIndex(), 1); + ASSERT_TRUE(DissolveSystemTest::checkPartialSet(rawGRBaseline, rawGRSimple, 1.0e-8)); + + // Test against cells method + ASSERT_TRUE(grNode->setOption("Method", GRNode::PartialsMethod::CellsMethod)); + ASSERT_EQ(grNode->run(), NodeConstants::ProcessResult::Success); + auto rawGRCells = *grNode->getOutputValue("RawGR"); + ASSERT_EQ(grNode->versionIndex(), 2); + ASSERT_TRUE(DissolveSystemTest::checkPartialSet(rawGRBaseline, rawGRCells, 1.0e-8)); +} + +TEST(GRNodeTest, Water) +{ + GraphTestData data; + createWater1000Graph(&data.graphRoot, CoordinateImportFileFormat("epsr25/water1000-neutron/waterbox.ato", + CoordinateImportFileFormat::CoordinateImportFormat::EPSR)); + + // Set GR options + auto grNode = data.graphRoot.findNode("GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(grNode->setOption("IntraBroadening", Function1DWrapper())); + ASSERT_TRUE(grNode->setOption("BinWidth", 0.03)); + + // Run the graph from the GR node + ASSERT_EQ(grNode->run(), NodeConstants::ProcessResult::Success); + ASSERT_EQ(grNode->versionIndex(), 0); + + // Get the raw GR + auto rawGR = grNode->getOutputValue("RawGR"); + + // Partial g(r) (unbound terms) + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "OW")), "OW-OW Unbound Partial", + {"epsr25/water1000-neutron/water.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 6.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "OW")), "HW-OW Unbound Partial", + {"epsr25/water1000-neutron/water.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 2.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "HW")), "HW-HW Unbound Partial", + {"epsr25/water1000-neutron/water.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 2.0e-2)); + + // Partial g(r) (intramolecular terms) + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "OW")), "OW-OW Bound Partial", + {"epsr25/water1000-neutron/water.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "OW")), "HW-OW Bound Partial", + {"epsr25/water1000-neutron/water.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.1)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "HW")), "HW-HW Bound Partial", + {"epsr25/water1000-neutron/water.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 1.5e-2)); +} + +TEST(GRNodeTest, WaterMethanol) +{ + GraphTestData data; + createWaterMethanolGraph(&data.graphRoot); + + // Set GR options + auto grNode = data.graphRoot.findNode("GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(grNode->setOption("IntraBroadening", Function1DWrapper())); + ASSERT_TRUE(grNode->setOption("BinWidth", 0.03)); + + // Run the graph from the GR node + ASSERT_EQ(grNode->run(), NodeConstants::ProcessResult::Success); + ASSERT_EQ(grNode->versionIndex(), 0); + + // Get the raw GR + auto rawGR = grNode->getOutputValue("RawGR"); + + /* + * Partial Radial Distribution Functions + * Order of partials in EPSR files is: + * 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 + * OW-OW OW-HW OW-CT OW-HC OW-OH OW-HO HW-HW HW-CT HW-HC HW-OH HW-HO CT-CT CT-HC CT-OH CT-HO HC-HC HC-OH + * 36 38 40 42 + * HC-HO OH-OH OH-HO HO-HO + */ + + // Partial g(r) (unbound terms + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "OW")), "OW-OW Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 1.0)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "HW")), "OW-HW Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.5)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "CT")), "OW-CT Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 0.2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "HC")), "OW-HC Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 8}, 7.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "OH")), "OW-OH Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 10}, 0.2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OW", "HO")), "OW-HO Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 12}, 0.3)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "HW")), "HW-HW Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 14}, 0.4)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "CT")), "HW-CT Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 16}, 0.1)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "HC")), "HW-HC Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 18}, 4.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "OH")), "HW-OH Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 20}, 0.2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HW", "HO")), "HW-HO Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 22}, 0.2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("CT", "CT")), "CT-CT Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 24}, 0.2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("CT", "HC")), "CT-HC Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 26}, 4.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("CT", "OH")), "CT-OH Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 28}, 0.1)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("CT", "HO")), "CT-HO Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 30}, 0.1)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HC", "HC")), "HC-HC Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 32}, 4.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HC", "OH")), "HC-OH Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 34}, 4.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HC", "HO")), "HC-HO Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 36}, 5.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OH", "OH")), "OH-OH Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 38}, 0.3)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("OH", "HO")), "OH-HO Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 40}, 0.1)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HO", "HO")), "HO-HO Unbound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 42}, 0.3)); + + // Partial g(r) (intramolecular terms) + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "HW")), "OW-HW Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.8)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "HW")), "HW-HW Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 14}, 0.5)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("CT", "HC")), "CT-HC Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 26}, 0.3)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("CT", "OH")), "CT-OH Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 28}, 0.5)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("CT", "HO")), "CT-HO Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 30}, 0.2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HC", "HC")), "HC-HC Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 32}, 0.06)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HC", "OH")), "HC-OH Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 34}, 0.08)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HC", "HO")), "HC-HO Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 36}, 0.5)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OH", "HO")), "OH-HO Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 40}, 0.5)); + + // Partial g(r) (intramolecular terms, zero) + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "OW")), "OW-OW Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "CT")), "OW-CT Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "HC")), "OW-HC Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 8}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "OH")), "OW-OH Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 10}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OW", "HO")), "OW-HO Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 12}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "CT")), "HW-CT Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 16}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "HC")), "HW-HC Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 18}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "OH")), "HW-OH Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 20}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HW", "HO")), "HW-HO Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 22}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("CT", "CT")), "CT-CT Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 24}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("OH", "OH")), "OH-OH Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 38}, 1.0e-5, + Error::ErrorType::RMSEError)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HO", "HO")), "HO-HO Bound Partial", + {"epsr25/water300methanol600/watermeth.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 42}, 1.0e-5, + Error::ErrorType::RMSEError)); +} + +TEST(GRNodeTest, Benzene) +{ + GraphTestData data; + createBenzeneGraph(&data.graphRoot); + + // Set GR options + auto grNode = data.graphRoot.findNode("GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(grNode->setOption("IntraBroadening", Function1DWrapper())); + ASSERT_TRUE(grNode->setOption("BinWidth", 0.03)); + + // Run the graph from the GR node + ASSERT_EQ(grNode->run(), NodeConstants::ProcessResult::Success); + ASSERT_EQ(grNode->versionIndex(), 0); + + // Get the raw GR + auto rawGR = grNode->getOutputValue("RawGR"); + + // Partial g(r) (unbound terms) + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("CA", "CA")), "CA-CA Unbound Partial", + {"epsr25/benzene200-neutron/benzene.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 3.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("CA", "HA")), "CA-HA Unbound Partial", + {"epsr25/benzene200-neutron/benzene.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 2.0e-2)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->unboundPartials().get(DoubleKeyedMapKey("HA", "HA")), "HA-HA Unbound Partial", + {"epsr25/benzene200-neutron/benzene.EPSR.g01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 4.0e-2)); + + // Partial g(r) (intramolecular terms) + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("CA", "CA")), "CA-CA Bound Partial", + {"epsr25/benzene200-neutron/benzene.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 2}, 0.12)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("CA", "HA")), "CA-HA Bound Partial", + {"epsr25/benzene200-neutron/benzene.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 4}, 0.18)); + EXPECT_TRUE(DissolveSystemTest::checkData1D( + rawGR->boundPartials().get(DoubleKeyedMapKey("HA", "HA")), "HA-HA Bound Partial", + {"epsr25/benzene200-neutron/benzene.EPSR.y01", Data1DImportFileFormat::Data1DImportFormat::XY, 1, 6}, 9.0e-2)); +} + +} // namespace UnitTest \ No newline at end of file diff --git a/tests/nodes/graph_argon.cpp b/tests/nodes/graph_argon.cpp index bf3c126bcd..444ab71108 100644 --- a/tests/nodes/graph_argon.cpp +++ b/tests/nodes/graph_argon.cpp @@ -3,172 +3,86 @@ #include "base/units.h" #include "data/structureFactors.h" -#include "io/import/coordinates.h" #include "io/import/data1D.h" -#include "nodes/atomicMC/atomicMC.h" -#include "nodes/atomicSpecies.h" -#include "nodes/configuration.h" -#include "nodes/data1DImport.h" -#include "nodes/dissolve.h" -#include "nodes/energy/energy.h" -#include "nodes/gr/gr.h" -#include "nodes/importConfigurationCoordinates.h" -#include "nodes/insert.h" -#include "nodes/md/md.h" -#include "nodes/neutronSQ/neutronSQ.h" -#include "nodes/number.h" -#include "nodes/sq/sq.h" +#include "tests/graphData.h" #include "tests/testData.h" #include namespace UnitTest { -class GraphArgonTest : public ::testing::Test +TEST(GraphArgonTest, InitSimulation) { - public: - GraphArgonTest() : dissolve_(coreData_), root_(dissolve_) { Node::echo_ = true; } + GraphTestData data; + createArgonGraph(&data.graphRoot); - // Create a graph for testing - void createGraph(bool advanced = false) - { - /* - * Configuration (Bulk) - * ------------------ - * - Configuration-o ---+ - * - | \ Insert (Insert) - * -----------------/ \ ------------------ - * +---- o-Configuration | - * AtomicSpecies (Ar) +------ o-Species | - * ------------------ / -----------------/ - * - Species-o --+ - * - | - * -----------------/ - */ - - // Create nodes - arNode_ = - dynamic_cast(root_.addNode(std::make_unique(&root_, Elements::Ar), "Ar")); - configurationNode_ = dynamic_cast(root_.createNode("Configuration", "Bulk")); - importConfigCoordsNode_ = - dynamic_cast(root_.createNode("ImportConfigurationCoordinates", "BulkXYZ")); - insertNode_ = dynamic_cast(root_.createNode("Insert", "Insert")); - - ASSERT_TRUE(arNode_); - ASSERT_TRUE(configurationNode_); - ASSERT_TRUE(importConfigCoordsNode_); - ASSERT_TRUE(insertNode_); - - ASSERT_TRUE(root_.addEdge({"Ar", "Species", "Insert", "Species"})); - - /* - * Set up configuration XYZ - */ - ASSERT_TRUE(root_.addEdge({"Bulk", "Configuration", "Insert", "Configuration"})); - ASSERT_TRUE(importConfigCoordsNode_->setOption("FilePath", "dissolve2/argon/Ar_bulk_step1000.xyz")); - ASSERT_TRUE(importConfigCoordsNode_->setOption( - "FileFormat", CoordinateImportFileFormat::CoordinateImportFormat::XYZ)); - ASSERT_TRUE(root_.addEdge({"Insert", "Configuration", "BulkXYZ", "Configuration"})); - - /* - * Set Insert node options - */ - auto population = insertNode_->findInput("Population"); - population->set(1000); - - auto density = insertNode_->findInput("Density"); - density->set(0.0213); - - ASSERT_TRUE(insertNode_->setOption("DensityUnits", Units::DensityUnits::AtomsPerAngstromUnits)); - - if (advanced) - { - grNode_ = dynamic_cast(root_.createNode("GR", "GR")); - sqNode_ = dynamic_cast(root_.createNode("SQ", "SQ")); - neutronSQNode_ = dynamic_cast(root_.createNode("NeutronSQ", "NeutronSQ")); - data1DImportNode_ = dynamic_cast(root_.createNode("Data1DImport", "ReferenceSQ")); - - ASSERT_TRUE(grNode_); - ASSERT_TRUE(sqNode_); - ASSERT_TRUE(neutronSQNode_); - ASSERT_TRUE(data1DImportNode_); - - /* - * Set up GR options - */ - ASSERT_TRUE(grNode_->setOption("BinWidth", 0.025)); - - /* - * Set up reference SQ data - */ - ASSERT_TRUE(data1DImportNode_->setOption("FilePath", "dissolve2/argon/yarnell.sq")); - ASSERT_TRUE(data1DImportNode_->setOption( - "ImportFormat", Data1DImportFileFormat::Data1DImportFormat::XY)); - ASSERT_TRUE(data1DImportNode_->setOption>("RemoveAverageFromX", 9.0)); - - /* - * Set up neutron SQ options - */ - ASSERT_TRUE(neutronSQNode_->setOption( - "ReferenceNormalisedTo", StructureFactors::SquareOfAverageNormalisation)); - - // Create nodes - ASSERT_TRUE(root_.addEdge({"BulkXYZ", "Configuration", "GR", "Configuration"})); - ASSERT_TRUE(root_.addEdge({"GR", "UnweightedGR", "SQ", "UnweightedGR"})); - ASSERT_TRUE(root_.addEdge({"ReferenceSQ", "Data", "NeutronSQ", "ReferenceData"})); - ASSERT_TRUE(root_.addEdge({"SQ", "UnweightedGR", "NeutronSQ", "UnweightedGR"})); - ASSERT_TRUE(root_.addEdge({"SQ", "UnweightedSQ", "NeutronSQ", "UnweightedSQ"})); - } - } - - protected: - // We need a CoreData and Dissolve definition to properly instantiate DissolveGraph at present. - CoreData coreData_; - Dissolve dissolve_; - DissolveGraph root_; - AtomicSpeciesNode *arNode_{nullptr}; - ConfigurationNode *configurationNode_{nullptr}; - ImportConfigurationCoordinatesNode *importConfigCoordsNode_{nullptr}; - InsertNode *insertNode_{nullptr}; - AtomicMCNode *atomicMCNode_{nullptr}; - MDNode *mdNode_{nullptr}; - EnergyNode *energyNode_{nullptr}; - GRNode *grNode_{nullptr}; - SQNode *sqNode_{nullptr}; - NeutronSQNode *neutronSQNode_{nullptr}; - Data1DImportNode *data1DImportNode_{nullptr}; -}; - -TEST_F(GraphArgonTest, InitSimulation) -{ - createGraph(); - - ASSERT_EQ(insertNode_->run(), NodeConstants::ProcessResult::Success); - auto *cfg = insertNode_->getOutputValue("Configuration"); + // Get the Insert node and run the graph + auto insertNode = data.graphRoot.findNode("Insert"); + ASSERT_TRUE(insertNode); + ASSERT_EQ(insertNode->run(), NodeConstants::ProcessResult::Success); // Check Configuration contents + auto *cfg = insertNode->getOutputValue("Configuration"); ASSERT_TRUE(cfg); - EXPECT_EQ(cfg->nMolecules(), insertNode_->getInputValue("Population").asInteger()); + EXPECT_EQ(cfg->nMolecules(), insertNode->getInputValue("Population").asInteger()); }; -TEST_F(GraphArgonTest, AdvancedSimulation) +TEST(GraphArgonTest, AllCorrelations) { - createGraph(true); - ASSERT_EQ(neutronSQNode_->run(), NodeConstants::ProcessResult::Success); + GraphTestData data; + createArgonGraph(&data.graphRoot); + + // Add in imported coordinates + auto importConfigCoordsNode = data.graphRoot.createNode("ImportConfigurationCoordinates", "BulkXYZ"); + ASSERT_TRUE(importConfigCoordsNode); + ASSERT_TRUE(importConfigCoordsNode->setOption("FilePath", "dissolve2/argon/Ar_bulk_step1000.xyz")); + ASSERT_TRUE(importConfigCoordsNode->setOption( + "FileFormat", CoordinateImportFileFormat::CoordinateImportFormat::XYZ)); + ASSERT_TRUE(data.graphRoot.addEdge({"Insert", "Configuration", "BulkXYZ", "Configuration"})); + + // Add GR node and link in configuration + auto grNode = data.graphRoot.createNode("GR", "GR"); + ASSERT_TRUE(grNode); + ASSERT_TRUE(grNode->setOption("BinWidth", 0.025)); + ASSERT_TRUE(data.graphRoot.addEdge({"BulkXYZ", "Configuration", "GR", "Configuration"})); + + // Add SQ node and link in GR + auto sqNode = data.graphRoot.createNode("SQ", "SQ"); + ASSERT_TRUE(sqNode); + ASSERT_TRUE(data.graphRoot.addEdge({"GR", "UnweightedGR", "SQ", "UnweightedGR"})); + + // Add NeutronSQ node and link in data + auto neutronSQNode = data.graphRoot.createNode("NeutronSQ", "NeutronSQ"); + ASSERT_TRUE(neutronSQNode); + ASSERT_TRUE(neutronSQNode->setOption("ReferenceNormalisedTo", + StructureFactors::SquareOfAverageNormalisation)); + ASSERT_TRUE(data.graphRoot.addEdge({"SQ", "UnweightedGR", "NeutronSQ", "UnweightedGR"})); + ASSERT_TRUE(data.graphRoot.addEdge({"SQ", "UnweightedSQ", "NeutronSQ", "UnweightedSQ"})); + + // Set reference F(Q) data + auto data1DImportNode = data.graphRoot.createNode("Data1DImport", "Yarnell"); + ASSERT_TRUE(data1DImportNode); + ASSERT_TRUE(data1DImportNode->setOption("FilePath", "dissolve2/argon/yarnell.sq")); + ASSERT_TRUE(data1DImportNode->setOption( + "ImportFormat", Data1DImportFileFormat::Data1DImportFormat::XY)); + ASSERT_TRUE(data1DImportNode->setOption>("RemoveAverageFromX", 9.0)); + ASSERT_TRUE(data.graphRoot.addEdge({"Yarnell", "Data", "NeutronSQ", "ReferenceData"})); /* - * Check total unweighted SQ + * Run the Graph */ - auto unweightedSQ = sqNode_->getOutputValue("UnweightedSQ"); - ASSERT_NO_THROW_VERBOSE(unweightedSQ); + + ASSERT_EQ(neutronSQNode->run(), NodeConstants::ProcessResult::Success); + + // Check total unweighted SQ + auto unweightedSQ = sqNode->getOutputValue("UnweightedSQ"); + ASSERT_TRUE(unweightedSQ); ASSERT_TRUE( DissolveSystemTest::checkData1D(unweightedSQ->total(), "UnweightedSQ", {"dissolve2/argon/SQ01-UnweightedSQ-total.sq"})); - /* - * Check neutron weighted SQ - */ + // Check neutron weighted SQ const auto tolerance = 0.025; - auto weightedSQ = neutronSQNode_->getOutputValue("WeightedSQ"); - ASSERT_NO_THROW_VERBOSE(weightedSQ); + auto weightedSQ = neutronSQNode->getOutputValue("WeightedSQ"); + ASSERT_TRUE(weightedSQ); ASSERT_TRUE(DissolveSystemTest::checkData1D(weightedSQ->total(), "WeightedSQ", {"dissolve2/argon/NeutronSQ01-WeightedSQ-total.sq"}, tolerance)); } diff --git a/tests/nodes/parameters.cpp b/tests/nodes/parameters.cpp index 57c2a60c7c..63e01e6964 100644 --- a/tests/nodes/parameters.cpp +++ b/tests/nodes/parameters.cpp @@ -4,28 +4,27 @@ #include "nodes/dissolve.h" #include "nodes/numberNode.h" #include "nodes/test.h" +#include "tests/graphData.h" #include namespace UnitTest { TEST(ParametersTest, OptionalPointerOutput) { - CoreData coreData_; - Dissolve dissolve_(coreData_); - DissolveGraph root_(dissolve_); + GraphTestData data; // Create a couple of TestNodes - auto *a = dynamic_cast(root_.addNode(std::make_unique(&root_), "TestA")); + auto *a = dynamic_cast(data.graphRoot.addNode(std::make_unique(&data.graphRoot), "TestA")); ASSERT_TRUE(a); auto createOptA = a->findInput("CreateConfiguration"); ASSERT_TRUE(createOptA); - auto *b = dynamic_cast(root_.addNode(std::make_unique(&root_), "TestB")); + auto *b = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "TestB"); ASSERT_TRUE(b); auto configInputB = b->findInput("ConfigurationInput"); ASSERT_TRUE(configInputB); // Create an edge between nodes - ASSERT_TRUE(root_.addEdge({"TestA", "OptionalConfiguration", "TestB", "ConfigurationInput"})); + ASSERT_TRUE(data.graphRoot.addEdge({"TestA", "OptionalConfiguration", "TestB", "ConfigurationInput"})); // Inputs to TestB should be valid (there is no optional data yet, but the Edge definitions are correct) EXPECT_TRUE(b->inputsAreValid()); @@ -45,12 +44,10 @@ TEST(ParametersTest, OptionalPointerOutput) TEST(ParametersTest, VectorParameter) { - CoreData coreData_; - Dissolve dissolve_(coreData_); - DissolveGraph root_(dissolve_); + GraphTestData data; // Create a couple of TestNodes - auto *a = dynamic_cast(root_.addNode(std::make_unique(&root_), "TestA")); + auto *a = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "TestA"); ASSERT_TRUE(a); auto numbersABase = a->findInput("NumberVector"); ASSERT_TRUE(numbersABase); @@ -64,7 +61,7 @@ TEST(ParametersTest, VectorParameter) EXPECT_ANY_THROW(numbersABase->set(Number{1.0})); // Create a Number node - auto n1 = dynamic_cast(root_.addNode(std::make_unique(&root_), "Number1")); + auto n1 = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "Number1"); ASSERT_TRUE(n1); auto number1 = n1->findOption("X"); number1->set(Number{1.0}); @@ -86,31 +83,31 @@ TEST(ParametersTest, VectorInputOutput) { CoreData coreData_; Dissolve dissolve_(coreData_); - DissolveGraph root_(dissolve_); + GraphTestData data; // Create a couple of TestNodes - auto *a = dynamic_cast(root_.addNode(std::make_unique(&root_), "TestA")); + auto *a = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "TestA"); ASSERT_TRUE(a); auto numbersA = a->findInput("NumberVector"); ASSERT_TRUE(numbersA); - auto *b = dynamic_cast(root_.addNode(std::make_unique(&root_), "TestB")); + auto *b = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "TestB"); ASSERT_TRUE(b); auto numbersB = b->findInput("NumberVector"); ASSERT_TRUE(numbersB); // Create an edge linking the vector output from A to the vector input of B - ASSERT_TRUE(root_.addEdge({"TestA", "NumberVector", "TestB", "NumberVector"})); + ASSERT_TRUE(data.graphRoot.addEdge({"TestA", "NumberVector", "TestB", "NumberVector"})); // Create three Number nodes as inputs for TestA's number vector - auto *n1 = dynamic_cast(root_.addNode(std::make_unique(&root_), "Number1")); + auto *n1 = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "Number1"); ASSERT_TRUE(n1); auto number1 = n1->findOption("X"); ASSERT_TRUE(number1); - auto *n2 = dynamic_cast(root_.addNode(std::make_unique(&root_), "Number2")); + auto *n2 = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "Number2"); ASSERT_TRUE(n2); auto number2 = n2->findOption("X"); ASSERT_TRUE(number2); - auto *n3 = dynamic_cast(root_.addNode(std::make_unique(&root_), "Number3")); + auto *n3 = data.graphRoot.addNode(std::make_unique(&data.graphRoot), "Number3"); ASSERT_TRUE(n3); auto number3 = n3->findOption("X"); ASSERT_TRUE(number3); @@ -121,9 +118,9 @@ TEST(ParametersTest, VectorInputOutput) number3->set(Number{8.0}); // Link all three numbers in to the TestA vector - ASSERT_TRUE(root_.addEdge({"Number1", "X", "TestA", "NumberVector"})); - ASSERT_TRUE(root_.addEdge({"Number2", "X", "TestA", "NumberVector"})); - ASSERT_TRUE(root_.addEdge({"Number3", "X", "TestA", "NumberVector"})); + ASSERT_TRUE(data.graphRoot.addEdge({"Number1", "X", "TestA", "NumberVector"})); + ASSERT_TRUE(data.graphRoot.addEdge({"Number2", "X", "TestA", "NumberVector"})); + ASSERT_TRUE(data.graphRoot.addEdge({"Number3", "X", "TestA", "NumberVector"})); // Run the TestB node to pull the number edge vector from TestA, using the numbers from the three number nodes EXPECT_TRUE(a->inputsAreValid()); @@ -157,7 +154,7 @@ TEST(ParametersTest, VectorInputOutput) EXPECT_EQ(numbersA->get>(), numbersB->get>()); // Remove a single edge - this should flag TestA and TestB as being out of date - EXPECT_TRUE(root_.removeEdge({"Number1", "X", "TestA", "NumberVector"})); + EXPECT_TRUE(data.graphRoot.removeEdge({"Number1", "X", "TestA", "NumberVector"})); EXPECT_TRUE(a->inputsAreValid()); EXPECT_TRUE(b->inputsAreValid()); EXPECT_FALSE(a->isUpToDate()); @@ -172,8 +169,8 @@ TEST(ParametersTest, VectorInputOutput) EXPECT_EQ(numbersA->get>(), numbersB->get>()); // Remove both of the other edges - EXPECT_TRUE(root_.removeEdge({"Number3", "X", "TestA", "NumberVector"})); - EXPECT_TRUE(root_.removeEdge({"Number2", "X", "TestA", "NumberVector"})); + EXPECT_TRUE(data.graphRoot.removeEdge({"Number3", "X", "TestA", "NumberVector"})); + EXPECT_TRUE(data.graphRoot.removeEdge({"Number2", "X", "TestA", "NumberVector"})); EXPECT_TRUE(a->inputsAreValid()); EXPECT_TRUE(b->inputsAreValid()); EXPECT_FALSE(a->isUpToDate()); diff --git a/tests/speciesData.h b/tests/speciesData.h new file mode 100644 index 0000000000..8a6605bacd --- /dev/null +++ b/tests/speciesData.h @@ -0,0 +1,185 @@ +// SPDX-License-Identifier: GPL-3.0-or-later +// Copyright (c) 2026 Team Dissolve and contributors + +#pragma once + +#include "classes/species.h" +#include "data/elements.h" +#include "math/mathFunc.h" +#include "nodes/graph.h" +#include "nodes/species.h" +#include + +namespace UnitTest +{ +// Create and return water test species in the specified graph +inline SpeciesNode *createWater(Graph *parentGraph) +{ + const auto name = "Water"; + + // Add species node + auto speciesNodeUniquePtr = std::make_unique(parentGraph); + auto speciesNodePtr = speciesNodeUniquePtr.get(); + auto species = &(speciesNodePtr->species()); + parentGraph->addNode(std::move(speciesNodeUniquePtr), name); + species->setName(name); + + // Set up atom types + auto oW = species->addAtomType(Elements::Element::O, "OW"); + oW->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJones, "epsilon=0.6503 sigma=3.165492"); + oW->setCharge(-0.82); + species->addAtom(Elements::Element::O, {}, -0.82, oW); + auto hW = species->addAtomType(Elements::Element::H, "HW"); + hW->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJones, "epsilon=0.0 sigma=0.0"); + hW->setCharge(0.41); + species->addAtom(Elements::Element::H, {1, 0, 0}, 0.41, hW); + species->addAtom(Elements::Element::H, {cos(DissolveMath::toRadians(113.24)), sin(DissolveMath::toRadians(113.24)), 0.0}, + 0.41, hW); + + // Apply intramolecular terms + species->addBond(0, 1).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=4431.53 eq=1.0"); + species->addBond(0, 2).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=4431.53 eq=1.0"); + species->addAngle(1, 0, 2).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=317.5656 eq=113.24"); + + // Create isotopologue + auto iso = species->addIsotopologue("D2O"); + iso->setAtomTypeIsotope(hW.get(), Sears91::H_2); + + return speciesNodePtr; +} + +// Create and return methanol test species in the specified graph +inline SpeciesNode *createMethanol(Graph *parentGraph) +{ + const auto name = "Methanol"; + + // Add methanol species node + auto speciesNodeUniquePtr = std::make_unique(parentGraph); + auto speciesNodePtr = speciesNodeUniquePtr.get(); + auto species = &(speciesNodePtr->species()); + parentGraph->addNode(std::move(speciesNodeUniquePtr), name); + species->setName(name); + + // Create atom types + auto CT = species->addAtomType(Elements::Element::C, "CT"); + CT->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJonesGeometric, + "epsilon=0.276 sigma=3.55"); + CT->setCharge(-0.18); + auto HC = species->addAtomType(Elements::Element::H, "HC"); + HC->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJonesGeometric, + "epsilon=0.126 sigma=2.5"); + HC->setCharge(0.06); + auto OH = species->addAtomType(Elements::Element::O, "OH"); + OH->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJonesGeometric, + "epsilon=0.711 sigma=3.12"); + OH->setCharge(-0.68); + auto HO = species->addAtomType(Elements::Element::H, "HO"); + HO->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJonesGeometric, + "epsilon=0.126 sigma=2.4"); + HO->setCharge(0.68); + + // Add atoms + species->addAtom(Elements::Element::C, {0.0, 0.0, 0.0}, -0.18, CT); + species->addAtom(Elements::Element::H, {1.1187, 0.0, 0.0}, 0.06, HC); + species->addAtom(Elements::Element::O, {-0.3683, 1.3617, 0.0}, -0.68, OH); + species->addAtom(Elements::Element::H, {-0.3834, -0.5181, -0.9144}, 0.06, HC); + species->addAtom(Elements::Element::H, {-0.3834, -0.5177, 0.9146}, 0.06, HC); + species->addAtom(Elements::Element::H, {-1.3318, 1.3955, -0.17}, 0.68, HO); + + // Apply intramolecular terms + species->addBond(0, 1).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3000.0 eq=1.12"); + species->addBond(0, 2).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3000.0 eq=1.41"); + species->addBond(0, 3).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3000.0 eq=1.12"); + species->addBond(0, 4).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3000.0 eq=1.12"); + species->addBond(2, 5).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3000.0 eq=0.964"); + species->addAngle(1, 0, 2).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + species->addAngle(1, 0, 3).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + species->addAngle(1, 0, 4).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + species->addAngle(2, 0, 3).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + species->addAngle(2, 0, 4).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + species->addAngle(3, 0, 4).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + species->addAngle(5, 2, 0).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=300.0 eq=109.5"); + + // Create isotopologues + auto D = species->addIsotopologue("Deuteriated"); + D->setAtomTypeIsotope(HO.get(), Sears91::H_2); + D->setAtomTypeIsotope(HC.get(), Sears91::H_2); + auto MeD = species->addIsotopologue("MethylD-OH"); + D->setAtomTypeIsotope(HC.get(), Sears91::H_2); + auto MeH = species->addIsotopologue("MethylH-OD"); + D->setAtomTypeIsotope(HO.get(), Sears91::H_2); + + return speciesNodePtr; +} + +// Create and return benzene test species in the specified graph +inline SpeciesNode *createBenzene(Graph *parentGraph) +{ + const auto name = "Benzene"; + + // Add species node + auto speciesNodeUniquePtr = std::make_unique(parentGraph); + auto speciesNodePtr = speciesNodeUniquePtr.get(); + auto species = &(speciesNodePtr->species()); + parentGraph->addNode(std::move(speciesNodeUniquePtr), name); + species->setName(name); + + // Set up atom types + auto CA = species->addAtomType(Elements::Element::C, "CA"); + CA->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJonesGeometric, + "epsilon=0.29288 sigma=3.55"); + CA->setCharge(-0.115); + auto HA = species->addAtomType(Elements::Element::H, "HA"); + HA->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJonesGeometric, + "epsilon=0.12552 sigma=2.42"); + HA->setCharge(0.115); + + // Add atoms + species->addAtom(Elements::Element::C, {-1.203775, 0.695, 0.0}, -0.115, CA); + species->addAtom(Elements::Element::H, {-2.069801, 1.195, 0.0}, -0.115, HA); + species->addAtom(Elements::Element::C, {-0.000000, 1.390, 0.0}, -0.115, CA); + species->addAtom(Elements::Element::H, {-0.000000, 2.390, 0.0}, -0.115, HA); + species->addAtom(Elements::Element::C, {1.203775, 0.695, 0.0}, -0.115, CA); + species->addAtom(Elements::Element::H, {2.069801, 1.195, 0.0}, -0.115, HA); + species->addAtom(Elements::Element::C, {1.203775, -0.695, 0.0}, -0.115, CA); + species->addAtom(Elements::Element::H, {2.069801, -1.195, 0.0}, -0.115, HA); + species->addAtom(Elements::Element::C, {-0.000000, -1.390, 0.0}, -0.115, CA); + species->addAtom(Elements::Element::H, {-0.000000, -2.390, 0.0}, -0.115, HA); + species->addAtom(Elements::Element::C, {-1.203775, -0.695, 0.0}, -0.115, CA); + species->addAtom(Elements::Element::H, {-2.069801, -1.195, 0.0}, -0.115, HA); + + // Add intramolecular terms + for (auto i = 0; i < 12; i += 2) + { + // C-H bond: (i , i+1) + species->addBond(i, i + 1).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3071.056 eq=1.08"); + // C-C bond: (i, i+2) + species->addBond(i, (i + 2) % 12).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=3924.592 eq=1.4"); + // H-C-C angles: (i+1, i, i+2) and (i+1, i, i-2) + species->addAngle(i + 1, i, (i + 2) % 12) + .setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=292.88 eq=120.0"); + species->addAngle(i + 1, i, DissolveMath::wrap(i - 2, 0, 11)) + .setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=292.88 eq=120.0"); + // C-C-C angle: (i-2, i, i+2) + species->addAngle(DissolveMath::wrap(i - 2, 0, 11), i, (i + 2) % 12) + .setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=527.184 eq=120.0"); + // H-C-C-H torsion: (i+1, i, i+2, i+3) + species->addTorsion(i + 1, i, DissolveMath::wrap(i + 2, 0, 11), DissolveMath::wrap(i + 3, 0, 11)) + .setInteractionFormAndParameters(TorsionFunctions::Form::Cos3, "k1=0.0 k2=30.334 k3=0.0"); + // H-C-C-C torsion: (i+1, i, i+2, i+4) + species->addTorsion(i + 1, i, DissolveMath::wrap(i + 2, 0, 11), DissolveMath::wrap(i + 4, 0, 11)) + .setInteractionFormAndParameters(TorsionFunctions::Form::Cos3, "k1=0.0 k2=30.334 k3=0.0"); + // C-C-C-C torsion: (i, i+2, i+4, i+6) + species + ->addTorsion(i, DissolveMath::wrap(i + 2, 0, 11), DissolveMath::wrap(i + 4, 0, 11), + DissolveMath::wrap(i + 6, 0, 11)) + .setInteractionFormAndParameters(TorsionFunctions::Form::Cos3, "k1=0.0 k2=30.334 k3=0.0"); + } + + // Create isotopologue + auto iso = species->addIsotopologue("C6D6"); + iso->setAtomTypeIsotope(HA.get(), Sears91::H_2); + + return speciesNodePtr; +} +} // namespace UnitTest diff --git a/tests/testData.h b/tests/testData.h index 5a27b0bc60..14a7838e73 100644 --- a/tests/testData.h +++ b/tests/testData.h @@ -493,6 +493,71 @@ class DissolveSystemTest checkIntramolecularTerms(std::format("improper {}", joinStrings(atoms, "-")), expectedParams, i->get().interactionPotential(), tolerance); } + // Test consistency between the two supplied double-keyed Data1D maps + static bool checkDoubleKeyedMap(std::string_view mapContents, const DoubleKeyedMap &mapA, + const DoubleKeyedMap &mapB, double testThreshold) + { + // Check map sizes + if (mapA.size() != mapB.size()) + { + std::cout << std::format("Maps containing {} data are of dissimilar size (A = {}, B = {})\n", mapContents, + mapA.size(), mapB.size()); + return false; + } + + // Check individual data + for (auto &[key, dataA] : mapA) + { + // Find same-keyed data in mapB + if (mapB.contains(key)) + { + auto errorReport = Error::percent(dataA, mapB.get(key)); + std::cout << Error::errorReportString(errorReport) << std::endl; + std::cout << std::format("{} '{}' in map B has {} error of {:7.3f}{} with data in map A and is " + "{} (threshold is {:6.3f}%)\n\n", + mapContents, key, Error::errorTypes().keyword(errorReport.errorType), + errorReport.error, errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", + errorReport.error <= testThreshold ? "OK" : "NOT OK", testThreshold); + if (errorReport.error > testThreshold) + return false; + } + else + { + std::cout << std::format("{} '{}' is present in map A but not in map B.\n", mapContents, key); + return false; + } + } + + return true; + } + // Test consistency, and error, between supplied partial sets + static bool checkPartialSet(const PartialSet &setA, const PartialSet &setB, double testThreshold) + { + // Full partials + if (!checkDoubleKeyedMap("Full Partials", setA.partials(), setB.partials(), testThreshold)) + return false; + + // Bound partials + if (!checkDoubleKeyedMap("Bound Partials", setA.boundPartials(), setB.boundPartials(), testThreshold)) + return false; + + // Unbound partials + if (!checkDoubleKeyedMap("Unbound Partials", setA.unboundPartials(), setB.unboundPartials(), testThreshold)) + return false; + + // Total + auto errorReport = Error::percent(setA.total(), setB.total()); + std::cout << Error::errorReportString(errorReport) << std::endl; + std::cout << std::format( + "Total in set B has {} error of {:7.3f}{} with data in set A and is {} (threshold is {:6.3f}%)\n\n", + Error::errorTypes().keyword(errorReport.errorType), errorReport.error, + errorReport.errorType == Error::ErrorType::PercentError ? "%" : "", + errorReport.error <= testThreshold ? "OK" : "NOT OK", testThreshold); + if (errorReport.error > testThreshold) + return false; + + return true; + } }; // Return argon test species @@ -705,43 +770,6 @@ const Species &benzeneSpecies() return benzene_; } -// Return Node-Graph water test species -Species *createWater(Graph *parentGraph) -{ - const auto name = "Water"; - - // Add water species node - auto speciesNode = std::make_unique(parentGraph); - auto species = &(speciesNode.get()->species()); - parentGraph->addNode(std::move(speciesNode), name); - - // Set up water species and atom types - species->clear(); - species->setName(name); - - auto oW = std::make_shared(Elements::Element::O, "OW"); - oW->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJones, "epsilon=0.6503 sigma=3.165492"); - oW->setCharge(-0.82); - species->addAtom(Elements::Element::O, {}, -0.82, oW); - auto hW = std::make_shared(Elements::Element::H, "HW"); - hW->interactionPotential().setFormAndParameters(ShortRangeFunctions::Form::LennardJones, "epsilon=0.0 sigma=0.0"); - hW->setCharge(0.41); - species->addAtom(Elements::Element::H, {1, 0, 0}, 0.41, hW); - species->addAtom(Elements::Element::H, {cos(DissolveMath::toRadians(113.24)), sin(DissolveMath::toRadians(113.24)), 0.0}, - 0.41, hW); - - assert(species->atom(0).Z() == Elements::Element::O); - assert(species->atom(1).Z() == Elements::Element::H); - assert(species->atom(2).Z() == Elements::Element::H); - - // Apply intramolecular terms - species->addBond(0, 1).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=4431.53 eq=1.0"); - species->addBond(0, 2).setInteractionFormAndParameters(BondFunctions::Form::Harmonic, "k=4431.53 eq=1.0"); - species->addAngle(1, 0, 2).setInteractionFormAndParameters(AngleFunctions::Form::Harmonic, "k=317.5656 eq=113.24"); - - return species; -} - // Return small molecules data class SmallMolecules {