From a076c628b8f3c085c0f0309bdafab8cf1b8685e0 Mon Sep 17 00:00:00 2001 From: Xander <102053371+fractalsbyx@users.noreply.github.com> Date: Fri, 20 Dec 2024 18:43:44 -0500 Subject: [PATCH] added more intuitive insertion functions for equation dependencies (#371) now the old functions just call the new ones --- include/core/variableAttributeLoader.h | 64 ++++++++++++++++-- src/core/variableAttributeLoader.cc | 91 ++++++++++++++++++-------- 2 files changed, 121 insertions(+), 34 deletions(-) diff --git a/include/core/variableAttributeLoader.h b/include/core/variableAttributeLoader.h index 4f05e627..503b7605 100644 --- a/include/core/variableAttributeLoader.h +++ b/include/core/variableAttributeLoader.h @@ -63,7 +63,7 @@ class variableAttributeLoader set_variable_equation_type(const unsigned int &index, const PDEType &var_eq_type) const; /** - * \brief Set the dependencies for the value term of the RHS equation of the variable at + * \brief Add dependencies for the value term of the RHS equation of the variable at * `index`. * * \param index Index of variable @@ -75,7 +75,7 @@ class variableAttributeLoader const std::string &dependencies); /** - * \brief Set the dependencies for the gradient term of the RHS equation of the variable + * \brief Add dependencies for the gradient term of the RHS equation of the variable * at `index`. * * \param index Index of variable @@ -87,7 +87,7 @@ class variableAttributeLoader const std::string &dependencies); /** - * \brief Set the dependencies for the value term of the LHS equation of the variable at + * \brief Add dependencies for the value term of the LHS equation of the variable at * `index`. * * \param index Index of variable @@ -99,7 +99,7 @@ class variableAttributeLoader const std::string &dependencies) const; /** - * \brief Set the dependencies for the gradient term of the LHS equation of the variable + * \brief Add dependencies for the gradient term of the LHS equation of the variable * at `index`. * * \param index Index of variable @@ -110,6 +110,58 @@ class variableAttributeLoader set_dependencies_gradient_term_LHS(const unsigned int &index, const std::string &dependencies) const; + /** + * \brief Insert dependencies for the value term of the RHS equation of the variable at + * `index`. + * + * \param index Index of variable + * \param dependencies Container containing list of dependency strings for + * variable at `index` Hint: {"variable", "grad(variable)", "hess(variable)"} + */ + template + void + insert_dependencies_value_term_RHS(const unsigned int &index, + const Iterable &dependencies); + + /** + * \brief Insert dependencies for the gradient term of the RHS equation of the variable + * at `index`. + * + * \param index Index of variable + * \param dependencies Container containing list of dependency strings for + * variable at `index` Hint: {"variable", "grad(variable)", "hess(variable)"} + */ + template + void + insert_dependencies_gradient_term_RHS(const unsigned int &index, + const Iterable &dependencies); + + /** + * \brief Insert dependencies for the value term of the LHS equation of the variable at + * `index`. + * + * \param index Index of variable + * \param dependencies Container containing list of dependency strings for + * variable at `index` Hint: {"variable", "grad(variable)", "hess(variable)"} + */ + template + void + insert_dependencies_value_term_LHS(const unsigned int &index, + const Iterable &dependencies) const; + + /** + * \brief Insert dependencies for the gradient term of the LHS equation of the variable + * at `index`. + * + * \param index Index of variable + * \param dependencies Container containing list of dependency strings for + * variable at `index` Hint: {"variable", "grad(variable)", "hess(variable)"} + */ + template + void + insert_dependencies_gradient_term_LHS(const unsigned int &index, + const Iterable &dependencies) const; + /** * \brief Flag whether the variable at `index` is needed to calculate the nucleation * probability. @@ -143,13 +195,13 @@ class variableAttributeLoader /** * \brief getter function for variable attributes list (copy). */ - AttributesList + [[nodiscard]] AttributesList get_var_attributes() const; /** * \brief getter function for postprocessing attributes list (copy). */ - AttributesList + [[nodiscard]] AttributesList get_pp_attributes() const; /** diff --git a/src/core/variableAttributeLoader.cc b/src/core/variableAttributeLoader.cc index 152cc430..9b42ac6c 100644 --- a/src/core/variableAttributeLoader.cc +++ b/src/core/variableAttributeLoader.cc @@ -109,18 +109,7 @@ variableAttributeLoader::set_dependencies_value_term_RHS(const unsigned int &ind { std::vector dependencies_set = dealii::Utilities::split_string_list(strip_whitespace(dependencies)); - /* (*relevant_attributes)[index].dependencies_value_RHS = - std::set(dependencies_set.begin(), dependencies_set.end()); */ - if (relevant_attributes != &pp_attributes) - { - var_attributes[index].dependencies_value_RHS = - std::set(dependencies_set.begin(), dependencies_set.end()); - } - else - { - pp_attributes[index].dependencies_value_PP = - std::set(dependencies_set.begin(), dependencies_set.end()); - } + insert_dependencies_value_term_RHS(index, dependencies_set); } void @@ -130,18 +119,7 @@ variableAttributeLoader::set_dependencies_gradient_term_RHS( { std::vector dependencies_set = dealii::Utilities::split_string_list(strip_whitespace(dependencies)); - /* (*relevant_attributes)[index].dependencies_gradient_RHS = - std::set(dependencies_set.begin(), dependencies_set.end()); */ - if (relevant_attributes != &pp_attributes) - { - var_attributes[index].dependencies_gradient_RHS = - std::set(dependencies_set.begin(), dependencies_set.end()); - } - else - { - pp_attributes[index].dependencies_gradient_PP = - std::set(dependencies_set.begin(), dependencies_set.end()); - } + insert_dependencies_gradient_term_RHS(index, dependencies_set); } void @@ -151,8 +129,7 @@ variableAttributeLoader::set_dependencies_value_term_LHS( { std::vector dependencies_set = dealii::Utilities::split_string_list(strip_whitespace(dependencies)); - (*relevant_attributes)[index].dependencies_value_LHS = - std::set(dependencies_set.begin(), dependencies_set.end()); + insert_dependencies_value_term_LHS(index, dependencies_set); } void @@ -162,8 +139,66 @@ variableAttributeLoader::set_dependencies_gradient_term_LHS( { std::vector dependencies_set = dealii::Utilities::split_string_list(strip_whitespace(dependencies)); - (*relevant_attributes)[index].dependencies_gradient_LHS = - std::set(dependencies_set.begin(), dependencies_set.end()); + insert_dependencies_gradient_term_LHS(index, dependencies_set); +} + +template +void +variableAttributeLoader::insert_dependencies_value_term_RHS(const unsigned int &index, + const Iterable &dependencies) +{ + /* (*relevant_attributes)[index].dependencies_value_RHS.insert(dependencies.begin(), + * dependencies.end()); */ + if (relevant_attributes != &pp_attributes) + { + var_attributes[index].dependencies_value_RHS.insert(dependencies.begin(), + dependencies.end()); + } + else + { + pp_attributes[index].dependencies_value_PP.insert(dependencies.begin(), + dependencies.end()); + } +} + +template +void +variableAttributeLoader::insert_dependencies_gradient_term_RHS( + const unsigned int &index, + const Iterable &dependencies) +{ + /* (*relevant_attributes)[index].dependencies_gradient_RHS.insert(dependencies.begin(), + * dependencies.end()); */ + if (relevant_attributes != &pp_attributes) + { + var_attributes[index].dependencies_gradient_RHS.insert(dependencies.begin(), + dependencies.end()); + } + else + { + pp_attributes[index].dependencies_gradient_PP.insert(dependencies.begin(), + dependencies.end()); + } +} + +template +void +variableAttributeLoader::insert_dependencies_value_term_LHS( + const unsigned int &index, + const Iterable &dependencies) const +{ + (*relevant_attributes)[index].dependencies_value_LHS.insert(dependencies.begin(), + dependencies.end()); +} + +template +void +variableAttributeLoader::insert_dependencies_gradient_term_LHS( + const unsigned int &index, + const Iterable &dependencies) const +{ + (*relevant_attributes)[index].dependencies_gradient_LHS.insert(dependencies.begin(), + dependencies.end()); } void