From 1748f12b32bcbd2332678e6b9a3c504c81e04044 Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Thu, 11 Dec 2025 17:26:49 +0000 Subject: [PATCH 1/7] Add ability for components to request additional components --- include/component.hxx | 23 ++++++++++++ src/component_scheduler.cxx | 19 ++++++++-- tests/unit/test_component_scheduler.cxx | 48 +++++++++++++++++++++++++ 3 files changed, 88 insertions(+), 2 deletions(-) diff --git a/include/component.hxx b/include/component.hxx index 372bf1461..c1a888898 100644 --- a/include/component.hxx +++ b/include/component.hxx @@ -76,6 +76,24 @@ private: } }; +/// Lightweight type to store information used to build a component. +struct ComponentInformation { + std::string name; + std::string type; + + ComponentInformation(const std::string& name_, const std::string& type_) : name(name_), type(type_) {} + + ComponentInformation(std::string&& name_, std::string&& type_) : name(std::move(name_)), type(std::move(type_)) {} + + bool operator<(const ComponentInformation& other) const { + return std::pair(name, type) < std::pair(other.name, other.type); + } + + bool operator==(const ComponentInformation& other) const { + return std::pair(name, type) == std::pair(other.name, other.type); + } +}; + /// Interface for a component of a simulation model /// /// The constructor of derived types should have signature @@ -92,6 +110,11 @@ struct Component { virtual ~Component() {} + /// Return a list of names/types of other components needed by this + /// component. All configurations for these components will take the + /// default value, unless set in the input file. + virtual std::vector additionalComponents() { return {}; } + /// Modify the given simulation state. This method will wrap the /// state in a GuardedOptions object and pass that to the private /// implementation of transform provided by each component. diff --git a/src/component_scheduler.cxx b/src/component_scheduler.cxx index 7df0bdcbf..1f3c7c6e8 100644 --- a/src/component_scheduler.cxx +++ b/src/component_scheduler.cxx @@ -435,6 +435,8 @@ ComponentScheduler::ComponentScheduler(Options& scheduler_options, .doc("Components in order of execution") .as(); + std::set required_components; + std::vector electrons; std::vector neutrals; std::vector positive_ions; @@ -483,9 +485,22 @@ ComponentScheduler::ComponentScheduler(Options& scheduler_options, continue; } - components.push_back( - Component::create(type_trimmed, name_trimmed, component_options, solver)); + required_components.emplace(name_trimmed, type_trimmed); + } + } + + std::set created_components; + for (auto it = required_components.begin(); it != required_components.end(); + it = required_components.begin()) { + auto comp = Component::create(it->type, it->name, component_options, solver); + for (const auto& sub_comp : comp->additionalComponents()) { + if (required_components.count(sub_comp) == 0 + and created_components.count(sub_comp) == 0) { + required_components.insert(sub_comp); + } } + components.push_back(std::move(comp)); + created_components.insert(required_components.extract(it)); } const SpeciesInformation species(electrons, neutrals, positive_ions, negative_ions); diff --git a/tests/unit/test_component_scheduler.cxx b/tests/unit/test_component_scheduler.cxx index 59c7c3437..90d34248f 100644 --- a/tests/unit/test_component_scheduler.cxx +++ b/tests/unit/test_component_scheduler.cxx @@ -29,6 +29,20 @@ struct TestMultiply : public NamedComponent { } }; +struct TestAdditionalComponent : public NamedComponent { + TestAdditionalComponent(const std::string&, Options&, Solver*) + : Component({}) {} + + void transform_impl(GuardedOptions&) override { + } + + std::vector additionalComponents() override { + return {{"TestComponent", "testcomponent"}, {"component2", "multiply"}}; + } + + static constexpr auto type = "testadditionalcomponent"; +}; + struct OrderChecker : public NamedComponent { OrderChecker(const std::string& name, Options& alloptions, Solver*) : NamedComponent(name, getPermissions(name, alloptions)) {} @@ -54,6 +68,7 @@ std::vector OrderChecker::execution_order; RegisterComponent registertestcomponent; RegisterComponent registertestcomponent2; +RegisterComponent registertestcomponent3; RegisterComponent registercomponentorderchecker; } // namespace @@ -92,6 +107,39 @@ TEST(SchedulerTest, SubComponents) { ASSERT_TRUE(options["answer"] == 42 * 2); } +TEST(SchedulerTest, AdditionalComponents) { + Options options; + options["components"] = "additionalcomponent"; + auto scheduler = ComponentScheduler::create(options, options, nullptr); + + EXPECT_FALSE(options.isSet("answer")); + scheduler->transform(options); + ASSERT_TRUE(options.isSet("answer")); + ASSERT_TRUE(options["answer"] == 42 * 2); +} + +TEST(SchedulerTest, AdditionalComponentsPredeclared) { + Options options; + options["components"] = "testcomponent, additionalcomponent"; + auto scheduler = ComponentScheduler::create(options, options, nullptr); + + EXPECT_FALSE(options.isSet("answer")); + scheduler->transform(options); + ASSERT_TRUE(options.isSet("answer")); + ASSERT_TRUE(options["answer"] == 42 * 2); +} + +TEST(SchedulerTest, AdditionalComponentsPredeclared2) { + Options options; + options["components"] = "additionalcomponent, testcomponent"; + auto scheduler = ComponentScheduler::create(options, options, nullptr); + + EXPECT_FALSE(options.isSet("answer")); + scheduler->transform(options); + ASSERT_TRUE(options.isSet("answer")); + ASSERT_TRUE(options["answer"] == 42 * 2); +} + using Parameter = std::pair>; class ComponentOrderTest : public testing::TestWithParam { From 996cf06857558e21b96efac57a231ebbce266475 Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Thu, 11 Dec 2025 17:28:14 +0000 Subject: [PATCH 2/7] Add Braginskii meta-component to activate all parts of closure --- CMakeLists.txt | 1 + hermes-3.cxx | 1 + include/braginskii_closure.hxx | 70 +++++++++++++++++++ tests/unit/test_braginskii_closure.cxx | 89 +++++++++++++++++++++++++ tests/unit/test_component_scheduler.cxx | 11 ++- 5 files changed, 166 insertions(+), 6 deletions(-) create mode 100644 include/braginskii_closure.hxx create mode 100644 tests/unit/test_braginskii_closure.cxx diff --git a/CMakeLists.txt b/CMakeLists.txt index f176f2986..c627babfe 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -125,6 +125,7 @@ set(HERMES_SOURCES include/classical_diffusion.hxx include/cx_reaction.hxx include/binormal_stpm.hxx + include/braginskii_closure.hxx include/braginskii_collisions.hxx include/braginskii_conduction.hxx include/braginskii_electron_viscosity.hxx diff --git a/hermes-3.cxx b/hermes-3.cxx index e734b03cb..406773900 100644 --- a/hermes-3.cxx +++ b/hermes-3.cxx @@ -29,6 +29,7 @@ #include "include/amjuel_data.hxx" #include "include/anomalous_diffusion.hxx" #include "include/binormal_stpm.hxx" +#include "include/braginskii_closure.hxx" #include "include/braginskii_collisions.hxx" #include "include/braginskii_conduction.hxx" #include "include/braginskii_electron_viscosity.hxx" diff --git a/include/braginskii_closure.hxx b/include/braginskii_closure.hxx new file mode 100644 index 000000000..eda9e3567 --- /dev/null +++ b/include/braginskii_closure.hxx @@ -0,0 +1,70 @@ +#pragma once +#ifndef BRAGINSKII_H +#define BRAGINSKII_H + +#include + +#include "component.hxx" + +/// Meta-component to set up all components necessary for the +/// Braginskii closure: `braginskii_collisions`, +/// `braginskii_friction`, `braginskii_heat_exchange`, +/// `braginskii_conduction`, `braginskii_electron_viscosity`, +/// `braginskii_ion_viscosity`, and `braginskii_thermal_force`. Each +/// of these components will have the same name as its type. +class BraginskiiClosure : public NamedComponent { +public: + /// @param alloptions Settings, which may include + /// - + /// - electron_viscosity : bool Include electron viscosity? (default: true) + /// - ion_viscosity : bool Include ion viscosity? (default: true) + /// - thermal_force : bool Inlucde thermal force between species? (default: true) + BraginskiiClosure(std::string name, Options& alloptions, Solver*) + : NamedComponent(name, {}) { + Options& options = alloptions[name]; + electron_viscosity = options["electron_viscosity"] + .doc("Include electron viscosity terms?") + .withDefault(true); + ion_viscosity = options["ion_viscosity"] + .doc("Include ion viscosity terms?") + .withDefault(true); + thermal_force = options["thermal_force"] + .doc("Include thermal force terms?") + .withDefault(true); + } + + virtual std::vector additionalComponents() override { + std::vector result = { + {"braginskii_collisions", "braginskii_collisions"}, + {"braginskii_friction", "braginskii_friction"}, + {"braginskii_heat_exchange", "braginskii_heat_exchange"}, + {"braginskii_conduction", "braginskii_conduction"}}; + if (electron_viscosity) { + result.emplace_back("braginskii_electron_viscosity", + "braginskii_electron_viscosity"); + } + if (ion_viscosity) { + result.emplace_back("braginskii_ion_viscosity", "braginskii_ion_viscosity"); + } + if (thermal_force) { + result.emplace_back("braginskii_thermal_force", "braginskii_thermal_force"); + } + return result; + } + + static constexpr auto type = "braginskii_closure"; + +private: + bool electron_viscosity; /// Whether to include electron viscosity terms + bool ion_viscosity; /// Whether to include ion viscosity terms + bool thermal_force; /// Whether to include thermal force terms + + /// Empty transform; all the work actually happens in the subcomponents. + void transform_impl(GuardedOptions&) override {} +}; + +namespace { +RegisterComponent registercomponentbraginskiiclosure; +} + +#endif // BRAGINSKII_H diff --git a/tests/unit/test_braginskii_closure.cxx b/tests/unit/test_braginskii_closure.cxx new file mode 100644 index 000000000..3f9da557f --- /dev/null +++ b/tests/unit/test_braginskii_closure.cxx @@ -0,0 +1,89 @@ + +#include "gtest/gtest.h" + +#include "fake_mesh_fixture.hxx" +#include "test_extras.hxx" // FakeMesh + +#include "../../include/braginskii_closure.hxx" + +/// Global mesh +namespace bout { +namespace globals { +extern Mesh* mesh; +} // namespace globals +} // namespace bout + +// The unit tests use the global mesh +using namespace bout::globals; + +// Reuse the "standard" fixture for FakeMesh +using BraginskiiClosureTest = FakeMeshFixture; + +std::set makeExpected(std::initializer_list names) { + std::set result; + for (const auto& name : names) { + result.emplace(name, name); + } + return result; +} + +std::set toSet(std::vector input) { + return std::set(input.begin(), input.end()); +} + +TEST_F(BraginskiiClosureTest, CreateDefault) { + Options options; + BraginskiiClosure component("test", options, nullptr); + EXPECT_EQ(toSet(component.additionalComponents()), + makeExpected({"braginskii_collisions", "braginskii_friction", + "braginskii_heat_exchange", "braginskii_conduction", + "braginskii_electron_viscosity", "braginskii_ion_viscosity", + "braginskii_thermal_force"})); +} + +TEST_F(BraginskiiClosureTest, CreateMinimal) { + Options options = {{"test", + {{"electron_viscosity", false}, + {"ion_viscosity", false}, + {"thermal_force", false}}}}; + BraginskiiClosure component("test", options, nullptr); + EXPECT_EQ(toSet(component.additionalComponents()), + makeExpected({"braginskii_collisions", "braginskii_friction", + "braginskii_heat_exchange", "braginskii_conduction"})); +} + +TEST_F(BraginskiiClosureTest, CreateThermalForce) { + Options options = {{"test", + {{"electron_viscosity", false}, + {"ion_viscosity", false}, + {"thermal_force", true}}}}; + BraginskiiClosure component("test", options, nullptr); + EXPECT_EQ(toSet(component.additionalComponents()), + makeExpected({"braginskii_collisions", "braginskii_friction", + "braginskii_heat_exchange", "braginskii_conduction", + "braginskii_thermal_force"})); +} + +TEST_F(BraginskiiClosureTest, CreateIonViscosity) { + Options options = {{"test", + {{"electron_viscosity", false}, + {"ion_viscosity", true}, + {"thermal_force", false}}}}; + BraginskiiClosure component("test", options, nullptr); + EXPECT_EQ(toSet(component.additionalComponents()), + makeExpected({"braginskii_collisions", "braginskii_friction", + "braginskii_heat_exchange", "braginskii_conduction", + "braginskii_ion_viscosity"})); +} + +TEST_F(BraginskiiClosureTest, CreateElectronViscosity) { + Options options = {{"test", + {{"electron_viscosity", true}, + {"ion_viscosity", false}, + {"thermal_force", false}}}}; + BraginskiiClosure component("test", options, nullptr); + EXPECT_EQ(toSet(component.additionalComponents()), + makeExpected({"braginskii_collisions", "braginskii_friction", + "braginskii_heat_exchange", "braginskii_conduction", + "braginskii_electron_viscosity"})); +} diff --git a/tests/unit/test_component_scheduler.cxx b/tests/unit/test_component_scheduler.cxx index 90d34248f..3f44819c3 100644 --- a/tests/unit/test_component_scheduler.cxx +++ b/tests/unit/test_component_scheduler.cxx @@ -30,17 +30,16 @@ struct TestMultiply : public NamedComponent { }; struct TestAdditionalComponent : public NamedComponent { - TestAdditionalComponent(const std::string&, Options&, Solver*) - : Component({}) {} + TestAdditionalComponent(const std::string& name, Options&, Solver*) + : NamedComponent(name, {}) {} - void transform_impl(GuardedOptions&) override { - } + void transform_impl(GuardedOptions&) override {} std::vector additionalComponents() override { - return {{"TestComponent", "testcomponent"}, {"component2", "multiply"}}; + return {{"TestComponent", "testcomponent"}, {"component2", "multiply"}}; } - static constexpr auto type = "testadditionalcomponent"; + static constexpr auto type = "additionalcomponent"; }; struct OrderChecker : public NamedComponent { From 7fa9766412074231a1fef5a233c100196c0654f7 Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Mon, 15 Dec 2025 14:19:38 +0000 Subject: [PATCH 3/7] Add documentation for creating sub-components. --- docs/sphinx/closure.rst | 5 +++++ docs/sphinx/developer.rst | 16 +++++++++++++++- 2 files changed, 20 insertions(+), 1 deletion(-) diff --git a/docs/sphinx/closure.rst b/docs/sphinx/closure.rst index 8794ed2f0..ffe5ff64f 100644 --- a/docs/sphinx/closure.rst +++ b/docs/sphinx/closure.rst @@ -40,6 +40,11 @@ is captured in the :ref:`sec-neutral_parallel_diffusion` top-level component, wh both parallel Braginskii transport and perpendicular pressure-diffusion for 2D/3D are captured in the :ref:`sec-neutral_mixed` species-level component. +A user can automatically activate all of these components at once +using the `BraginskiiClosure` component. + +.. doxygenclass:: BraginskiiClosure + :members: Collision frequency selection diff --git a/docs/sphinx/developer.rst b/docs/sphinx/developer.rst index 0534b616a..8407e2fb0 100644 --- a/docs/sphinx/developer.rst +++ b/docs/sphinx/developer.rst @@ -617,6 +617,17 @@ species and density of ions:: See the documentation for `Component::declareAllSpecies` for a list of all substitutions that will be performed. +Components may request the creation of additional components, upon +which they depend. This is done by overriding the +`Component::additionalComponents` method, which returns a list of +`ComponentInformation` structs that specify the names and types of +components required. These extra components will use the default +settings for components of this type, unless other values are +specified in a section of the input file with the component name. + +.. doxygenstruct:: ComponentInformation + :members: + Component scheduler ~~~~~~~~~~~~~~~~~~~ @@ -631,7 +642,10 @@ and then in `Hermes::rhs` the components are run by a call:: The call to `ComponentScheduler::create` treats the "components" option as a comma-separated list of names. For each name in the list, -the scheduler looks up the options under the section of that name. The +the scheduler looks up the options under the section of that name. If +any of the listed components request further components +(via the `Component::AdditionalComponents` method) then these will be +created too. The ``ComponentScheduler`` will use permission information stored by each component in `Component::state_variable_access` to work out the order to execute components. It will ensure that all writes to a variable From 07a5b1fa58bad73b50bd8bda78021731ca0bd5c4 Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Fri, 4 Sep 2026 18:14:49 +0100 Subject: [PATCH 4/7] Add ability to print component order and deactivate sorting --- include/component.hxx | 20 +++++++++++-- src/component.cxx | 6 ++++ src/component_scheduler.cxx | 58 ++++++++++++++++++++++++++++++------- 3 files changed, 71 insertions(+), 13 deletions(-) diff --git a/include/component.hxx b/include/component.hxx index c1a888898..fa683e16a 100644 --- a/include/component.hxx +++ b/include/component.hxx @@ -81,9 +81,13 @@ struct ComponentInformation { std::string name; std::string type; - ComponentInformation(const std::string& name_, const std::string& type_) : name(name_), type(type_) {} + ComponentInformation() {}; - ComponentInformation(std::string&& name_, std::string&& type_) : name(std::move(name_)), type(std::move(type_)) {} + ComponentInformation(const std::string& name_, const std::string& type_) + : name(name_), type(type_) {} + + ComponentInformation(std::string&& name_, std::string&& type_) + : name(std::move(name_)), type(std::move(type_)) {} bool operator<(const ComponentInformation& other) const { return std::pair(name, type) < std::pair(other.name, other.type); @@ -94,6 +98,18 @@ struct ComponentInformation { } }; +/// Format `ComponentInformation` to string. Format string specification is the +/// same as when formatting a string. +/// See https://fmt.dev/12.0/syntax/#format-specification-mini-language. +/// +/// TODO: provide custom formatting to configure exactly how the +/// component name and type are displayed. +template <> +struct fmt::formatter : formatter { + auto format(const ComponentInformation& ci, format_context& ctx) const + -> format_context::iterator; +}; + /// Interface for a component of a simulation model /// /// The constructor of derived types should have signature diff --git a/src/component.cxx b/src/component.cxx index cd523e286..8e63ae111 100644 --- a/src/component.cxx +++ b/src/component.cxx @@ -9,6 +9,12 @@ #include "../include/guarded_options.hxx" #include "../include/permissions.hxx" +auto fmt::formatter::format(const ComponentInformation& ci, + format_context& ctx) const + -> format_context::iterator { + return formatter::format(fmt::format("{} ({})", ci.name, ci.type), ctx); +} + std::unique_ptr Component::create(const std::string& type, const std::string& name, Options& alloptions, Solver* solver) { diff --git a/src/component_scheduler.cxx b/src/component_scheduler.cxx index 1f3c7c6e8..676ca75ca 100644 --- a/src/component_scheduler.cxx +++ b/src/component_scheduler.cxx @@ -11,6 +11,7 @@ #include #include #include +#include #include // for trim, strsplit #include #include @@ -308,6 +309,15 @@ setReadDependencies(const std::vector>& components, return missing; } +void printComponents(const std::vector>& components) { + if (!components.empty()) { + output_info << "Components will be executed in the following order:\n"; + } + for (const auto& comp : components) { + output_info << fmt::format("\t{}\n", *comp); + } +} + /// Topologically sorts the list of components to ensure variables are /// written and read in the right order. /// @@ -425,6 +435,7 @@ void sortComponents(std::vector>& components) { } components = std::move(result); + printComponents(components); } } // namespace @@ -434,8 +445,12 @@ ComponentScheduler::ComponentScheduler(Options& scheduler_options, const std::string component_names = scheduler_options["components"] .doc("Components in order of execution") .as(); + const bool autosort = scheduler_options["autosort"] + .doc("Perform a topological sort to ensure components " + "executed in the right order?") + .withDefault(true); - std::set required_components; + std::list required_components; std::vector electrons; std::vector neutrals; @@ -485,22 +500,39 @@ ComponentScheduler::ComponentScheduler(Options& scheduler_options, continue; } - required_components.emplace(name_trimmed, type_trimmed); + required_components.emplace_back(name_trimmed, type_trimmed); } } + // Use sets for efficient lookup of whether a component is already + // in the queue to be created or already has been created. + // + // We use a vector to decide the order in which to create + // components, to keep this close to what is in the file. + std::set unbuilt_components(required_components.begin(), + required_components.end()); std::set created_components; - for (auto it = required_components.begin(); it != required_components.end(); - it = required_components.begin()) { - auto comp = Component::create(it->type, it->name, component_options, solver); - for (const auto& sub_comp : comp->additionalComponents()) { - if (required_components.count(sub_comp) == 0 - and created_components.count(sub_comp) == 0) { - required_components.insert(sub_comp); + // We use a vector to keep track of the actual order in which + // components are created + std::vector component_order; + + while (required_components.size() > 0) { + const ComponentInformation component = required_components.front(); + required_components.pop_front(); + auto comp = + Component::create(component.type, component.name, component_options, solver); + std::vector sub_components = comp->additionalComponents(); + for (auto sub_comp = sub_components.rbegin(); sub_comp != sub_components.rend(); + ++sub_comp) { + if (unbuilt_components.count(*sub_comp) == 0 + and created_components.count(*sub_comp) == 0) { + required_components.push_front(*sub_comp); + unbuilt_components.insert(*sub_comp); } } + component_order.push_back(component); components.push_back(std::move(comp)); - created_components.insert(required_components.extract(it)); + created_components.insert(unbuilt_components.extract(component)); } const SpeciesInformation species(electrons, neutrals, positive_ions, negative_ions); @@ -509,7 +541,11 @@ ComponentScheduler::ComponentScheduler(Options& scheduler_options, component->declareAllSpecies(species); } - ::sortComponents(components); + if (autosort) { + ::sortComponents(components); + } else { + printComponents(components); + } } std::unique_ptr ComponentScheduler::create(Options& scheduler_options, From b636bf4afa07770314231d184f7d68d19498a7e5 Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Fri, 4 Sep 2026 18:21:43 +0100 Subject: [PATCH 5/7] Tidy up permission --- src/braginskii_collisions.cxx | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/braginskii_collisions.cxx b/src/braginskii_collisions.cxx index 3168b61c5..5321ccada 100644 --- a/src/braginskii_collisions.cxx +++ b/src/braginskii_collisions.cxx @@ -25,7 +25,8 @@ BraginskiiCollisions::BraginskiiCollisions(const std::string& name, Options& alloptions, Solver*) : NamedComponent(name, - {readOnly("species:{non_electrons}:density", Regions::Interior), + {readOnly("species:{all_species}:density", Regions::Interior), + readOnly("species:{electrons}:temperature", Regions::Interior), readIfSet("species:{non_electrons}:charge"), readIfSet("species:{negative_ions}:temperature", Regions::Interior), readOnly("species:{all_species}:AA")}) { @@ -65,8 +66,6 @@ BraginskiiCollisions::BraginskiiCollisions(const std::string& name, Options& all diagnose = options["diagnose"].doc("Output additional diagnostics?").withDefault(false); - setPermissions(readOnly("species:{electrons}:temperature", Regions::Interior)); - setPermissions(readOnly("species:{electrons}:density", Regions::Interior)); if (electron_electron) { setPermissions(readWrite( "species:{electrons}:collision_frequencies:{electrons}_{electrons2}_coll")); From 84f53b4499023af0f929433723ab5657a0300f15 Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Fri, 4 Sep 2026 18:40:20 +0100 Subject: [PATCH 6/7] Debugging ordering issue with upstream_density and viscosity --- src/braginskii_ion_viscosity.cxx | 21 +++++++++++++++------ 1 file changed, 15 insertions(+), 6 deletions(-) diff --git a/src/braginskii_ion_viscosity.cxx b/src/braginskii_ion_viscosity.cxx index b104a3184..3c4b2a18a 100644 --- a/src/braginskii_ion_viscosity.cxx +++ b/src/braginskii_ion_viscosity.cxx @@ -28,6 +28,8 @@ #include "../include/div_ops.hxx" #include "../include/hermes_utils.hxx" +// FIXME: Why does it matter if upstream_density_feedback executed after this? Is there some wierd interaction of the readIfSet variables? + using bout::globals::mesh; BraginskiiIonViscosity::BraginskiiIonViscosity(const std::string& name, @@ -36,10 +38,11 @@ BraginskiiIonViscosity::BraginskiiIonViscosity(const std::string& name, name, { readIfSet("species:{non_electrons}:pressure"), - readIfSet("species:{non_electrons}:temperature"), - readIfSet("species:{non_electrons}:density"), - readIfSet("species:{non_electrons}:velocity"), readIfSet("species:{non_electrons}:charge"), + // FIXME: The rest of these (except DivJextra) only + // apply for species in which the above are set + readOnly("species:{non_electrons}:temperature"), + readOnly("species:{non_electrons}:density"), readIfSet("species:{non_electrons}:collision_frequencies:{coll_type}"), readWrite("species:{non_electrons}:momentum_source"), readWrite("species:{non_electrons}:energy_source"), @@ -121,11 +124,20 @@ BraginskiiIonViscosity::BraginskiiIonViscosity(const std::string& name, Curlb_B.z *= SQ(Lnorm); Curlb_B *= 2. / coord->Bxy; + + setPermissions(readOnly("fields:phi")); + } + + if (parallel) { + setPermissions(readOnly("species:{non_electrons}:velocity")); } + if (bounce_frequency) { const Options& units = alloptions["units"]; const BoutReal Lnorm = units["meters"]; bounce_frequency_R /= Lnorm; + + setPermissions(readOnly("species:{non_electrons}:pressure")); } std::vector coll_types; @@ -135,9 +147,6 @@ BraginskiiIonViscosity::BraginskiiIonViscosity(const std::string& name, coll_types.push_back("{non_electrons}_{all_species}_coll"); coll_types.push_back("{non_electrons}_{all_species}_cx"); } - if (perpendicular) { - setPermissions(readOnly("fields:phi")); - } substitutePermissions("coll_type", coll_types); } From f82ce21003e0b76f97b06a1457c891f0d175845c Mon Sep 17 00:00:00 2001 From: Chris MacMackin Date: Tue, 8 Sep 2026 10:50:09 +0100 Subject: [PATCH 7/7] Fix permission error for electron momentum source --- include/electron_force_balance.hxx | 18 ++++++++---------- .../integrated/1D-recycling-dthe/data/BOUT.inp | 5 +++-- 2 files changed, 11 insertions(+), 12 deletions(-) diff --git a/include/electron_force_balance.hxx b/include/electron_force_balance.hxx index e96273d23..815e8384f 100644 --- a/include/electron_force_balance.hxx +++ b/include/electron_force_balance.hxx @@ -20,16 +20,14 @@ /// struct ElectronForceBalance : public NamedComponent { ElectronForceBalance(std::string name, Options& alloptions, Solver*) - : NamedComponent(name, - {readOnly("species:e:pressure"), - readOnly("species:e:density", Regions::Interior), - readOnly("species:e:charge"), - // FIXME: Only writes if already exists - readWrite("species:e:momentum_source"), - readIfSet("species:{non_electrons}:density", Regions::Interior), - readIfSet("species:{non_electrons}:charge"), - // FIXME: Only written if density and charge have been set. - readWrite("species:{non_electrons}:momentum_source")}) { + : NamedComponent( + name, {readOnly("species:e:pressure"), + readOnly("species:e:density", Regions::Interior), + readOnly("species:e:charge"), readIfSet("species:e:momentum_source"), + readIfSet("species:{non_electrons}:density", Regions::Interior), + readIfSet("species:{non_electrons}:charge"), + // FIXME: Only written if density and charge have been set. + readWrite("species:{non_electrons}:momentum_source")}) { auto& options = alloptions[name]; diagnose = options["diagnose"] .doc("Save additional output diagnostics") diff --git a/tests/integrated/1D-recycling-dthe/data/BOUT.inp b/tests/integrated/1D-recycling-dthe/data/BOUT.inp index 55df01856..e02a5bd24 100644 --- a/tests/integrated/1D-recycling-dthe/data/BOUT.inp +++ b/tests/integrated/1D-recycling-dthe/data/BOUT.inp @@ -44,9 +44,10 @@ ixseps2 = -1 # Evolve ion density, ion and electron pressure, then calculate force on ions due # to electron pressure by using electron force balance. components = (d+, d, t+, t, he+, he, e, - sheath_boundary, braginskii_collisions, braginskii_friction, braginskii_heat_exchange, + sheath_boundary, braginskii_collisions, braginskii_heat_exchange, recycling, reactions, braginskii_thermal_force, electron_force_balance, - neutral_parallel_diffusion, braginskii_ion_viscosity, braginskii_conduction) + neutral_parallel_diffusion, braginskii_ion_viscosity, braginskii_conduction, + braginskii_friction) Nnorm = 1e19 Bnorm = 1